08-ilang-propagation / 08-tilt
08tiltverified
Determines whether a push on one cell orients rotation-symmetric exchange cells far away, how small deviations match the anchored oscillator field, and which sources have a long-range mutual cost.
ed_push.py6 KB · raw
"""Check (item 1): exact diagonalization of exchange cells with a push.
C = J sum_<xy> S_x.S_y - h_p S^3_h (J = 1, n = e_3), spin S = 1/2 and S = 1.
Clusters: rings of even size L (S=1/2: 8,12,16,20; S=1: 6,8,10,12) and the
16-object periodic cubic cluster Z^3 / Lambda, Lambda spanned by (2,2,0),(2,0,2),(0,2,2).
C commutes with M = sum_x S^3_x; the lowest eigenvector is searched in M = 0,1,2.
Output: output/ed_push.txt
"""
from pathlib import Path
from collections import deque
import itertools
import numpy as np
import scipy.sparse as sp
from scipy.sparse.linalg import eigsh
OUT = Path(__file__).parent / "output"
OUT.mkdir(exist_ok=True)
HPS = [0.01, 0.1, 0.3, 1.0, 3.0, 10.0]
def ring(L):
links = [(i, (i + 1) % L) for i in range(L)]
eps = np.array([(-1) ** i for i in range(L)])
return links, eps
def cubic16():
A = np.array([[2, 2, 0], [2, 0, 2], [0, 2, 2]]).T # columns are the periods
Ainv = np.linalg.inv(A)
def same(u, v):
c = Ainv @ (np.array(u) - np.array(v))
return np.allclose(c, np.round(c))
reps = []
for p in itertools.product(range(4), repeat=3):
if not any(same(p, q) for q in reps):
reps.append(p)
assert len(reps) == 16
def idx(p):
return [k for k, q in enumerate(reps) if same(p, q)][0]
links = set()
for k, p in enumerate(reps):
nb = set()
for d in np.eye(3, dtype=int):
for s in (1, -1):
j = idx(tuple(np.array(p) + s * d))
nb.add(j)
links.add(tuple(sorted((k, j))))
assert len(nb) == 6
links = sorted(links)
assert len(links) == 48
eps = np.array([(-1) ** (sum(p) % 2) for p in reps])
for a, b in links:
assert eps[a] != eps[b]
return links, eps, reps
def distances(links, n, h=0):
adj = [[] for _ in range(n)]
for a, b in links:
adj[a].append(b)
adj[b].append(a)
d = [-1] * n
d[h] = 0
q = deque([h])
while q:
u = q.popleft()
for v in adj[u]:
if d[v] < 0:
d[v] = d[u] + 1
q.append(v)
return np.array(d)
def sector(n, S2, M2):
"""Basis states (digit arrays, digit = m + S) with 2*sum(m) = M2; S2 = 2S."""
d = S2 + 1
codes = np.arange(d ** n, dtype=np.int64)
digits = np.zeros((codes.size, n), dtype=np.int8)
c = codes.copy()
for i in range(n):
digits[:, i] = c % d
c //= d
tot = 2 * digits.sum(axis=1).astype(np.int64) - n * S2
keep = tot == M2
return codes[keep], digits[keep]
def build(n, S2, M2, links, J=1.0):
S = S2 / 2.0
d = S2 + 1
codes, digits = sector(n, S2, M2)
m = digits.astype(float) - S
diag = np.zeros(codes.size)
rows, cols, vals = [], [], []
for a, b in links:
diag += J * m[:, a] * m[:, b]
for i, j in ((a, b), (b, a)): # S^+_i S^-_j
ok = (digits[:, i] < S2) & (digits[:, j] > 0)
src = np.nonzero(ok)[0]
new = codes[src] + d ** i - d ** j
dst = np.searchsorted(codes, new)
amp = 0.5 * J * np.sqrt(S * (S + 1) - m[src, i] * (m[src, i] + 1)) \
* np.sqrt(S * (S + 1) - m[src, j] * (m[src, j] - 1))
rows.append(dst); cols.append(src); vals.append(amp)
off = sp.csr_matrix((np.concatenate(vals), (np.concatenate(rows), np.concatenate(cols))),
shape=(codes.size, codes.size))
return off, diag, m
def lowest(H, k=1):
if H.shape[0] <= 1500:
w, v = np.linalg.eigh(H.toarray())
return w[:k], v[:, :k]
w, v = eigsh(H, k=k, which="SA", tol=1e-12)
o = np.argsort(w)
return w[o], v[:, o]
def run(name, n, S2, links, eps, out):
dist = distances(links, n)
far = np.nonzero(dist == dist.max())[0]
out.append(f"\n## {name}: N={n}, S={S2}/2, links={len(links)}, far sites {list(far)} at distance {dist.max()}")
secs = {}
for M2 in (0, 2, 4):
if M2 // 2 > n * S2 // 2:
continue
secs[M2] = build(n, S2, M2, links)
# item 1(a): no push
w0, _ = lowest(secs[0][0] + sp.diags(secs[0][1]), k=2)
w1, _ = lowest(secs[2][0] + sp.diags(secs[2][1]), k=1)
out.append(f"no push: E0(M=0)={w0[0]:.10f}, E1(M=0)={w0[1]:.10f}, E0(M=1)={w1[0]:.10f}"
f" -> unique singlet: {w0[0] < w0[1] - 1e-9 and w0[0] < w1[0] - 1e-9}")
out.append("h_p/J M_low E_low m(far)=eps*<S3_far> m(h)")
for hp in HPS:
best = None
for M2, (off, diag, m) in secs.items():
H = off + sp.diags(diag - hp * m[:, 0])
w, v = lowest(H, k=1)
if best is None or w[0] < best[0] - 1e-12:
best = (w[0], v[:, 0], m, M2)
E, psi, m, M2 = best
p = np.abs(psi) ** 2
sz = p @ m
mfar = np.mean(eps[far] * sz[far])
out.append(f"{hp:6.2f} {M2/2:5.1f} {E:.10f} {mfar: .6f} {sz[0]: .6f}")
def main():
out = ["# ED check, item 1: exchange cells, J=1, push -h_p S^3_h at h=0 (class +)"]
for L in (8, 12, 16, 20):
links, eps = ring(L)
run(f"ring L={L}", L, 1, links, eps, out)
links, eps, reps = cubic16()
out.append("\n16-object cubic cluster sites (representatives mod Lambda): " + str(reps))
out.append("links: each site x to x+e_i, i=1,2,3 (mod Lambda); 48 links, 6 distinct neighbours each")
run("cubic cluster 16", 16, 1, links, eps, out)
for L in (6, 8, 10, 12):
links, eps = ring(L)
run(f"ring L={L}", L, 2, links, eps, out)
(OUT / "ed_push.txt").write_text("\n".join(out) + "\n")
print("\n".join(out))
if __name__ == "__main__":
main()
spinwave_mutual.py10 KB · raw
"""Check (item 3, leading order in 1/S): mutual cost of two held sources on the cubic torus.
Units J S = 1. Quadratic form C_2 = 1/2 P^T T P + 1/2 Q^T V Q, T_k = z(1+gamma_k), V_k = z(1-gamma_k).
(1) Transverse direction sources (linear in Q): E_12 = -S eps1 eps2 g1 g2 G_a(x1,x2); the separation part is
J E_r/(g1 g2) = -eps1 eps2 ghat(r), ghat = Lap^+ (single objects); -eps1 eps2 sum_{a,b} ghat(a-b) (two-object cells).
(2) Quadratic sources dH_i = 1/2 Q^T A_i Q + 1/2 P^T B_i P (along-n sources, bound held spin as a spin-length defect,
changed link): E_12 = -2 int_0^inf dtau <dH_1(tau) dH_2(0)>_c, Wick with G_QQ, G_PP and G_QP = (i/2) F (tau > 0);
far-from-push limit: the two zero modes (k = 0, k = Q) are left out.
Output: output/spinwave_mutual.txt
"""
from pathlib import Path
import numpy as np
from scipy.integrate import quad
from scipy.special import ive
OUT = Path(__file__).parent / "output"
OUT.mkdir(exist_ok=True)
Z = 6
def kgrid(L):
k = 2 * np.pi * np.arange(L) / L
c = np.cos(k)
g = (c[:, None, None] + c[None, :, None] + c[None, None, :]) / 3.0
return g
def ghat(L):
k = 2 * np.pi * np.arange(L) / L
lap = 2 * ((1 - np.cos(k))[:, None, None] + (1 - np.cos(k))[None, :, None] + (1 - np.cos(k))[None, None, :])
inv = np.zeros_like(lap)
inv[lap > 1e-12] = 1.0 / lap[lap > 1e-12]
return np.real(np.fft.ifftn(inv))
def ghat_inf(x):
"""Infinite cubic lattice: (Lap^{-1})(x) = int_0^inf prod_i e^{-2t} I_{x_i}(2t) dt."""
f = lambda t: np.prod([ive(abs(int(c)), 2 * t) for c in x])
return quad(f, 0, np.inf, limit=400, epsabs=1e-13, epsrel=1e-11)[0]
def eps(x):
return (-1) ** (int(sum(x)) % 2)
def check_58(L=6, eta=0.7):
n = L ** 3
co = np.array(np.unravel_index(np.arange(n), (L,) * 3)).T
adj = np.zeros((n, n))
for d in range(3):
sh = co.copy(); sh[:, d] = (sh[:, d] + 1) % L
j = np.ravel_multi_index(sh.T, (L,) * 3)
adj[np.arange(n), j] += 1; adj[j, np.arange(n)] += 1
V = Z * np.eye(n) - adj
V[0, 0] += eta
Ga = np.linalg.inv(V)
g = ghat(L).ravel()
dif = lambda a, b: np.ravel_multi_index(tuple((co[a] - co[b]) % L), (L,) * 3)
err = 0.0
for a, b in [(5, 17), (40, 101), (3, 3), (0, 77)]:
f = 1 / eta + g[dif(a, b)] - g[dif(a, 0)] - g[dif(b, 0)] + g[0]
err = max(err, abs(f - Ga[a, b]))
return err
def pert(kind, x, e=None):
"""Return (A, B) as lists of (site, site, value); per unit strength (see header of output)."""
x = tuple(x)
if kind == "along":
A = [(x, x, -eps(x))]
return A, A
if kind == "along_cell":
x2 = (x[0], x[1] + 1, x[2])
A = [(x, x, -eps(x)), (x2, x2, -eps(x))]
return A, A
if kind == "defect":
nb = []
for d in range(3):
for s in (1, -1):
y = list(x); y[d] += s; nb.append(tuple(y))
A = [(x, x, -Z)] + [(x, y, 0.5) for y in nb] + [(y, x, 0.5) for y in nb]
B = [(x, x, -Z)] + [(x, y, -0.5) for y in nb] + [(y, x, -0.5) for y in nb]
return A, B
if kind == "link":
y = tuple(np.array(x) - np.array(e))
A = [(x, x, 1), (y, y, 1), (x, y, -1), (y, x, -1)]
B = [(x, x, 1), (y, y, 1), (x, y, 1), (y, x, 1)]
return A, B
raise ValueError(kind)
def contraction(M1, M2, L):
w, i1, i2 = [], [], []
for a, b, u in M1:
for c, d, v in M2:
w.append(u * v)
i1.append(np.ravel_multi_index(tuple(np.mod(np.array(a) - np.array(c), L)), (L,) * 3))
i2.append(np.ravel_multi_index(tuple(np.mod(np.array(b) - np.array(d), L)), (L,) * 3))
return np.array(w), np.array(i1), np.array(i2)
def mutual_quadratic(L, configs, nodes=400):
g = kgrid(L)
om = Z * np.sqrt(np.clip(1 - g ** 2, 0, None))
mask = om > 1e-9
Tk = np.where(mask, Z * (1 + g), 0.0)
Vk = np.where(mask, Z * (1 - g), 0.0)
omm = np.where(mask, om, 1.0)
pre = []
for (A1, B1), (A2, B2) in configs:
pre.append((contraction(A1, A2, L), contraction(B1, B2, L),
contraction(A1, B2, L), contraction(B1, A2, L)))
s, ws = np.polynomial.legendre.leggauss(nodes)
lo, hi = np.log(1e-7), np.log(200.0)
sv = 0.5 * (hi - lo) * s + 0.5 * (hi + lo)
wv = 0.5 * (hi - lo) * ws
acc = np.zeros(len(configs))
for sj, wj in zip(sv, wv):
tau = np.exp(sj)
ex = np.where(mask, np.exp(-omm * tau), 0.0)
GQ = np.real(np.fft.ifftn(Tk / (2 * omm) * ex)).ravel()
GP = np.real(np.fft.ifftn(Vk / (2 * omm) * ex)).ravel()
F = np.real(np.fft.ifftn(ex)).ravel()
for n, (qq, pp, qp, pq) in enumerate(pre):
val = 0.5 * np.sum(qq[0] * GQ[qq[1]] * GQ[qq[2]]) + 0.5 * np.sum(pp[0] * GP[pp[1]] * GP[pp[2]]) \
- 0.125 * (np.sum(qp[0] * F[qp[1]] * F[qp[2]]) + np.sum(pq[0] * F[pq[1]] * F[pq[2]]))
acc[n] += wj * tau * val
return -2 * acc
def cont_link(e1, e2, e12):
"""Continuum: J E/(delta1 delta2 S) * r^7 for two changed links (unit vectors e1, e2), r along axis 1."""
def C(t):
X = np.array([t, 1.0, 0, 0]); rho = np.linalg.norm(X); n = X / rho
I = np.eye(4) - 4 * np.outer(n, n)
E1 = np.r_[0, e1]; E2 = np.r_[0, e2]; t0 = np.r_[1.0, 0, 0, 0]
M1 = np.outer(E1, E1) - (2 / Z) * np.outer(t0, t0)
M2 = np.outer(E2, E2) - (2 / Z) * np.outer(t0, t0)
cross = I[0, 0] * (E1 @ I @ E2) - (I @ E1)[0] * (I @ E2)[0]
return (4 * np.trace(M1 @ I @ M2 @ I) - (16 / Z) * e12 * cross) / rho ** 8
val = quad(C, -np.inf, np.inf, epsabs=1e-13, epsrel=1e-12)[0]
return -np.sqrt(2 * Z) / (16 * np.pi ** 4) * val
def main():
out = ["# Mutual cost of two held sources, leading order in 1/S, cubic torus (units as stated per block)"]
out.append(f"\n(5.8) with kappa_h = J S, kappa_a = h_p against direct inverse (L=6): max error {check_58():.2e}")
out.append("\n## (1) Transverse sources: separation part J E_r/(g1 g2 (m1perp.m2perp))")
out.append("single objects: -eps1 eps2 ghat(r); cells (x, x+e2), + on the first object: -eps1 eps2 sum ghat")
out.append(f"ghat_inf(0) = {ghat_inf((0, 0, 0)):.7f} (gamma_3 = 0.2527310)")
out.append("r eps1eps2 single: J E_r/(g1g2) 4 pi r ghat_inf(r) cells: J E_r/(g1g2) pi r sum ghat_inf"
" torus check: ghat_L(r)-ghat_L(8) for L=16, 32, 64 vs ghat_inf(r)-ghat_inf(8)")
gl = {L: ghat(L) for L in (16, 32, 64)}
g8 = ghat_inf((8, 0, 0))
for r in range(2, 9):
e12 = (-1) ** r
gs = ghat_inf((r, 0, 0))
gc = 2 * gs + ghat_inf((r, 1, 0)) + ghat_inf((r, -1, 0))
tor = " ".join(f"{gl[L][r, 0, 0] - gl[L][8, 0, 0]:.6f}" for L in (16, 32, 64))
out.append(f"{r:2d} {e12:+d} {-e12 * gs: .6f} {4 * np.pi * r * gs:.4f}"
f" {-e12 * gc: .6f} {np.pi * r * gc:.4f} {tor} {gs - g8:.6f}")
kap = np.sqrt(2 * Z) / (32 * np.pi ** 3)
lam = Z ** 1.5 / (16 * np.sqrt(2) * np.pi ** 3)
cA = {1: cont_link([1, 0, 0], [1, 0, 0], 1), -1: cont_link([1, 0, 0], [1, 0, 0], -1)}
cB = {1: cont_link([0, 1, 0], [0, 1, 0], 1), -1: cont_link([0, 1, 0], [0, 1, 0], -1)}
pre = np.sqrt(3) / (8 * np.pi ** 3)
out.append("\n## Continuum coefficients")
out.append(f"along n, single: J S E/(g1 g2) = -eps1 eps2 kappa_par / r^3, kappa_par = sqrt(2z)/(32 pi^3) = {kap:.7f}")
out.append(f"held spin (defect), r >> l_b: S E/(J sig1 sig2) = +eps1 eps2 lambda / r^5, lambda = {lam:.7f}")
out.append(f"changed links, J E/(d1 d2 S) * r^7: parallel to r: eps=+1 {cA[1]:.7f} (formula {-pre * (43/6 + 2.5):.7f}),"
f" eps=-1 {cA[-1]:.7f} (formula {-pre * (43/6 - 2.5):.7f})")
out.append(f" perpendicular: eps=+1 {cB[1]:.7f} (formula {-pre * (4/3 - 5/12):.7f}),"
f" eps=-1 {cB[-1]:.7f} (formula {-pre * (4/3 + 5/12):.7f})")
out.append("\n## (2) Quadratic sources, lattice (zero modes left out); ratio = lattice / continuum")
out.append("columns: r, eps1eps2, along single [J S E/(g1g2), ratio], along cell [.., ratio vs 4x], "
"defect [S E/(J s1 s2), ratio], link par [J E/(d1d2S), ratio], link perp [.., ratio]")
rs = list(range(2, 13))
ratios = {}
for L in (32, 64, 128):
configs = []
for r in rs:
x2 = (r, 0, 0)
configs += [(pert("along", (0, 0, 0)), pert("along", x2)),
(pert("along_cell", (0, 0, 0)), pert("along_cell", x2)),
(pert("defect", (0, 0, 0)), pert("defect", x2)),
(pert("link", (0, 0, 0), (1, 0, 0)), pert("link", x2, (1, 0, 0))),
(pert("link", (0, 0, 0), (0, 1, 0)), pert("link", x2, (0, 1, 0)))]
res = mutual_quadratic(L, configs, nodes=300).reshape(len(rs), 5)
out.append(f"L = {L}")
ratios[L] = np.zeros_like(res)
for i, r in enumerate(rs):
e = (-1) ** r
p = [-e * kap / r ** 3, -4 * e * kap / r ** 3, e * lam / r ** 5, cA[e] / r ** 7, cB[e] / r ** 7]
ratios[L][i] = res[i] / np.array(p)
row = " ".join(f"{res[i, j]: .3e} {res[i, j] / p[j]:6.3f}" for j in range(5))
out.append(f"{r:2d} {e:+d} {row}")
out.append("ratio extrapolated to L -> infinity: Aitken from L = 32, 64, 128 (geometric finite-size error);"
" the L = 128 value where the change from 64 to 128 is below 1e-3 or not geometric")
for i, r in enumerate(rs):
x1, x2, x3 = ratios[32][i], ratios[64][i], ratios[128][i]
d1, d2 = x2 - x1, x3 - x2
geo = (np.abs(d2) > 1e-3) & (d1 * d2 > 0) & (np.abs(d2) < np.abs(d1))
ex = np.where(geo, x3 + d2 / np.where(geo, d1 / np.where(geo, d2, 1.0) - 1, 1.0), x3)
out.append(f"{r:2d} " + " ".join(f"{v:6.3f}" for v in ex))
(OUT / "spinwave_mutual.txt").write_text("\n".join(out) + "\n")
print("\n".join(out))
if __name__ == "__main__":
main()
spinwave_tilt.py3 KB · raw
"""Check (item 1(b), leading order in 1/S): first correction to m(x) far from the push.
Quadratic form (units J S = 1): C_2 = 1/2 P^T T P + 1/2 Q^T V Q with
V = Lap + eta P_h, T = E V E = Lap^+ + eta P_h, eta = h_p/(J S), E = diag(eps).
Ground-state covariances: <QQ^T> = 1/2 T^(1/2) (T^(1/2) V T^(1/2))^(-1/2) T^(1/2), same for P with T<->V.
First correction: m(x) = eps_h (S - dS(x)), dS(x) = (<Q_x^2> + <P_x^2> - 1)/2.
Output: output/spinwave_tilt.txt
"""
from pathlib import Path
import numpy as np
OUT = Path(__file__).parent / "output"
OUT.mkdir(exist_ok=True)
def lattice(L, D):
n = L ** D
coords = np.array(np.unravel_index(np.arange(n), (L,) * D)).T
adj = np.zeros((n, n))
for d in range(D):
sh = coords.copy()
sh[:, d] = (sh[:, d] + 1) % L
j = np.ravel_multi_index(sh.T, (L,) * D)
adj[np.arange(n), j] += 1
adj[j, np.arange(n)] += 1
eps = (-1.0) ** coords.sum(axis=1)
return coords, adj, eps
def sqrtm_sym(A, power):
w, U = np.linalg.eigh(A)
return (U * w ** power) @ U.T
def dS_push(L, D, eta):
coords, adj, eps = lattice(L, D)
z = 2 * D
n = L ** D
V = z * np.eye(n) - adj
V[0, 0] += eta
T = (eps[:, None] * V) * eps[None, :]
Vh = sqrtm_sym(V, 0.5)
Th = (eps[:, None] * Vh) * eps[None, :]
covQ = 0.5 * Th @ sqrtm_sym(Th @ V @ Th, -0.5) @ Th
covP = 0.5 * Vh @ sqrtm_sym(Vh @ T @ Vh, -0.5) @ Vh
far = np.ravel_multi_index(tuple(np.full(D, L // 2)), (L,) * D)
dS = 0.5 * (np.diag(covQ) + np.diag(covP) - 1)
return dS[far], dS[0], np.max(np.abs(np.diag(covQ) - np.diag(covP)))
def dS3_bulk(L):
k = 2 * np.pi * (np.arange(L) + 0.5) / L
c = np.cos(k)
tot = 0.0
for c1 in c:
g = (c1 + c[:, None] + c[None, :]) / 3.0
tot += np.sum(1.0 / np.sqrt(1.0 - g ** 2))
W = tot / L ** 3
return 0.5 * (W - 1.0)
def main():
out = ["# Leading order in 1/S: first correction dS, m = eps_h (S - dS) + O(1/S)"]
out.append("\n## Cubic torus, L -> infinity, no push (midpoint sums): dS_3 = (W - 1)/2")
vals = {}
for L in (64, 128, 256):
vals[L] = dS3_bulk(L)
out.append(f"L={L:4d}: dS_3 = {vals[L]:.6f}")
rich = (4 * vals[256] - vals[128]) / 3
out.append(f"Richardson (error ~ L^-2): dS_3 = {rich:.5f}")
out.append("\n## Cubic torus with push (exact Gaussian), far site (L/2,L/2,L/2)")
out.append("L eta=h_p/(JS) dS(far) dS(h) zero-mode estimate sqrt(2z/(N eta))/2 max|<Q^2>-<P^2>|")
for L in (8, 12, 16):
for eta in (0.1, 1.0, 10.0):
f, h, sym = dS_push(L, 3, eta)
est = 0.5 * np.sqrt(12.0 / (L ** 3 * eta))
out.append(f"{L:3d} {eta:6.2f} {f:.5f} {h:.5f} {est:.5f} {sym:.1e}")
out.append("\n## Ring with push (exact Gaussian), eta = 1, far site L/2")
out.append("L dS(L/2) dS(L/2) - (1/pi) ln(L/2)")
ring = []
for L in (64, 128, 256, 512, 1024):
f, h, _ = dS_push(L, 1, 1.0)
ring.append((L, f))
out.append(f"{L:5d} {f:.5f} {f - np.log(L / 2) / np.pi:.5f}")
sl = np.polyfit(np.log([r[0] for r in ring[-3:]]), [r[1] for r in ring[-3:]], 1)[0]
out.append(f"slope d dS / d ln L (last three sizes) = {sl:.5f}; 1/pi = {1 / np.pi:.5f} (Dirichlet image at the push)")
(OUT / "spinwave_tilt.txt").write_text("\n".join(out) + "\n")
print("\n".join(out))
if __name__ == "__main__":
main()
output/ed_push.txt6 KB · raw
# ED check, item 1: exchange cells, J=1, push -h_p S^3_h at h=0 (class +)
## ring L=8: N=8, S=1/2, links=8, far sites [np.int64(4)] at distance 4
no push: E0(M=0)=-3.6510934089, E1(M=0)=-3.1284190638, E0(M=1)=-3.1284190638 -> unique singlet: True
h_p/J M_low E_low m(far)=eps*<S3_far> m(h)
0.01 0.0 -3.6511215930 0.002915 0.005637
0.10 0.0 -3.6539014094 0.028913 0.055951
0.30 0.0 -3.6756584654 0.081547 0.158772
1.00 0.0 -3.8706747638 0.178327 0.361393
3.00 0.0 -4.7406497893 0.217820 0.469512
10.00 0.0 -8.1656749710 0.223623 0.496231
## ring L=12: N=12, S=1/2, links=12, far sites [np.int64(6)] at distance 6
no push: E0(M=0)=-5.3873909174, E1(M=0)=-5.0315434037, E0(M=1)=-5.0315434037 -> unique singlet: True
h_p/J M_low E_low m(far)=eps*<S3_far> m(h)
0.01 0.0 -5.3874232606 0.003078 0.006468
0.10 0.0 -5.3906079037 0.030399 0.063994
0.30 0.0 -5.4152047132 0.083233 0.177593
1.00 0.0 -5.6229268437 0.162829 0.374293
3.00 0.0 -6.5029771989 0.187370 0.471246
10.00 0.0 -9.9314563940 0.190367 0.496371
## ring L=16: N=16, S=1/2, links=16, far sites [np.int64(8)] at distance 8
no push: E0(M=0)=-7.1422963606, E1(M=0)=-6.8721066784, E0(M=1)=-6.8721066784 -> unique singlet: True
h_p/J M_low E_low m(far)=eps*<S3_far> m(h)
0.01 0.0 -7.1423316575 0.003176 0.007059
0.10 0.0 -7.1458016524 0.031226 0.069621
0.30 0.0 -7.1722882312 0.083171 0.189514
1.00 0.0 -7.3866985161 0.150071 0.379989
3.00 0.0 -8.2707749474 0.167171 0.471891
10.00 0.0 -11.7005556368 0.169000 0.496425
## ring L=20: N=20, S=1/2, links=20, far sites [np.int64(10)] at distance 10
no push: E0(M=0)=-8.9043865299, E1(M=0)=-8.6864409862, E0(M=1)=-8.6864409862 -> unique singlet: True
h_p/J M_low E_low m(far)=eps*<S3_far> m(h)
0.01 0.0 -8.9044241355 0.003245 0.007520
0.10 0.0 -8.9081154951 0.031764 0.073953
0.30 0.0 -8.9359905145 0.082444 0.197877
1.00 0.0 -9.1544799698 0.139769 0.383032
3.00 0.0 -10.0405891867 0.152531 0.472206
10.00 0.0 -13.4710121778 0.153753 0.496452
16-object cubic cluster sites (representatives mod Lambda): [(0, 0, 0), (0, 0, 1), (0, 0, 2), (0, 0, 3), (0, 1, 0), (0, 1, 1), (0, 1, 2), (0, 1, 3), (1, 0, 0), (1, 0, 1), (1, 0, 2), (1, 0, 3), (1, 1, 0), (1, 1, 1), (1, 1, 2), (1, 1, 3)]
links: each site x to x+e_i, i=1,2,3 (mod Lambda); 48 links, 6 distinct neighbours each
## cubic cluster 16: N=16, S=1/2, links=48, far sites [np.int64(13), np.int64(15)] at distance 3
no push: E0(M=0)=-15.5783288133, E1(M=0)=-14.7965997285, E0(M=1)=-14.7965997285 -> unique singlet: True
h_p/J M_low E_low m(far)=eps*<S3_far> m(h)
0.01 0.0 -15.5783455887 0.002435 0.003355
0.10 0.0 -15.5800039442 0.024280 0.033454
0.30 0.0 -15.5932340835 0.071115 0.098118
1.00 0.0 -15.7272156212 0.192088 0.268207
3.00 0.0 -16.4641224768 0.294607 0.424892
10.00 0.0 -19.7476521880 0.330312 0.487059
## ring L=6: N=6, S=2/2, links=6, far sites [np.int64(3)] at distance 3
no push: E0(M=0)=-8.6174231818, E1(M=0)=-7.8967958192, E0(M=1)=-7.8967958192 -> unique singlet: True
h_p/J M_low E_low m(far)=eps*<S3_far> m(h)
0.01 0.0 -8.6174901294 0.010352 0.013389
0.10 0.0 -8.6240752905 0.102132 0.132193
0.30 0.0 -8.6745349327 0.278218 0.362050
1.00 0.0 -9.0869706859 0.546592 0.732523
3.00 0.0 -10.7994879963 0.660086 0.919464
10.00 0.0 -17.5670872881 0.694926 0.986112
## ring L=8: N=8, S=2/2, links=8, far sites [np.int64(4)] at distance 4
no push: E0(M=0)=-11.3369560779, E1(M=0)=-10.7434008235, E0(M=1)=-10.7434008235 -> unique singlet: True
h_p/J M_low E_low m(far)=eps*<S3_far> m(h)
0.01 0.0 -11.3370291952 0.010221 0.014622
0.10 0.0 -11.3442097178 0.100476 0.143920
0.30 0.0 -11.3985873076 0.267565 0.386678
1.00 0.0 -11.8253272598 0.494463 0.745935
3.00 0.0 -13.5497976377 0.582260 0.922190
10.00 0.0 -20.3241940732 0.609254 0.986461
## ring L=10: N=10, S=2/2, links=10, far sites [np.int64(5)] at distance 5
no push: E0(M=0)=-14.0941299549, E1(M=0)=-13.5693220045, E0(M=1)=-13.5693220045 -> unique singlet: True
h_p/J M_low E_low m(far)=eps*<S3_far> m(h)
0.01 0.0 -14.0942062020 0.009913 0.015248
0.10 0.0 -14.1016869890 0.097223 0.149801
0.30 0.0 -14.1579598525 0.255412 0.398181
1.00 0.0 -14.5907289613 0.456752 0.751132
3.00 0.0 -16.3197109015 0.529673 0.923214
10.00 0.0 -23.0967344366 0.550821 0.986603
## ring L=12: N=12, S=2/2, links=12, far sites [np.int64(6)] at distance 6
no push: E0(M=0)=-16.8695561395, E1(M=0)=-16.3853596696, E0(M=1)=-16.3853596696 -> unique singlet: True
h_p/J M_low E_low m(far)=eps*<S3_far> m(h)
0.01 0.0 -16.8696337457 0.009055 0.015520
0.10 0.0 -16.8772443685 0.088705 0.152333
0.30 0.0 -16.9343124777 0.231425 0.402886
1.00 0.0 -17.3693927041 0.407686 0.753039
3.00 0.0 -19.1000782201 0.470341 0.923626
10.00 0.0 -25.8782395163 0.488416 0.986670
output/spinwave_mutual.txt6 KB · raw
# Mutual cost of two held sources, leading order in 1/S, cubic torus (units as stated per block)
(5.8) with kappa_h = J S, kappa_a = h_p against direct inverse (L=6): max error 4.88e-14
## (1) Transverse sources: separation part J E_r/(g1 g2 (m1perp.m2perp))
single objects: -eps1 eps2 ghat(r); cells (x, x+e2), + on the first object: -eps1 eps2 sum ghat
ghat_inf(0) = 0.2527310 (gamma_3 = 0.2527310)
r eps1eps2 single: J E_r/(g1g2) 4 pi r ghat_inf(r) cells: J E_r/(g1g2) pi r sum ghat_inf torus check: ghat_L(r)-ghat_L(8) for L=16, 32, 64 vs ghat_inf(r)-ghat_inf(8)
2 +1 -0.042889 1.0779 -0.157642 0.9905 0.029359 0.032566 0.032862 0.032901
3 -1 0.027545 1.0384 0.106137 1.0003 0.014233 0.017247 0.017521 0.017557
4 +1 -0.020289 1.0198 -0.079621 1.0005 0.007303 0.010028 0.010269 0.010301
5 -1 0.016101 1.0117 0.063679 1.0003 0.003573 0.005889 0.006087 0.006113
6 +1 -0.013365 1.0077 -0.053058 1.0001 0.001456 0.003213 0.003358 0.003377
7 -1 0.011430 1.0055 0.045475 1.0001 0.000347 0.001353 0.001432 0.001442
8 +1 -0.009988 1.0041 -0.039790 1.0000 0.000000 0.000000 0.000000 0.000000
## Continuum coefficients
along n, single: J S E/(g1 g2) = -eps1 eps2 kappa_par / r^3, kappa_par = sqrt(2z)/(32 pi^3) = 0.0034913
held spin (defect), r >> l_b: S E/(J sig1 sig2) = +eps1 eps2 lambda / r^5, lambda = 0.0209480
changed links, J E/(d1 d2 S) * r^7: parallel to r: eps=+1 -0.0674991 (formula -0.0674991), eps=-1 -0.0325858 (formula -0.0325858)
perpendicular: eps=+1 -0.0064008 (formula -0.0064008), eps=-1 -0.0122197 (formula -0.0122197)
## (2) Quadratic sources, lattice (zero modes left out); ratio = lattice / continuum
columns: r, eps1eps2, along single [J S E/(g1g2), ratio], along cell [.., ratio vs 4x], defect [S E/(J s1 s2), ratio], link par [J E/(d1d2S), ratio], link perp [.., ratio]
L = 32
2 +1 -4.624e-04 1.060 -1.491e-03 0.854 5.671e-04 0.866 -1.274e-03 2.417 -1.174e-04 2.347
3 -1 1.198e-04 0.927 4.142e-04 0.801 -1.211e-04 1.405 -2.430e-05 1.631 -1.119e-05 2.002
4 +1 -3.922e-05 0.719 -1.480e-04 0.678 1.949e-05 0.952 -6.880e-06 1.670 -6.874e-07 1.759
5 -1 1.640e-05 0.587 6.200e-05 0.555 -6.361e-06 0.949 -5.446e-07 1.306 -2.119e-07 1.355
6 +1 -7.170e-06 0.444 -2.790e-05 0.432 2.055e-06 0.763 -3.006e-07 1.247 -2.761e-08 1.207
7 -1 3.351e-06 0.329 1.293e-05 0.318 -8.651e-07 0.694 -4.524e-08 1.143 -1.714e-08 1.155
8 +1 -1.510e-06 0.221 -5.918e-06 0.217 3.521e-07 0.551 -3.768e-08 1.171 -3.017e-09 0.988
9 -1 6.659e-07 0.139 2.578e-06 0.135 -1.646e-07 0.464 -7.892e-09 1.158 -2.875e-09 1.125
10 +1 -2.585e-07 0.074 -1.014e-06 0.073 7.468e-08 0.356 -8.418e-09 1.247 -5.211e-10 0.814
11 -1 8.809e-08 0.034 3.369e-07 0.032 -3.849e-08 0.296 -2.283e-09 1.366 -7.612e-10 1.214
12 +1 -2.467e-08 0.012 -9.778e-08 0.012 2.099e-08 0.249 -2.906e-09 1.543 -1.184e-10 0.663
L = 64
2 +1 -4.970e-04 1.139 -1.616e-03 0.926 5.756e-04 0.879 -1.274e-03 2.416 -1.174e-04 2.348
3 -1 1.393e-04 1.077 4.876e-04 0.943 -1.251e-04 1.451 -2.430e-05 1.631 -1.119e-05 2.002
4 +1 -5.146e-05 0.943 -1.957e-04 0.897 2.119e-05 1.036 -6.868e-06 1.667 -6.922e-07 1.772
5 -1 2.483e-05 0.889 9.493e-05 0.850 -7.344e-06 1.096 -5.434e-07 1.303 -2.117e-07 1.353
6 +1 -1.312e-05 0.812 -5.144e-05 0.796 2.620e-06 0.972 -2.969e-07 1.231 -2.873e-08 1.256
7 -1 7.695e-06 0.756 3.007e-05 0.738 -1.231e-06 0.987 -4.446e-08 1.124 -1.699e-08 1.145
8 +1 -4.679e-06 0.686 -1.850e-05 0.678 5.840e-07 0.914 -3.586e-08 1.114 -3.397e-09 1.113
9 -1 2.998e-06 0.626 1.181e-05 0.616 -3.205e-07 0.904 -7.274e-09 1.068 -2.759e-09 1.080
10 +1 -1.950e-06 0.559 -7.737e-06 0.554 1.761e-07 0.841 -7.248e-09 1.074 -6.713e-10 1.049
11 -1 1.306e-06 0.498 5.163e-06 0.492 -1.058e-07 0.814 -1.750e-09 1.047 -6.603e-10 1.053
12 +1 -8.771e-07 0.434 -3.487e-06 0.431 6.334e-08 0.752 -1.997e-09 1.060 -1.788e-10 1.001
L = 128
2 +1 -5.071e-04 1.162 -1.653e-03 0.947 5.768e-04 0.881 -1.274e-03 2.416 -1.174e-04 2.348
3 -1 1.453e-04 1.124 5.105e-04 0.987 -1.256e-04 1.457 -2.430e-05 1.631 -1.119e-05 2.002
4 +1 -5.554e-05 1.018 -2.116e-04 0.970 2.146e-05 1.049 -6.867e-06 1.667 -6.925e-07 1.773
5 -1 2.784e-05 0.997 1.068e-04 0.956 -7.510e-06 1.120 -5.434e-07 1.303 -2.117e-07 1.353
6 +1 -1.545e-05 0.956 -6.063e-05 0.938 2.724e-06 1.011 -2.967e-07 1.230 -2.881e-08 1.260
7 -1 9.546e-06 0.938 3.740e-05 0.919 -1.304e-06 1.046 -4.445e-08 1.123 -1.699e-08 1.145
8 +1 -6.180e-06 0.906 -2.447e-05 0.897 6.354e-07 0.994 -3.578e-08 1.112 -3.429e-09 1.123
9 -1 4.239e-06 0.885 1.674e-05 0.874 -3.590e-07 1.012 -7.265e-09 1.066 -2.757e-09 1.079
10 +1 -2.985e-06 0.855 -1.186e-05 0.849 2.049e-07 0.978 -7.205e-09 1.067 -6.865e-10 1.072
11 -1 2.178e-06 0.830 8.638e-06 0.823 -1.282e-07 0.986 -1.743e-09 1.042 -6.589e-10 1.051
12 +1 -1.616e-06 0.800 -6.433e-06 0.796 8.072e-08 0.959 -1.970e-09 1.046 -1.869e-10 1.046
ratio extrapolated to L -> infinity: Aitken from L = 32, 64, 128 (geometric finite-size error); the L = 128 value where the change from 64 to 128 is below 1e-3 or not geometric
2 1.172 0.956 0.881 2.416 2.348
3 1.144 1.007 1.459 1.631 2.002
4 1.055 1.007 1.051 1.667 1.773
5 1.057 1.016 1.125 1.303 1.353
6 1.047 1.029 1.020 1.230 1.260
7 1.073 1.053 1.060 1.123 1.145
8 1.104 1.095 1.017 1.112 1.124
9 1.180 1.170 1.047 1.066 1.079
10 1.322 1.318 1.032 1.067 1.075
11 1.672 1.675 1.071 1.042 1.051
12 3.173 3.225 1.103 1.046 1.053
output/spinwave_tilt.txt2 KB · raw
# Leading order in 1/S: first correction dS, m = eps_h (S - dS) + O(1/S)
## Cubic torus, L -> infinity, no push (midpoint sums): dS_3 = (W - 1)/2
L= 64: dS_3 = 0.078304
L= 128: dS_3 = 0.078344
L= 256: dS_3 = 0.078354
Richardson (error ~ L^-2): dS_3 = 0.07836
## Cubic torus with push (exact Gaussian), far site (L/2,L/2,L/2)
L eta=h_p/(JS) dS(far) dS(h) zero-mode estimate sqrt(2z/(N eta))/2 max|<Q^2>-<P^2>|
8 0.10 0.31135 0.29896 0.24206 0.0e+00
8 1.00 0.15197 0.11065 0.07655 0.0e+00
8 10.00 0.11169 0.01999 0.02421 0.0e+00
12 0.10 0.20637 0.19853 0.13176 0.0e+00
12 1.00 0.11960 0.08897 0.04167 0.0e+00
12 10.00 0.09770 0.01863 0.01318 0.0e+00
16 0.10 0.16197 0.15613 0.08558 0.0e+00
16 1.00 0.10561 0.07974 0.02706 0.0e+00
16 10.00 0.09139 0.01805 0.00856 0.0e+00
## Ring with push (exact Gaussian), eta = 1, far site L/2
L dS(L/2) dS(L/2) - (1/pi) ln(L/2)
64 1.10378 0.00061
128 1.31488 -0.00893
256 1.53065 -0.01380
512 1.74882 -0.01626
1024 1.96822 -0.01750
slope d dS / d ln L (last three sizes) = 0.31564; 1/pi = 0.31831 (Dirichlet image at the push)
nDot Physics Research