Module lattice 

Lattice sums of outgoing spherical and cylindrical waves, and lattice geometry.

Upstream: treams.lattice (lattice/_esum.pyx, _dsum.pyx and _misc.pyx of tfp-photonics/treams, MIT), whose formulae and conventions the sums follow, and the Brillouin-zone helpers of treams.misc.

Definition

A lattice of dimension d holds the points R = n_1 a_1 + ... + n_d a_d with integer n_i. For a wavenumber k, a real Bloch vector kpar (the wavevector components along the lattice) and a shift r, sum returns

D_lm(k, kpar, r) = Σ'_R h_l^(1)(k |r + R|) Y_lm(-r - R) exp(i kpar · R)        spherical
D_m(k, kpar, r)  = Σ'_R H_m^(1)(k |r + R|) exp(i m φ(-r - R)) exp(i kpar · R)   cylindrical

as in treams.lattice. Y_lm(-r - R) is the spherical harmonic in the direction of -r - R and φ its azimuth. The prime drops the point with r + R = 0, if there is one. Family selects the sum and its labels: degree l and order m for spherical sums (d <= 3), order m for cylindrical sums (d <= 2).

Frame

The lattice lies along fixed Cartesian axes:

Sum	Lattice axes
1D spherical	z
1D cylindrical	x
2D spherical and cylindrical	x and y (the plane z = 0)
3D spherical	x, y and z

BlochLattice takes the lattice vectors and kpar as components along these axes. The shift is a Cartesian 3-vector; cylindrical shifts lie in the xy plane.

Branch and Bloch sign

The waves are outgoing: h_l^(1) and H_m^(1), which decay for Im k > 0. The Bloch factor is exp(+i kpar · R), so the sum at r + R is exp(-i kpar · R) times the sum at r. The sum depends on the lattice and on kpar modulo the reciprocal lattice, but not on the primitive basis or on the Bloch cell. Exact diffraction thresholds (|kpar + G| = k for a reciprocal lattice vector G) are singular.

Ewald split

The Ewald split parameter eta divides each sum into a real-space part, whose terms fall like exp(-Re((k eta)^2) |R|^2 / 2), and a reciprocal-space part, whose terms fall like exp(-Re(1 / (k eta)^2) |kpar + G|^2 / 2). The exact sum does not depend on eta. eta = 0 selects the automatic split

eta = sqrt(2 pi) / (k L) max(|k L| / 8, 1),    L = V^(1/d)

with the cell measure V (length, area or volume). Derivatives hold the split fixed. The automatic split has the modulus sqrt(2 pi) / 8 = 0.31 from |k L| = 8 on, and more below. A split is small when it lies below every automatic one, Re(1 / (2 eta^2)) > 16 / pi: there the two parts grow like exp(Re(1 / (2 eta^2))) and cancel.

Ewald sums require Im k >= 0, and spherical sums additionally require Re k >= 0. Within this domain every sum has the same value at every split where its Ewald parts converge. Inputs outside it give InvalidInput: the gain-side continuation is not supported, and negative-real spherical wavenumbers select the wrong Kambe branch. Finite direct shells accept every finite nonzero complex wavenumber.

Crosswalk

Each treams.lattice function is one call of sum_part (or sum for Full). r is the shift, i the shell of a direct sum, and dim the lattice dimension of BlochLattice.

treams	Family	dim	Shift
lsumsw1d(l, k, kpar, a, r, eta)	Spherical { l, m: 0 }	1	[0, 0, r]
lsumsw1d_shift(l, m, ...)	Spherical { l, m }	1	r
lsumsw2d(l, m, ...)	Spherical { l, m }	2	[r_x, r_y, 0]
lsumsw2d_shift(l, m, ...)	Spherical { l, m }	2	r
lsumsw3d(l, m, ...)	Spherical { l, m }	3	r
lsumcw1d(m, ...)	Cylindrical { m }	1	[r, 0, 0]
lsumcw1d_shift(m, ...)	Cylindrical { m }	1	[r_x, r_y, 0]
lsumcw2d(m, ...)	Cylindrical { m }	2	[r_x, r_y, 0]

The prefix selects the SumPart:

treams prefix	SumPart
lsum*	Full
realsum*	Real
recsum*, with zero2d or zero3d added at r = 0	Reciprocal, which adds the self term at r = 0 for every degree
dsum*(..., i)	Direct(i); eta is not used

The geometry helpers map one to one:

treams	treams-rs
lattice.volume, lattice.area	volume
lattice.reciprocal	reciprocal()
lattice.cube(d, n), lattice.cubeedge(d, n)	cube with edge false or true
lattice.diffr_orders_circle(b, rmax)	diffraction_orders
misc.firstbrillouin1d	first_brillouin_1d
misc.firstbrillouin2d, misc.firstbrillouin3d	first_brillouin

