Freezer

Freezer#

What’s in this notebook? This notebook demonstrates the module freezer — a framework for freezing out heavy moduli from the EFT. The ConifoldFreezer integrates out the conifold modulus \(z_{\text{cf}}\) by solving its leading-order EOM analytically.

In this notebook, you will learn:

  • How to set up a ConifoldFreezer from an existing FluxEFT or rather FluxVacuaFinder model

  • How the conifold modulus is solved analytically and frozen out

  • How to compute the reduced F-terms \(D_\alpha W_{\text{bulk}}\) in the light-field EFT

  • How to write a custom Freezer subclass for other moduli

Prerequisites: NB10: coniLCS pipeline for conifold geometry setup.

We use the \(h^{1,1}=99\), \(h^{1,2}=3\) example from 2009.03312 throughout.

(Created: March 2026)

Outline#

  1. Setup

  2. Model setup

  3. freezer — freezing out heavy moduli

  4. Theory: \(z_{\rm cf}\) EOM and the \(\widetilde W_1\) corrections

  5. Take-aways

  6. Further reading

Setup#

# General imports
import warnings
import numpy as np
from scipy.optimize import root

# JAX imports
import jax
import jax.numpy as jnp
jax.config.update("jax_enable_x64", True)

# JAXVacua
import jaxvacua as jvc
from jaxvacua.flux_utils import flux_to_pfv, pfv_to_flux, pfv_to_moduli
from jaxvacua.freezer import Freezer, ConifoldFreezer

warnings.filterwarnings('ignore')

Model setup#

We use the \(h^{1,1}=99\), \(h^{1,2}=3\) coniLCS model from 2009.03312, with the conifold curve \([-1,1,0]\).

from cytools import Polytope

poly = Polytope([[-1,3,-2,-1],[1,-1,0,0],[-1,0,0,1],[-1,0,0,0],[-1,0,1,1],[-1,0,2,0],[-1,0,1,0]])
cy = poly.triangulate().get_cy()

basis_matrix = np.array([[0, 1, 1], [1, 1, 0], [0, 0, 1]])
conifold_curve = np.array([-1, 1, 0])

model = jvc.FluxEFT(
    h12=cy.h11(), Q=cy.h11()+cy.h12()+2, 
    use_cytools=True, mirror_cy=cy, ncf=2, use_gvs=True,
    maximum_degree=6, basis_change=basis_matrix, conifold_curve=conifold_curve,
    limit="coniLCS", prange=10, conifold_basis=True
)
model

freezer — freezing out heavy moduli#

The idea#

When a modulus is parametrically heavy, its leading-order EOM determines it as a function of the remaining light fields. We can substitute this solution back into the superpotential and work with an EFT that has fewer degrees of freedom.

The Freezer base class defines the general interface. Subclasses implement:

  • heavy_indices — which moduli to freeze out

  • solve_heavy — how to solve for them given light fields

  • _real_light_to_full — how to convert real coordinates

The base class then provides reconstruct_full_moduli, superpotential, DW_light, DW_x_light, and dDW_x_light for free.

ConifoldFreezer: freezing out \(z_{\text{cf}}\)#

Near the conifold locus, \(z_{\text{cf}}\) acquires a parametrically large mass. Its leading-order EOM gives

\[z_{\text{cf}} = -\frac{1}{2\pi i}\exp\!\left(-\frac{2\pi i\,\widetilde{W}_1}{n_{\text{cf}}(M_1 - \tau H_1)}\right)\]

Documentation note. The setup cell only needs the leading formula. The appendix below fixes the notation used by conifold.zcf_solver, including the logarithmic prefactor, W_log_coeff, and the Kähler-covariant correction. Future notation changes should update the notebook and solver documentation together.

The ConifoldFreezer implements this formula. Let’s set it up:

freezer = ConifoldFreezer(model, conifold_index=0)

print(f"Heavy indices: {freezer.heavy_indices}")
print(f"Light indices: {freezer.light_indices}")
print(f"Freezing out {freezer.n_heavy} modulus, keeping {freezer.n_light} light moduli + tau")
print(f"Conifold degree (sourced from lcs_tree.conifold.ncf): {freezer.ncf}")

Solving for \(z_{\text{cf}}\) from the bulk moduli#

