1import numpy as np
 2from degree_phase_monitor import jacobian, integrate, order_parameter
 3
 4def field(m,k,p,g):
 5    u=(p*k)@m/(p@k)
 6    return -m+np.tanh(g*k*u)
 7
 8# Prediction 1: at zero, exact collective eigenvalue crosses at g_c=<k>/<k^2>.
 9# Prediction 2: equal degrees have no degree separation, while heterogeneity lowers g_c.
10# Prediction 3: degree penalty lowers S with less mean-scale distortion than global penalty.
11for d in [0.,1.,2.,4.]:
12    k=np.array([3-d,3+d]); p=np.array([.5,.5]); gc=(p@k)/(p@(k*k))
13    A=jacobian(np.zeros(2),k,p,.37)
14    eps=1e-6; N=2; FD=np.column_stack([(field(eps*np.eye(N)[j],k,p,.37)-field(-eps*np.eye(N)[j],k,p,.37))/(2*eps) for j in range(N)])
15    err=np.max(np.abs(A-FD))
16    print('d',d,'pred_gc',gc,'fd_jac_err',err,'eig_at_gc',max(np.linalg.eigvals(jacobian(np.zeros(2),k,p,gc)).real))