derivatives, derivatives_part and SumResidual have no upstream counterpart: they give the derivatives and gradients of the sums.

Differences from treams.lattice

• Automatic split in 3D. The cell length is L = V^(1/3). Upstream _check_eta divides by k V^(2/3). Both choices give the same sum wherever the parts converge; they differ in cost and in which sums converge.
• Convergence. Each Ewald part adds integer cube shells until two consecutive shells add less than 2e-13 max(|P|, 1) of the part P, and fails after 200, 32 or 16 shells in one, two or three dimensions (up to four times that where a large |k| L or a small |k eta| needs it). Upstream stops once two shells change a part by less than 1e-10, or after 200, 20 or 10 shells, and reports no failure.
• Reduced bases. Both Ewald parts sum over shells of reduced bases, the shift is reduced into a cell of the direct lattice and kpar into a cell of the reciprocal lattice. Upstream sums over the given basis, whose shells can miss short lattice vectors of a skewed basis.
• Outgoing sheets. For Im k >= 0 the reciprocal integrals at a complex split take the sheet that continues the sum from the automatic split, so the value does not depend on the split. Upstream takes the principal branches.
• Small splits. Sums at small splits predict the rounding their cancelling parts leave and fail where it exceeds 1e-3 of the scale of a component. Upstream returns the cancelled sum.
• Exact zeros. For odd l (spherical) or odd m (cylindrical), kpar = 0 and a shift on the lattice line or plane with 2 r a lattice vector, inversion maps the sum onto its negative. The sum is then exactly zero at every split. Upstream returns exact zeros only at r = 0 (and, in 1D, at r = ±a / 2); at the other such shifts it returns the rounding of its parts.
• Spectral series. Off the axis the Ewald parts of 1D spherical sums cancel by up to exp(w^2) (see the glossary). There each component whose rounding exceeds 1e-12 max(|S|, 1) comes from the spectral series, a sum over diffraction orders of cylindrical waves, where that is more accurate.
• Derivatives. derivatives returns the analytic derivatives in k, eta, the shift, kpar and the lattice vectors, from jets (see the glossary) and the derivative identities of the integrals.

Glossary

Term	Meaning
shell	The integer points n with max_i |n_i| = s of shell s; sums add shell by shell
split	The Ewald split parameter eta
small split	A split below every automatic one, Re(1 / (2 eta^2)) > 16 / pi
diffraction order	q = kpar + G for a reciprocal lattice vector G
v	(q · q / k^2 - 1) / (2 eta^2) of a diffraction order q
w	k rho eta, with the distance rho of the shift from the lattice plane or axis
t	w^2
F_n(v, t)	The reduced integral (-1)^n ∫_1^∞ s^(n-1) exp(-v s - t / (2 s)) ds of the reciprocal terms; F_n(v, 0) = Γ(n, v) / (-v)^n (upstream _redincgamma, _redintkambe)
I_n(z, eta)	The Kambe integral ∫_eta^∞ t^n exp(-z^2 t^2 / 2 + 1 / (2 t^2)) dt of the real-space terms, at z = k |r + R| (upstream intkambe)
sheet	The branch of ln v on which a reciprocal term takes F_n
self term	The correction at r = 0 that removes the missing R = 0 term from the reciprocal part; upstream adds it as zero2d or zero3d
jet	A value carried through the arithmetic together with its first derivatives
settled	A real-space shell sum that reaches its shell limit while the Bloch phases cancel each of its last two shells; it is kept only if it agrees with the sum at the automatic split

Failure messages

Message	Error	Cause
Ewald sums require Im(k) >= 0; gain media are unsupported	InvalidInput	An amplifying wavenumber
spherical Ewald sums require Re(k) >= 0	InvalidInput	A spherical wavenumber on the unsupported Kambe branch
Ewald sum did not converge within the shell limit	NotConverged	A part needs more shells than its limit
… ; use a larger split (eta = 0 selects one)	NotConverged	The same, at a small split
… , neither at this split nor at the automatic one	NotConverged	A settled sum whose check at the automatic split does not converge
Ewald split too small: its cancelling parts lose about ...; use a larger split (eta = 0 selects one)	NotConverged	The parts at a small split cancel by more than about 1e-3 of the scale of the value or a derivative
Ewald sum lost its accuracy to cancelling Kambe integrals; reduce the split parameter	NotConverged	A 1D spherical sum off the axis rounds by more than 1e-3 max(|S|, 1), and its spectral series does not help
non-finite Ewald summand; change the split parameter or reduce the order	NonFinite	A term is not finite
non-finite direct lattice term	NonFinite	A term of a direct sum is not finite
lattice sum is at a diffraction threshold; supply a limiting complex wavenumber	InvalidInput	A diffraction order lies on its threshold, |v| <= 1e-14
component adjoints require an explicit nonzero Ewald split	InvalidInput	Derivatives of the real or reciprocal part at eta = 0
direct half-cell shell grouping changes discontinuously; use the Ewald sum for adjoints	InvalidInput	Derivatives of a direct 1D shell at a half-cell shift