Given the bulk moduli and \(\tau\), the freezer solves for \(z_{\text{cf}}\). The available mode values mirror the underlying compute_zcf dispatcher:

  • mode="manual" is the full closed-form coni-LCS expression. It uses the actual bulk moduli supplied to the freezer and assembles the logarithmic coefficient from the intersection data, \(a\)-matrix, \(b\)-vector, and included worldsheet-instanton terms.

  • mode="pfv" is the perturbatively flat-vacuum / linear-racetrack approximation. It replaces the bulk equations by the analytic PFV racetrack relation expressed directly in the integer flux data.

The ConifoldFreezer is the reduced-EFT wrapper. The direct model.compute_zcf_x(...) calls below are the lower-level dispatcher used by older notebooks and by the freezer internally. They should agree mode by mode; the difference between the full/manual and PFV numbers is a physics approximation, not a wrapper change.

# Use the first solution from Table 1 of arXiv:2009.03312 as the running example.
M = np.array([4, -8, 8])
K = np.array([-8, 3, -6])
tau0 = 1j / 0.38

# PFV initial guess for the bulk moduli.
z0 = pfv_to_moduli(model, M, K, tau0)
zbulk = z0[1:]  # light (bulk) moduli
flux = pfv_to_flux(model, M, K)

# Full/manual route: closed-form coni-LCS expression evaluated at zbulk.
zcf_full = freezer.solve_heavy(zbulk, tau0, flux, mode="manual")

# PFV route: linear-racetrack approximation evaluated from the flux data.
zcf_pfv = freezer.solve_heavy(zbulk, tau0, flux, mode="pfv")

# Cross-check directly against the model's `compute_zcf_x` dispatcher.
# `compute_zcf_x` takes the bulk-only real vector (length 2*h12, no z_cf
# direction); slicing `x[2:]` strips Re(z_cf) / Im(z_cf) from the full real x.
x_full = model._convert_complex_to_real(z0, jnp.conj(z0), tau0, jnp.conj(tau0))
x_bulk = x_full[2:]

zcf_old_pfv = model.compute_zcf_x(x_bulk, flux, mode="pfv")
zcf_old_full = model.compute_zcf_x(x_bulk, flux, mode="manual")

print("Full/manual coni-LCS solve")
print("  ConifoldFreezer:       ", np.abs(zcf_full[0]))
print("  model.compute_zcf_x:   ", np.abs(zcf_old_full))
print()
print("PFV / linear-racetrack approximation")
print("  ConifoldFreezer:       ", np.abs(zcf_pfv[0]))
print("  model.compute_zcf_x:   ", np.abs(zcf_old_pfv))

Reconstructing full moduli and computing the superpotential#

The freezer can reconstruct the full moduli vector (inserting the solved \(z_{\text{cf}}\)) and evaluate EFT quantities on the reduced field space:

# Reconstruct full moduli from light ones
z_full = freezer.reconstruct_full_moduli(zbulk, tau0, flux)
print(f"Full moduli (freezer):  {z_full}")
print(f"Full moduli (PFV):      {z0}")
print()

# Superpotential evaluated with z_cf on-shell
W_reduced = freezer.superpotential(zbulk, tau0, flux)
W_full = model.superpotential(z0, tau0, flux)
print(f"|W| (reduced EFT): {np.abs(W_reduced):.6e}")
print(f"|W| (full PFV):    {np.abs(W_full):.6e}")

Using the freezer for optimisation#

The key use case: solve for the true vacuum using only the bulk moduli as free variables. We first find the full vacuum, then show the freezer reproduces \(z_{\text{cf}}\) at the solution.

# Find the true vacuum using the full model (all moduli free)
x0 = model._convert_complex_to_real(z0, jnp.conj(z0), tau0, jnp.conj(tau0))
res = root(model.DW_x, x0, args=(flux,), jac=model.dDW_x, tol=1e-10, method="hybr")

x_sol = res.x
z_sol, _, tau_sol, _ = model._convert_real_to_complex(x_sol)

print(f"Full minimisation converged: {res.success}")
print(f"sum|DW|: {np.sum(np.abs(model.DW_x(x_sol, flux))):.2e}")
print(f"z_cf (numerical):  {z_sol[0]}")
print(f"|z_cf| (numerical): {np.abs(z_sol[0]):.6e}")

Now we use the freezer to predict \(z_{\text{cf}}\) from the bulk moduli at the solution and compare with the numerics:

