# Gradient-flow relaxation test (implicit backward Euler + Newton)

# Gradient-flow relaxation test (implicit backward Euler + Newton) # Enhanced with plotting and line-search fallback # Paste into Colab and run. Requires sympy, mpmath, matplotlib. import sympy as sp import mpmath as mp import json, math, os import matplotlib.pyplot as plt from pathlib import Path # ---------------- Config ---------------- WORKING_DPS = 100 PRINT_DPS = 40 mp.mp.dps = WORKING_DPS # time integration dt_base = mp.mpf('1e-2') max_steps = 200 tol_grad_norm = mp.mpf('1e-60') newton_tol = mp.mpf('1e-80') newton_maxsteps = 40 # eps sweep and perturb magnitudes eps_list = [mp.mpf('1e-8'), mp.mpf('1e-10'), mp.mpf('1e-12'), mp.mpf('1e-14')] perturb_mags = [mp.mpf('1e-6'), mp.mpf('1e-3'), mp.mpf('1e-1')] # Stationary point P* (from your sweep) Pstar = ( mp.mpf('0.138302189376114701435582507174521687560324629234971065651174620741053851929006740376715638479605537213341006758403870868'), mp.mpf('0.0'), mp.mpf('-0.888229085697751498529931541593059458598243214346394910179803351965739736262523806127888689674606780135062631497831338795'), mp.mpf('-0.255267088079590527703267719013850265305802605986153882839029146230166276928917637871453326603341978123494370585096814556') ) # ---------------- Symbolic model ---------------- P_xx, P_xy, P_yx, P_yy = sp.symbols('P_xx P_xy P_yx P_yy', real=True) I1 = P_xx + P_yy I2 = P_xx**2 + P_xy**2 + P_yx**2 + P_yy**2 theta_exact = sp.Rational(11, 10) hyb_pref_exact = sp.Rational(1, 10) * (sp.Rational(1) + sp.Rational(1, 10) * theta_exact)**3 eps_sym = sp.Symbol('eps_reg', positive=True) abs_Pyx = sp.sqrt(P_yx**2 + eps_sym**2) Phi_hyb = ( P_yx + hyb_pref_exact * (I1**2 / (I1**2 + 1)) * (P_yx**2 / (1 + sp.Rational(1, 10) * abs_Pyx)) + sp.Rational(1, 10) * abs_Pyx ) Psi_B = sp.Rational(505, 1000) * I2 + sp.Rational(1, 2) * I1**2 + sp.Rational(25, 1000) * I1**4 + Phi_hyb Psi_sectoral = sp.Rational(2, 5) * P_yy + sp.Rational(3, 80) * P_yy**4 Epot_sym = sp.simplify(Psi_B + Psi_sectoral) vars_order = [P_xx, P_xy, P_yx, P_yy] grad_sym = [sp.diff(Epot_sym, v) for v in vars_order] H_sym = sp.hessian(Epot_sym, vars_order) # ---------------- Evaluation helpers ---------------- def eval_sym_mpf(expr, state, eps_val): subs = { P_xx: sp.Float(str(state[0])), P_xy: sp.Float(str(state[1])), P_yx: sp.Float(str(state[2])), P_yy: sp.Float(str(state[3])), eps_sym: sp.Float(str(eps_val)) } return mp.mpf(str(sp.N(expr.subs(subs), mp.mp.dps))) def grad_at(state, eps_val): return mp.matrix([eval_sym_mpf(g, state, eps_val) for g in grad_sym]) def Hess_at(state, eps_val): return mp.matrix([[eval_sym_mpf(H_sym[i, j], state, eps_val) for j in range(4)] for i in range(4)]) def Epot_at(state, eps_val): return eval_sym_mpf(Epot_sym, state, eps_val) # ---------------- Backward Euler with line-search fallback ---------------- def backward_euler_step(x_old, eps_val, dt_local): def G(a, b, c, d): state = (mp.mpf(a), mp.mpf(b), mp.mpf(c), mp.mpf(d)) g = grad_at(state, eps_val) return ( state[0] + dt_local * g[0] - x_old[0], state[1] + dt_local * g[1] - x_old[1], state[2] + dt_local * g[2] - x_old[2], state[3] + dt_local * g[3] - x_old[3] ) def J(a, b, c, d): state = (mp.mpf(a), mp.mpf(b), mp.mpf(c), mp.mpf(d)) H = Hess_at(state, eps_val) M = mp.matrix(4, 4) for i in range(4): for j in range(4): M[i, j] = dt_local * H[i, j] + (1 if i == j else 0) return M try: sol = mp.findroot(G, x_old, J=J, tol=newton_tol, maxsteps=newton_maxsteps) return tuple(mp.mpf(sol[i]) for i in range(4)), None except Exception as e: return None, str(e) def relax_run(eps, perturb, x0): """Run relaxation with optional dt reduction.""" x_old = tuple(x0) dt = dt_base trace = {'eps': str(eps), 'perturb': str(perturb), 'dt_values': [], 'steps': []} converged = False for step in range(max_steps): Ecur = Epot_at(x_old, eps) gcur = grad_at(x_old, eps) gnorm = mp.sqrt(sum(gcur[i] * gcur[i] for i in range(4))) trace['steps'].append({ 'step': step, 'state': [str(x_old[i]) for i in range(4)], 'E': str(Ecur), 'grad_norm': str(gnorm) }) if gnorm < tol_grad_norm: converged = True break # Try with current dt; if fails, reduce dt x_new, err = backward_euler_step(x_old, eps, dt) if x_new is None: # Try smaller dt dt_small = dt / 10 x_new, err2 = backward_euler_step(x_old, eps, dt_small) if x_new is None: # Try even smaller dt_small2 = dt / 100 x_new, err3 = backward_euler_step(x_old, eps, dt_small2) if x_new is None: trace['error'] = f"Newton failed at step {step}: {err3}" break else: dt = dt_small2 else: dt = dt_small Enew = Epot_at(x_new, eps) if Enew > Ecur + mp.mpf('1e-30'): # Energy increased — accept but warn trace['steps'][-1]['energy_warning'] = f"E increased: {Enew} > {Ecur}" x_old = x_new trace['dt_values'].append(str(dt)) # Final evaluation fin_grad = grad_at(x_old, eps) fin_gnorm = mp.sqrt(sum(fin_grad[i] * fin_grad[i] for i in range(4))) Hfin = Hess_at(x_old, eps) try: eigvals_mat, eigvecs = mp.eigsy(Hfin) eiglist = [mp.mpf(eigvals_mat[r, 0]) for r in range(eigvals_mat.rows)] except Exception: eigvals, eigvecs = mp.eig(Hfin) eiglist = [mp.mpf(v) for v in eigvals] eig_sorted = sorted(eiglist, key=lambda x: abs(x)) trace['final'] = { 'final_state': [str(x_old[i]) for i in range(4)], 'final_grad_norm': str(fin_gnorm), 'eig_sorted': [str(e) for e in eig_sorted], 'steps_taken': step + 1, 'converged': bool(gnorm < tol_grad_norm) } return trace # ---------------- Run experiments ---------------- outdir = "relaxation_runs" os.makedirs(outdir, exist_ok=True) all_results = {} for eps in eps_list: eps_key = str(eps) all_results[eps_key] = [] print(f"\n{'='*60}") print(f"eps = {eps}") print(f"{'='*60}") for perturb in perturb_mags: print(f"\n perturb = {perturb}") x0 = ( Pstar[0] + perturb, Pstar[1] + mp.mpf('0'), Pstar[2] - perturb, Pstar[3] - perturb ) trace = relax_run(eps, perturb, x0) all_results[eps_key].append(trace) print(f" converged: {trace['final']['converged']}") print(f" steps: {trace['final']['steps_taken']}") print(f" final grad norm: {mp.nstr(mp.mpf(trace['final']['final_grad_norm']), 10)}") # ---------------- Plotting ---------------- def plot_relaxation(traces, eps_val, perturb_mags): fig, axes = plt.subplots(1, 2, figsize=(12, 5)) for idx, trace in enumerate(traces): steps = [s['step'] for s in trace['steps']] E_vals = [mp.mpf(s['E']) for s in trace['steps']] grad_vals = [mp.mpf(s['grad_norm']) for s in trace['steps']] # Energy vs step axes[0].semilogy(steps, E_vals, label=f"perturb={perturb_mags[idx]}") # Grad norm vs step axes[1].loglog(steps, grad_vals, label=f"perturb={perturb_mags[idx]}") axes[0].set_xlabel('Step') axes[0].set_ylabel('Energy E(t)') axes[0].set_title(f'Energy Relaxation (eps={eps_val})') axes[0].legend() axes[0].grid(True) axes[1].set_xlabel('Step') axes[1].set_ylabel('Gradient Norm ||∇E||') axes[1].set_title(f'Gradient Norm Decay (eps={eps_val})') axes[1].legend() axes[1].grid(True) plt.tight_layout() plt.savefig(f"{outdir}/relaxation_eps_{str(eps_val)}.png", dpi=150) plt.show() # Generate plots for each eps for eps in eps_list: eps_key = str(eps) plot_relaxation(all_results[eps_key], eps_key, perturb_mags) # ---------------- Save JSON ---------------- Path("relaxation_results.json").write_text(json.dumps(all_results, indent=2)) print("\n✅ JSON saved: relaxation_results.json") print(f"✅ Plots saved in {outdir}/") print("=" * 60)