Invalid labels, dimensions and non-finite inputs give InvalidInput with a message that names the input. A value-only call can return a sum whose derivatives fail.

Tuning constants

Name	Value	Meaning	Why
resolve_split	1 / 8, sqrt(2 pi)	Automatic split sqrt(2 pi) / (k L) max(|k L| / 8, 1)	As upstream _check_eta, except for L in 3D
SHELL_TOLERANCE	2e-13	Part of max(|P|, 1) that each of the last two shells may add	About 900 ulps of the scale
shells	200, 32, 16	Shell limit in 1D, 2D, 3D, raised up to four times where the terms need it (shells_within)	200 in 1D as upstream; above upstream’s 20 and 10 in 2D and 3D
TERM_DECAY	50	Shell limits reach every term above exp(-50)	exp(-50) = 2e-22, far below SHELL_TOLERANCE
AUTOMATIC_GROWTH	16 / pi (1 + 1e-12)	Largest Re(1 / (2 eta^2)) of the automatic split; larger values mark a small split	Exact, with a margin for rounding
MAX_LOSS	1e-3	Largest predicted loss of a sum, relative to the scale of a component, before it fails	Calibrated: the bounds exceed the errors by a median factor of 17 to 87
EARLY_FAILURE_LOSS	2e-3	Loss at which a complete sum at a small split fails before its far shells	Calibrated: without it, a 3D sum at the split 0.03 takes 29 s to fail
TERM_ULPS	8	Rounding of each Ewald term and of the self term, in ulps	The MAX_LOSS calibration compares the bounds built from it with the errors
REAL_ULPS	64	Rounding of the real-space part of a 1D spherical sum, in ulps	Its Kambe integrals of even order keep their relative accuracy
SPECTRAL_SWITCH_LOSS	1e-12	Rounding of a 1D spherical component, relative to its scale, from which the spectral series replaces it	Far below MAX_LOSS, so the series takes over long before a sum fails
SPECTRAL_FIRST_W	2.5	w = |k| rho |eta| from which a 1D spherical sum tries the spectral series first	The reciprocal part cancels by exp(w^2) = 5e2 there
SERIES_SHELLS	512	Diffraction orders the spectral series may take on each side	Caps the cost of one series; beyond it the Ewald sum is used
SERIES_DECAY	45	The spectral series runs only where its terms fall by exp(-45) within SERIES_SHELLS	exp(-45) = 3e-20, below QUIET_SHELL
QUIET_SHELL	1e-18	Part of the summed moduli below which a shell of the spectral series is quiet	The series stops after two quiet shells past its peak; 1e-18 is below the rounding of the sum
SERIES_T	1e-2	Largest |t| at which integer orders of F_n come from their series in t	Below about |w| = 1e-50 the Kambe integrals at eta = -i / w leave the range of f64; the series loses only a few ulps
SERIES_REACH	1	The same, for the orders below the Kambe base pair	Each of these orders would need a gamma-function series of its own
HALF_INTEGER_SERIES_T	1	Largest |t| at which half-integer orders come from their series in t	Calibrated: within 1e-13 of 30-digit references for |t| from 0.3 to 1
SERIES_TERMS	64	Terms of the series in t at most	|t| <= 1 stops within 20
NEAR_THRESHOLD	1 / 16	Bound on the real and imaginary parts of q · q / k^2 - 1 below which the zeroth order takes compensated sums	Rounded sums lose 50 ulps of v there
THRESHOLD_DISTANCE	1e-14	|v| (reciprocal terms) or |k_q^2 / k^2| (spectral series) at or below which a diffraction order counts as on its threshold	The sums are singular on the threshold |q| = k
SHEET_MARGIN	0.5	Angle margin within which a diffraction order needs the argument of its sheet	Far above the rounding of the angle (1e-2 next to a threshold)
axial_cylindrical_reciprocal	4	|v| above which the gamma functions of value-only 1D cylindrical sums on the axis recur upward	A downward step there subtracts the leading asymptotic term
kambe_arguments	1e-100	Imaginary part that puts a real Kambe argument on the outgoing side	Any positive part selects the upper half plane
reduce_basis	1e-8, 2^20, 64	Relative shortening a reduction step needs, largest coefficient, most passes	Only strict decreases count, so ties such as those of a hexagonal basis keep it as given
first_brillouin	1e-13, 1e-14	Relative and absolute tolerance on squared lengths of a reduced 2D basis	The absolute one is upstream’s
PARALLEL	8	Batched sums run in parallel from this many outputs	Fewer outputs run on the calling thread