# Use the freezer to predict z_cf from the solved bulk moduli
zbulk_sol = z_sol[1:]
zcf_frozen = freezer.solve_heavy(zbulk_sol, tau_sol, flux)

print(f"|z_cf| (numerical):         {np.abs(z_sol[0]):.10e}")
print(f"|z_cf| (freezer, full):     {np.abs(zcf_frozen[0]):.10e}")
print(f"|z_cf| (freezer, PFV):      {np.abs(freezer.solve_heavy(zbulk_sol, tau_sol, flux, mode='pfv')[0]):.10e}")
print()
print(f"Relative error (full):  {np.abs(np.abs(z_sol[0]) - np.abs(zcf_frozen[0])) / np.abs(z_sol[0]):.2e}")

Reduced DW: F-terms for the light moduli only#

The reduced F-terms are the light-modulus components of the full covariant F-term vector after substituting the on-shell heavy modulus. To make the comparison honest, the code below evaluates both the full and reduced quantities at the same frozen point z_frozen_sol, rather than mixing the numerical full root with a separately reconstructed analytic \(z_{\text{cf}}\).

# Full DW at the numerical solution, shown as a reference.
DW_full_numeric = model.DW(z_sol, jnp.conj(z_sol), tau_sol, jnp.conj(tau_sol), flux)

# Reconstruct the frozen full point from the solved light fields and evaluate DW there.
z_frozen_sol = freezer.reconstruct_full_moduli(zbulk_sol, tau_sol, flux, mode="manual")
DW_full_frozen = model.DW(z_frozen_sol, jnp.conj(z_frozen_sol), tau_sol, jnp.conj(tau_sol), flux)

# Reduced DW: keep only the light-modulus components and D_tau W.
light_idx = jnp.array(freezer.light_indices)
DW_light = jnp.append(DW_full_frozen[light_idx], DW_full_frozen[-1])

print("Full DW at the numerical solution (all moduli + tau):")
print(f"  |DW| = {np.abs(DW_full_numeric)}")
print()
print("Full DW at the frozen point (all moduli + tau):")
print(f"  |DW| = {np.abs(DW_full_frozen)}")
print()
print("Reduced DW at the same frozen point (light moduli + tau):")
print(f"  |DW_light| = {np.abs(DW_light)}")
print(f"  Dimensions: {len(DW_full_frozen)} (full) -> {len(DW_light)} (reduced)")

# Quantify why the frozen-point residual is larger than at the full numerical root.
delta_zcf = complex(z_frozen_sol[0] - z_sol[0])
delta_z = np.asarray(z_frozen_sol - z_sol)
delta_DW = np.asarray(DW_full_frozen - DW_full_numeric)

print("Frozen point versus full numerical root:")
print(f"  |Δz_cf|              = {abs(delta_zcf):.3e}")
print(f"  max |Δz|             = {float(np.max(np.abs(delta_z))):.3e}")
print(f"  max |ΔDW|            = {float(np.max(np.abs(delta_DW))):.3e}")
print(f"  max |DW| frozen      = {float(jnp.max(jnp.abs(DW_full_frozen))):.3e}")
print(f"  max |DW| numerical   = {float(jnp.max(jnp.abs(DW_full_numeric))):.3e}")

The reduced system is evaluated at the analytic frozen point, not at the fully numerical root. The analytic solution for \(z_{ m cf}\) is a controlled local approximation, so it shifts the full point by a tiny but non-zero amount. Since the exact numerical root has \(D_IW\) tuned to machine precision, even a very small displacement can lift some light-field F-terms from around \(10^{-15}\) to around \(10^{-8}\). The relevant consistency check is therefore the agreement between DW_x_light and the corresponding light slice of the full DW_x at the same frozen point, not equality with the separate full numerical root.

Real-coordinate interface: DW_x_light and dDW_x_light#

For use with scipy.optimize.root or similar solvers, the freezer also provides real-coordinate versions. These operate on the vector of light fields only, but internally reconstruct the full real vector by solving for \(z_{\text{cf}}\).

The comparison must again be made at the same reconstructed full point. The entries need not match the light slice of model.DW_x(x_sol, flux) if x_sol uses the numerical full root and DW_x_light uses the analytic frozen root.