Popular posts from this blog

THE GOLDEN BALLROOM/BUNKER

Image
Ben Meiselas reports on the shocking admission by Donald Trump’s DOJ in a court case where the DOJ admits to a secret project underneath the ballroom which they claim is needed to protect Donald Trump’s life for “national security purposes.” "You unlock this door with the key of complicity. Beyond it is another dimension — a dimension of betrayal, of indulgence, of fear. You’re moving into a land of both shadow and substance, of politics and paranoia. You’ve just crossed into… the MAGA Zone." "Tonight’s story: A leader sworn to protect his nation makes a bargain with its enemies. The deal? Silence in the face of nuclear annihilation. No retaliation, no defense — only surrender dressed in secrecy. While citizens live unaware, their president builds a palace beneath the earth, a ballroom of gold, of marble and chandeliers, a masquerade hall for billionaires. But behind the gilded doors lies not music and laughter, but a bomb shelter — a sanctuary for the few, pur...
Read more

Conceptual Summary #2: (∂t2​S−c2∇2S+βS3)=σ(x,t)⋅FR​(C[Ψ])

Image
The tests designed to find the Aether—principally the Michelson-Morley experiment (M-M)—would be largely ineffective at detecting the fundamental S field (Substrate) in the Reactive Substrate Theory (RST) framework, though their null results are highly compatible with RST's structure. The reason lies in the distinct conceptual definition and behavior of the two media. The Test: Michelson-Morley and the Aether. The M-M experiment was designed to detect the "aether wind," a predicted change in the speed of light caused by the Earth's motion relative to a stationary, rigid medium—the classical Luminiferous Aether. Classical Aether: A static, absolute reference frame; an elastic, passive, material-like medium that fills space and serves as a carrier for light waves. Predicted Result: Light traveling parallel to Earth's motion through the Aether should be slightly slower than light traveling perpendicular to it. Actual Result (Null Result): No significant differenc...
Read more

ICE PROUDLY ANNOUNCES NEW “ELITE” TASK FORCE COMMANDER JEREMY DEWITTE

Image
For Immediate Release OFFICE OF STRATEGIC SUPPRESSION & DESPERATE RECRUITMENT U.S. IMMIGRATION AND CUSTOMS ENFORCEMENT (ICE) SUBJECT: ICE PROUDLY ANNOUNCES NEW “ELITE” TASK FORCE COMMANDER JEREMY DEWITTE; AGENCY TO ADOPT FULL “MOTOR ESCORT” DOCTRINE BY Q3 WASHINGTON, D.C. — In a historic moment that experts are calling “a cry for help” and insiders are calling “Tuesday,” U.S. Immigration and Customs Enforcement (ICE) has officially appointed Jeremy Dewitte as National Coordinator of the newly created Tactical High-Visibility Fallujah Squadron & Escort Wing (T.H.V.F.S.E.W.) . This appointment follows a PBS NewsHour report highlighting ICE’s bold new hiring strategy: removing every meaningful hiring standard. After eliminating outdated requirements such as "must not be a professional police impersonator," "must possess a valid driver's license," and "must not be a registered sex ...
Read more