References

• K. Kambe, Z. Naturforsch. A 22 (1967), https://doi.org/10.1515/zna-1967-0305
• K. Kambe, Z. Naturforsch. A 23 (1968), https://doi.org/10.1515/zna-1968-0908
• C. M. Linton, SIAM Rev. 52, 630 (2010), https://doi.org/10.1137/09075130X
• D. Beutel et al., J. Opt. Soc. Am. B 38 (2021), https://doi.org/10.1364/JOSAB.419645
• NIST DLMF, chapter 8 (incomplete gamma functions), https://dlmf.nist.gov/8

Modules

accuracy 🔒

Rounding-error bounds, and the failure of sums whose Ewald parts cancel.

batch 🔒

Lattice sums over arrays of inputs, with their gradients.

cell 🔒

The lattice and Bloch vector of a sum, and the automatic Ewald split.

direct 🔒

Direct sums over cube shells, without Ewald acceleration.

ewald 🔒

The Ewald sum: real-space part, reciprocal-space part and self term.

geometry 🔒

Lattice geometry: cell volumes, reciprocal and reduced bases, integer cube points, diffraction orders and reduction to the first Brillouin zone.

inputs 🔒

The inputs of one sum as forward-mode jets, and their jet slots.

real 🔒

Real-space terms of the Ewald sum.

reciprocal 🔒

Reciprocal-space terms of the five sums (family and lattice dimension: sw3d, sw2d, sw1d, cw2d, cw1d) and the self term.

reduced 🔒

Reduced integrals F_n(v, t) of the reciprocal terms.

saved 🔒

Canonical numerical lattice geometry, including optional reduced bases. treams-rs extension.

sheets 🔒

Continuation of the reciprocal integrals to the sheet of the outgoing sum.

shells 🔒

Sums over cube shells and when they stop.

spectral 🔒

The spectral series of 1D spherical sums off the axis, and when to take it.

wave 🔒

The family of a sum, the scalar wave it adds up, and its angular factors.

Structs

BlochLattice

Row lattice vectors, a real Bloch vector kpar and their reduced bases in one, two or three dimensions.

Derivatives

Local derivatives of one scalar Ewald sum.

SumGradient

Cotangents of the continuous Ewald parameters under the real Hermitian pairing.

SumResidual

What sum_array saves for a batch of lattice sums: the original inputs. Every input holds one shared value or one value per output.

SumTangent

A direction in the continuous inputs of a lattice sum.

Enums

Family

The family of a lattice sum: spherical or cylindrical.

SumPart

Part of an Ewald sum, or one unaccelerated integer cube shell: the lsum*, realsum*, recsum* and dsum* functions of treams.lattice.

Functions

cube

All integer cube points, or just its boundary, in lexicographic row order.

derivatives

Differentiate an Ewald sum using analytic integral identities and Cartesian polynomials.

derivatives_part

Differentiate the selected component, holding a direct shell’s index fixed.

diffraction_orders

Diffraction orders inside a reciprocal-space circle, with adjacent opposite pairs.

evaluate 🔒

Validates the inputs and evaluates what evaluation selects: a direct shell, or the Ewald sum at the resolved split. N is the number of jet slots: 0 for a value, 1 for the split derivative of the real-space part, and for all derivatives 6, 10 or 16 in one, two or three dimensions (k, three shift components, dim Bloch components and dim * dim lattice-vector components).

evaluate_shared 🔒

Evaluates spherical orders sharing one degree’s real-space samples. Other callers pass no cache and retain the ordinary scalar path.

first_brillouin

Iterative nearest-cell reduction of a 2D or 3D wavevector.

first_brillouin_1d

Reduce a scalar wavevector to (-b/2, b/2], for positive reciprocal pitch b.

larger_split 🔒

error, advising a larger split where a sum at a split eta below every automatic one does not converge: its real-space shells converge sooner at a larger one. A sum whose shells settle but whose automatic split, which verifies it, does not converge either reports that instead (see ewald).

reciprocal

Reciprocal row vectors, padded by unused identity axes for dimensions below three.

spherical_degree 🔒

Consecutive orders of one spherical degree, sharing their real-space radial integrals while retaining each harmonic’s convergence policy.

spherical_degree_derivatives 🔒

Derivatives of the requested orders of one spherical degree. The jet cache shares the radial derivatives as well as the values, with no change to any order’s shell stopping rule or rounding checks.

sum

The lattice sum of wave at the wavenumber k, the shift shift and the Ewald split parameter eta, as treams.lattice.lsum* returns it.

sum_array

sum_part over arrays, with the residual for the pullback.

sum_part

Evaluate an Ewald component or one direct-summation shell.

volume

Signed determinant for a one-, two- or three-dimensional cell.