# Reduced real DW: only bulk moduli + tau (6 components instead of 8).
x_bulk = x_sol[2:]  # drop Re(z_cf), Im(z_cf)
DW_x_light = freezer.DW_x_light(x_bulk, flux, mode="manual")

# Full real DW evaluated at the same frozen point used by DW_x_light.
x_frozen = freezer._real_light_to_full(x_bulk, flux, mode="manual")
DW_x_full_frozen = model.DW_x(x_frozen, flux)
DW_x_expected = DW_x_full_frozen[freezer._real_light_slice]

# Full numerical root, shown only to indicate the nearby-but-different point.
DW_x_full_numeric = model.DW_x(x_sol, flux)

print(f"Full DW_x at numerical solution ({len(DW_x_full_numeric)} components):")
print(f"  |DW_x| = {np.abs(DW_x_full_numeric)}")
print()
print(f"Reduced DW_x_light at frozen point ({len(DW_x_light)} components):")
print(f"  |DW_x_light| = {np.abs(DW_x_light)}")
print(f"  max |DW_x_light - full frozen light slice| = {float(jnp.max(jnp.abs(DW_x_light - DW_x_expected))):.2e}")

Solving the reduced F-term conditions with scipy#

We can directly use DW_x_light and dDW_x_light as objective and Jacobian for scipy.optimize.root, solving only for the 6 bulk + tau real variables instead of 8:

# Initial guess: PFV values for bulk moduli + tau (drop z_cf)
x0_bulk = model._convert_complex_to_real(z0, jnp.conj(z0), tau0, jnp.conj(tau0))[2:]

# Solve the reduced system
res_reduced = root(
    lambda x: freezer.DW_x_light(x, flux),
    x0_bulk,
    jac=lambda x: freezer.dDW_x_light(x, flux),
    tol=1e-10,
    method="hybr"
)

print(f"Reduced minimisation converged: {res_reduced.success}")
print(f"Variables: {len(x0_bulk)} (reduced) vs {len(x0)} (full)")
print(f"sum|DW_x|: {np.sum(np.abs(res_reduced.fun)):.2e}")
print(f"nfev: {res_reduced.nfev}")

Compare the solution from the reduced system with the full minimisation:

# Recover full moduli from the reduced solution
x_red = res_reduced.x
x_full_reconstructed = np.array(freezer._real_light_to_full(x_red, flux))
z_red, _, tau_red, _ = model._convert_real_to_complex(x_full_reconstructed)

print("Bulk moduli comparison (full vs reduced):")
for i in range(1, model.h12):
    print(f"  z[{i}]: {z_sol[i]:.10f}  vs  {z_red[i]:.10f}")

print(f"\ntau: {tau_sol:.10f}  vs  {tau_red:.10f}")
print(f"\n|z_cf| (full):    {np.abs(z_sol[0]):.10e}")
print(f"|z_cf| (reduced): {np.abs(z_red[0]):.10e}")

# Check W0
W_red = model.superpotential(z_red, tau_red, flux, normalise=True)
W_full = model.superpotential(z_sol, tau_sol, flux, normalise=True)
print(f"\n|W0| (full):    {np.abs(W_full):.10e}")
print(f"|W0| (reduced): {np.abs(W_red):.10e}")

Validation across multiple flux choices#

Let’s verify the freezer against all five examples from Table 1 in 2009.03312:

Mlist = np.array([[4,-8,8], [4,-8,10], [8,-12,6], [-8,4,12], [-14,6,27]])
Klist = np.array([[-8,3,-6], [-6,3,-4], [-5,1,-2], [5,1,-4], [4,1,-2]])
gslist = np.array([0.38, 0.15, 0.125, 0.35, 0.0643])

print(f"{'#':>2} | {'converged':>9} | {'|z_cf| num':>14} | {'|z_cf| frozen':>14} | {'rel err':>10} | {'|W0|':>12}")
print("-" * 85)

# Construct FVF model
model_fvf = jvc.FluxVacuaFinder.from_model(model=model)

