1"""Finite-grid polarization-gap and spherical-triangle Chern monitor."""
 2import numpy as np
 3
 4
 5def spherical_area(a, b, c):
 6    num = np.einsum('...i,...i->...', a, np.cross(b, c))
 7    den = 1 + np.einsum('...i,...i->...', a, b) + np.einsum('...i,...i->...', b, c) + np.einsum('...i,...i->...', c, a)
 8    return 2 * np.arctan2(num, den)
 9
10
11def monitor_field(m, return_triangles=False):
12    """Return (minimum norm, Chern estimate) for periodic KxK vector field."""
13    m = np.asarray(m, dtype=float)
14    norms = np.linalg.norm(m, axis=-1)
15    gap = float(norms.min())
16    if gap <= 1e-14:
17        return (gap, np.nan) if not return_triangles else (gap, np.nan, [])
18    n = m / norms[..., None]
19    K = m.shape[0]
20    total, areas = 0.0, []
21    for i in range(K):
22        for j in range(K):
23            a, b = n[i,j], n[(i+1)%K,j]
24            c, d = n[(i+1)%K,(j+1)%K], n[i,(j+1)%K]
25            for tri in ((a,b,c), (a,c,d)):
26                ar = float(spherical_area(*tri)); total += ar; areas.append(ar)
27    result = (gap, total/(4*np.pi))
28    return result + (areas,) if return_triangles else result
29
30
31def qwz_field(K, mass, anisotropy=1.0):
32    """Periodic response field with gap closings at mass=-2, 0, 2."""
33    x = 2*np.pi*np.arange(K)/K
34    a, t = np.meshgrid(x, x, indexing='ij')
35    return np.stack([np.sin(a), anisotropy*np.sin(t), mass + np.cos(a) + np.cos(t)], axis=-1)
36
37
38def fixed_point_parity(mass):
39    # z signs at (0,0),(pi,0),(0,pi),(pi,pi), with nonzero perturbation avoided.
40    zs = [mass+2, mass, mass, mass-2]
41    return int(np.sign(np.prod(zs)))