for i in range(len(Mlist)):
    Mi, Ki, gsi = Mlist[i], Klist[i], gslist[i]
    
    # Build flux and initial guess via flux_utils
    fluxi = model.pfv_to_flux(Mi, Ki)
    tau_i = 1j / gsi
    z_i = model.pfv_to_moduli(Mi, Ki, tau_i)
    
    # Full minimisation
    x_i = model._convert_complex_to_real(z_i, jnp.conj(z_i), tau_i, jnp.conj(tau_i))
    res_i = root(model.DW_x, x_i, args=(fluxi,), jac=model.dDW_x, tol=1e-10, method="hybr")
    
    if not res_i.success:
        # Fallback: Newton + retry
        z_tmp, _, tau_tmp, _ = model._convert_real_to_complex(res_i.x)
        z_tmp, tau_tmp, _ = model_fvf.newton_method_flux_vacua(
            z_tmp, tau_tmp, fluxi, step_size_Newton=0.1, tol=1e-12, max_iters=100, mode="SUSY", solver_mode="real"
        )
        x_i = model._convert_complex_to_real(z_tmp, jnp.conj(z_tmp), tau_tmp, jnp.conj(tau_tmp))
        res_i = root(model.DW_x, x_i, args=(fluxi,), jac=model.dDW_x, tol=1e-8, method="hybr")
    
    z_i_sol, _, tau_i_sol, _ = model._convert_real_to_complex(res_i.x)
    
    # Freezer prediction
    zcf_frozen_i = freezer.solve_heavy(z_i_sol[1:], tau_i_sol, fluxi)
    
    zcf_num = np.abs(z_i_sol[0])
    zcf_frz = np.abs(zcf_frozen_i[0])
    rel_err = np.abs(zcf_num - zcf_frz) / zcf_num
    W0 = np.abs(model.superpotential(z_i_sol, tau_i_sol, fluxi, normalise=True))
    
    print(f"{i:>2} | {str(res_i.success):>9} | {zcf_num:>14.6e} | {zcf_frz:>14.6e} | {rel_err:>10.2e} | {W0:>12.6e}")

Writing a custom Freezer subclass#

The Freezer base class is designed to be extended. To freeze out a different set of heavy moduli, subclass Freezer and implement heavy_indices, solve_heavy, and _real_light_to_full. Here is a minimal skeleton:

from jaxvacua.freezer import Freezer

class MyCustomFreezer(Freezer):
    
    @property
    def heavy_indices(self):
        return (0, 1)  # freeze out the first two moduli
    
    def solve_heavy(self, z_light, tau, fluxes, **kwargs):
        # Your analytic (or numerical) solution for the heavy moduli
        # as functions of z_light, tau, and fluxes
        z_heavy = ...  
        return z_heavy
    
    def _real_light_to_full(self, x_light, fluxes, **kwargs):
        # Convert real light coordinates to full real array
        # by solving for heavy moduli and inserting them
        ...
        return x_full

The base class then automatically provides superpotential, DW_light, DW_x_light, and dDW_x_light for the reduced theory.

Theory: \(z_{\rm cf}\) EOM and the \(\widetilde W_1\) corrections#

This appendix documents the analytic solution that ConifoldFreezer.solve_heavy implements internally. The conventions follow the conifold-PFV setup used in arXiv:2009.03312 and arXiv:2004.10740.

Leading-order conifold \(F\)-term#

Near \(z_{\rm cf} = 0\), the conifold-modulus derivative of the GVW superpotential takes the local form

\[ \partial_{z_{\rm cf}} W_{\rm coni} = (M^1 - \tau H^1)\,\frac{n_{\rm cf}}{2\pi i}\, \ln(-2\pi i\,z_{\rm cf}) + \widetilde W_1 + \mathcal O(z_{\rm cf})\,, \]

where the logarithm originates from the conifold-period monodromy and \(\widetilde W_1\) packages all non-singular dependence on the bulk moduli, the axio-dilaton, and the flux quanta. In the code the two ingredients are accessible separately as

Quantity in code

Math symbol

Method

log_prefactor(tau, flux, conj)

\((M^1 - \tau H^1)\,n_{\rm cf}/(2\pi i)\)

conifold/zcf_solver

W_log_coeff(z, tau, flux, mode, conj)

\(\widetilde W_1\)

conifold/zcf_solver

W_log_coeff has three implementation modes: "manual" (closed-form Li sums), "autodiff" (F_coniLCS_exp derivative), and "pfv" (PFV/linear-racetrack approximation).

Leading-order solution for \(z_{\rm cf}\)#

Dropping the Kähler-covariant term and equating the leading singular + constant pieces, the EOM \(\partial_{z_{\rm cf}} W \simeq 0\) has the closed-form solution

\[ z_{\rm cf}^{(0)} = -\frac{1}{2\pi i}\, \exp\!\left( -\frac{2\pi i\,\widetilde W_1}{(M^1 - \tau H^1)\,n_{\rm cf}} \right)\,. \]

The logarithmic monodromy of the conifold period thus drives the vacuum exponentially close to \(z_{\rm cf} = 0\). The role of \(\widetilde W_1\) is only to set the exponent.

Kähler-covariant correction \(\delta\widetilde W_1\)#

Because \(\partial_{z_{\rm cf}} K\) is generically finite (not zero) at \(z_{\rm cf} = 0\), the full \(F\)-term is

\[ D_{z_{\rm cf}} W = (M^1 - \tau H^1)\,\frac{n_{\rm cf}}{2\pi i}\, \ln(-2\pi i\,z_{\rm cf}) + \widetilde W_1 + \delta\widetilde W_1 + \mathcal O(z_{\rm cf})\,, \]

where

\[ \delta\widetilde W_1 \coloneqq \left.(\partial_{z_{\rm cf}} K)\right|_{z_{\rm cf}=0}\, W_{\rm bulk}\,. \]

The corrected solution then reads

\[ z_{\rm cf} = -\frac{1}{2\pi i}\, \exp\!\left( -\frac{2\pi i\,(\widetilde W_1 + \delta\widetilde W_1)} {(M^1 - \tau H^1)\,n_{\rm cf}} \right)\,. \]

In code: the additive piece \(\delta\widetilde W_1\) is exposed as log_coeff_K_corr(z, cz, tau, ctau, flux, conj) and is included by passing apply_correction=True to compute_zcf* and solve_heavy. For coni-PFVs with tuned \(W \ll 1\), \(\delta\widetilde W_1\) is parametrically suppressed — the correction is small but non-zero. For generic AFVs with \(W \sim \mathcal O(1)\), \(\delta\widetilde W_1\) should not be neglected.

What freezer.solve_heavy(mode=...) does internally#

  z_bulk, tau, fluxes
          │
          ▼
  log_prefactor(tau, flux)                       ←  conifold/zcf_solver
          │
          │
          ▼
  W_log_coeff(z_bulk, tau, flux, mode=...)       ←  3 modes available
          │                                          •  "manual"   — closed-form Li sums
          │                                          •  "autodiff" — F_coniLCS_exp gradient
          │                                          •  "pfv"      — racetrack approximation
          │
          ▼
  [if apply_correction]:
      W_log_coeff += log_coeff_K_corr(z, cz, tau, ctau, flux)
          │
          ▼
  _zcf_from_log_coeff(W_log_coeff, tau, flux)
          │
          │      z_cf = -1/(2πi) · exp(-W_log_coeff / log_prefactor)
          ▼
  return z_cf

All three mode flavours give the same \(z_{\rm cf}\) to machine precision (manual ≡ autodiff are exact reformulations; "pfv" differs by deliberate approximations and is the cheapest of the three).

Take-aways#

  • The Freezer framework provides a uniform interface for integrating out heavy moduli analytically. ConifoldFreezer is the canonical implementation; the base class is designed for subclassing.

  • The conifold modulus EOM reduces to a closed form: \(z_{\rm cf} = -(2\pi i)^{-1} \exp\!\left(-2\pi i\,\widetilde W_1 / [(M^1 - \tau H^1) n_{\rm cf}]\right)\). The logarithmic monodromy of the conifold period is what makes \(z_{\rm cf}\) exponentially small in \(\widetilde W_1\).

  • \(\widetilde W_1\) packages the non-singular flux/moduli dependence. The Kähler-covariant correction \(\delta\widetilde W_1 = (\partial_{z_{\rm cf}} K) W_{\rm bulk}\) is small when \(W\) is tuned (coni-PFVs) and \(\mathcal O(1)\) otherwise.

  • The Freezer.DW_x_light / dDW_x_light interface lets you solve for SUSY vacua in the reduced bulk EFT (one fewer modulus). Newton refinement converges in fewer iterations because \(z_{\rm cf}\) is already at its true value.

  • solve_heavy(mode="manual" | "autodiff" | "pfv") controls the implementation of \(\widetilde W_1\); the first two agree to machine precision, the third is a deliberately cheaper approximation.

Further reading#