Sampling module

Sampling module#

What’s in this notebook? This notebook illustrates the sampling infrastructure within JAXVacua. The data_sampler class provides flexible methods for generating initial guesses for the moduli and fluxes, which serve as starting points for the numerical minimisation described in notebooks NB06 and NB07.

In this notebook, you will learn:

How to construct a data_sampler and understand its three independent bound specs:

moduli_bounds — per-direction box on \(\operatorname{Im}(z_i)\) (scalar or length-\(h^{1,2}\) vectors).

stretching + exclude_walls — Kähler-cone wall distance, with optional per-wall exemptions (useful for coniLCS setups).

cone_cutoff — scalar \(L^2\) cap used by cone / ray sampling modes.

The seven moduli sampling modes: "box", "cone", "stretched_cone", "tip_ray", "random_ray", "random_rays", "sphere".

The one-shot helper initial_guesses plus its filter knob filter_moduli=True (physics-aware pruning via filter_by_instantons + filter_by_km).

The component samplers (get_moduli, get_axions, get_dilaton, get_fluxes) and composed helpers (get_complex_moduli, get_complex_tau).

The standalone filter API (filter_by_instantons, filter_by_km, filter_moduli, filter_points) and cone-geometry helpers (find_interior_points, sample_ray, rescale_points).

How ISD sampling works: replacing a subset of the fluxes by solving the ISD condition \(\star G_3 = \mathrm{i}\, G_3\) for given moduli, plus the four ISD modes ("H", "F", "ISD+", "ISD-").

Prerequisites: Basic JAXVacua usage (NB02), flux vacua finding (NB05).

Related notebooks: NB06: ISD sampling principle (single-vacuum ISD example), NB07: ISD sampling workflows (large-scale moduli-first and flux-first sampling).

(Created: Andreas Schachner, June 25, 2024; major expansion 2026-04-24 to cover the direction-aware bounds refactor.)

How to use this notebook#

The sampler has many knobs, but most workflows use one of three entry points:

Goal	Recommended method	Main controls
Draw complete starting data	`initial_guesses(...)`	`moduli_sampling_mode`, `fluxes_sampling_mode`, `filter_moduli`, `include_fluxes`
Inspect one component	`get_moduli`, `get_complex_tau`, `get_fluxes`	Bounds, cone/ray options, tadpole or sphere constraints
Seed ISD-based searches	`initial_guesses_ISD(...)`	ISD mode, integer rounding, flux bounds and physical filters

Read §3 first if you want the high-level API, then use §4-§6 as a reference for individual samplers and filters. The later sections are intentionally more detailed because they document the edge cases needed by the vacuum-finding notebooks.

1. Imports#

# General imports
import warnings
import time
import numpy as np
from tqdm.auto import tqdm

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

# Plotting
import seaborn as sn
import matplotlib.pyplot as plt
cmap = sn.color_palette("viridis", as_cmap=True)

# JAXVacua
import jaxvacua as jvc

warnings.filterwarnings('ignore')

2. Model and sampler setup#

All the examples below use a small local model at \(h^{1,2}=2\) so you can run the notebook end-to-end on a laptop. Replace the model-loading cell with any jvc.FluxVacuaFinder / jvc.FluxEFT of your choosing — the sampler interface is identical.

2.1 Scalar `moduli_bounds` (backward-compatible)#

The simplest way to constrain the sampler is to pass a scalar pair (lo, hi). Each complex-structure modulus direction is then sampled from the same box \(\operatorname{Im}(z_i)\in[\texttt{lo}, \texttt{hi}]\) (in "box" mode), and the same pair also sets the default \(L^2\) cap for cone sampling (via cone_cutoff=max(moduli_upper)=hi auto-resolution).

Load model from files

#if False:
if True:
    h12=2
    model = jvc.FluxVacuaFinder(h12 = h12, Q=276, model_ID = 1, maximum_degree = 2,limit="LCS",model_type="KS")

Alternatively, use CYTools

from cytools import Polytope, Cone, fetch_polytopes, read_polytopes

if False:
    p = fetch_polytopes(h11=2,h12=272,limit=5,lattice="N",as_list=True)[0]
    cy = p.triangulate().get_cy()
    mcap = cy.mori_cone_cap(in_basis=True)
    Kcup = mcap.dual_cone()
    basis_change = Kcup.extremal_rays()

    model = jvc.FluxEFT(h12=cy.h11(), Q=cy.h11()+cy.h12()+2, model_type="KS", maximum_degree=2, 
                                      use_cytools=True, mirror_cy = cy, basis_change=basis_change)

As objective function, we take the \(F\)-terms for the moduli as computed in model.DW. For later purposes, we vectorise this function by using jax.vmap:

DW_v = vmap(model.DW)

It is also convenient to introduce a data sampler that constrains our sampling procedure to specific regions in moduli and flux space:

sampler = jvc.data_sampler(model,flux_bounds=[-5,5],moduli_bounds=[0,5],axion_bounds=[-0.5,0.5])

2.2 Direction-aware `moduli_bounds`#

In realistic CY geometries, different moduli directions often have very different natural scales — one direction might want \(\operatorname{Im}(z_1)\in(0.5, 2)\) while another wants \(\operatorname{Im}(z_2)\in(10, 100)\). A single scalar pair cannot express this. moduli_bounds therefore accepts a 2-tuple of length-\(h^{1,2}\) vectors as an alternative:

moduli_bounds = (lower_vec, upper_vec)   # element-wise per-direction bounds

After construction, sampler.moduli_lower / sampler.moduli_upper are always JAX arrays of shape (h12,), regardless of whether the user passed scalars (broadcast) or per-direction vectors. Legacy user code that does np.all(im_z >= sampler.moduli_lower) continues to work transparently via broadcasting.

# Scalar form — every direction in [2, 5]
s_scalar = jvc.data_sampler(
    model, moduli_bounds=(2., 5.),
    axion_bounds=(-0.5, 0.5),
    dilaton_bounds=(0.87, 10.),
)
print("scalar:        lower =", s_scalar.moduli_lower,
      "  upper =", s_scalar.moduli_upper)

# Per-direction form — z_1 in [1, 3], z_2 in [2, 5]
s_dir = jvc.data_sampler(
    model, moduli_bounds=([1., 2.], [3., 5.]),
    axion_bounds=(-0.5, 0.5),
    dilaton_bounds=(0.87, 10.),
)
print("per-direction: lower =", s_dir.moduli_lower,
      "  upper =", s_dir.moduli_upper)

2.3 Cone knobs: `stretching`, `exclude_walls`, `cone_cutoff`#

Three additional kwargs on the constructor control cone-mode sampling. They have no effect in "box" mode — the sampler simply doesn’t consult them.

Knob	Effect	Default
`stretching`	Distance from every Kähler-cone wall enforced in cone / ray sampling. Implemented as a rejection filter (`H @ Im(z) ≥ stretching`).	`0.0` (bare cone membership).
`exclude_walls`	Per-wall exemption. Accepts a list of integer indices (into `lcs_tree.hyperplanes`) or a list of row vectors. Excluded walls require only bare cone membership; all others enforce the full `stretching` distance. Designed for coniLCS setups where you want to be close to one specific wall.	`None` (all walls obey `stretching`).
`cone_cutoff`	Scalar \(L^2\) cap on sampled cone points and on ray-coefficient ranges. Replaces the implicit `moduli_upper` cap used before the 2026-04-24 refactor.	`float(max(moduli_upper))` (keeps back-compat).

An internal attribute sampler.cone_opted_in is auto-computed from stretching != 0 OR any excluded wall, and controls which checks the downstream bounded_fluxes._in_patch applies. It’s exposed read-only so the user can verify the configuration.

# coniLCS-style sampler: stay 1 away from every wall except wall 0
# (which we allow to be approached arbitrarily closely).
s_coni = jvc.data_sampler(
    model, moduli_bounds=(0., 3.),
    stretching=1,
    exclude_walls=[0],        # wall 0 gets bare cone membership only
    cone_cutoff=10.,          # scalar L² cap on sampled points
)
print(f"stretching    = {s_coni.stretching}")
print(f"exclude_walls = {s_coni.exclude_walls}")
print(f"cone_cutoff   = {s_coni.cone_cutoff}")
print(f"cone_opted_in = {s_coni.cone_opted_in}")

# Visualise three cone configurations side by side.
import numpy as _np
H = _np.asarray(model.lcs_tree.hyperplanes)

fig, ax = plt.subplots(1, 3, figsize=(12, 3.5), dpi=200, sharex=True, sharey=True)

for a, (s, title) in zip(ax, [
    (jvc.data_sampler(model, moduli_bounds=(0., 3.)),
     "stretching=0, no exclusions"),
    (jvc.data_sampler(model, moduli_bounds=(0., 3.), stretching=1.),
     "stretching=1, all walls"),
    (jvc.data_sampler(model, moduli_bounds=(0., 3.), stretching=1.,
                      exclude_walls=[0]),
     "stretching=1, excl. wall 0"),
]):
    
    mod = s.get_moduli(5000, sampling_mode="cone",
                        stretching=s.stretching,
                        n_rays=2, use_rays=True)
    a.scatter(_np.asarray(mod[:, 0]), _np.asarray(mod[:, 1]),
                s=4, alpha=0.6, color="steelblue")
    
    # Overlay Kähler-cone walls as dashed lines through the origin
    a.set_title(title, fontsize=9)
    a.set_xlabel(r"Im$(z^1)$")
    a.set_ylabel(r"Im$(z^2)$")
    a.set_aspect("equal")


plt.show()

3. `initial_guesses` — the one-shot API#

3.1 Minimal example#

sampler.initial_guesses(N) returns a 3-tuple (moduli, tau, fluxes) of \(N\) ready-to-go starting points:

moduli — complex array of shape (N, h12)
tau — complex array of shape (N,)
fluxes — integer array of shape (N, 2·n_fluxes)

Pass include_fluxes=False to return just (moduli, tau). For reproducibility pass a JAX PRNG key via rns_key= or simply np.random.seed(42) at the top of the cell (the default NumPy backend of data_sampler honours this).

np.random.seed(42)  # reproducible starting points

The example below draws 5 starting points from the sampler built in §2.1:

moduli,tau,fluxes = sampler.initial_guesses(5)
moduli,tau,fluxes

3.2 Keyword overrides#

Every stored bound can be overridden per call via minval_* / maxval_* kwargs. This is useful when you want most of the sampler’s defaults but need to zoom into a tight sub-region for one particular experiment (e.g. a parameter scan). minval_moduli / maxval_moduli accept scalars or length-h12 sequences — same flexible form as the constructor.

# Zoom into a tight per-direction sub-region for this call only.
np.random.seed(0)
mod, tau = sampler.initial_guesses(
    5,
    minval_moduli=[1., 2.],   # per-direction lower
    maxval_moduli=[2., 3.],   # per-direction upper
    moduli_sampling_mode="box",
    include_fluxes=False,
)
print("moduli imag:")
print(jnp.imag(mod))

3.3 Filter-aware sampling#

Pass filter_moduli=True to drop sampled moduli points that fail two physical requirements:

Kähler-metric positivity — the Kähler-metric eigenvalues at \((z, \tau)\) must be positive-definite so the kinetic terms of the moduli have the right sign. Implemented as filter_by_km (see §5.2).
Instanton control — the instanton-corrected prepotential must stay below a ratio cutoff (inst_cutoff, default 1e-1) relative to the perturbative part. This keeps the samples in the perturbative LCS regime where jaxvacua’s analytic prepotential is a good approximation. Implemented as filter_by_instantons (see §5.1).

The two filters are combined into filter_moduli (§5.3). Rejection rates for a small h12=2 geometry are typically in the 30–60% range; larger h12 usually sees higher retention because the cone has more interior volume relative to the walls.

# Retention-rate demo: compare filtered vs. unfiltered sampling.
np.random.seed(1)
N_req = 2000

mod_raw, _ = sampler.initial_guesses(N_req, include_fluxes=False)
mod_filtered, _ = sampler.initial_guesses(
    N_req, include_fluxes=False, filter_moduli=True,
)

print(f"Unfiltered:   {mod_raw.shape[0]:>5d} points")
print(f"filter_moduli={True}: {mod_filtered.shape[0]:>5d} points")
print(f"(Note: initial_guesses always returns exactly N points by "
      f"re-sampling until the target is met — retention shows up as "
      f"wall-clock cost, not as a smaller array.)")

sampler.update_interior_points(num_pts = 2000)

4. Component samplers#

initial_guesses is a thin convenience wrapper that composes five finer-grained samplers:

get_moduli — the imaginary part of \(z\).
get_axions — the real part of \(z\).
get_dilaton — \(\tau\) (in the standard SL(2,ℤ) fundamental domain).
get_fluxes — integer flux vectors.
get_complex_moduli / get_complex_tau — glue the pieces back into complex scalars.

If you want finer control over any individual component — e.g. sample moduli from the Kähler cone but fluxes from an L²-sphere — you can bypass initial_guesses and call these directly.

4.1 `get_moduli` and its 7 sampling modes#

Mode	Samples	Remarks
`"box"`	\(\operatorname{Im}(z_i) \sim U[\texttt{moduli\_lower}[i],\texttt{moduli\_upper}[i]]\) per direction	One-shot, no filter. Does not enforce cone membership.
`"cone"`	Positive combination of cone generators (or interior-point perturbation)	Filter applied only when `stretching > 0` or walls are excluded.
`"stretched_cone"`	`cone` + offset `stretching · tip_skc`	Equivalent to `cone` with the same stretching for walls that satisfy `H_k · tip ≥ 1`.
`"random_rays"`	Like `"cone"` but sampling random rays from the full ray set	Only available when `lcs_tree.rays_kahler_cone` is non-empty.
`"tip_ray"`	Coefficient times `tip_skc` (single direction)	Useful as a deterministic probe.
`"random_ray"`	Coefficient times a single randomly-chosen extremal ray	Rejection-filtered to stay inside the cone.
`"sphere"`	Uniform on an \(L^2\)-sphere of radius `cone_cutoff`	Useful for pure norm-constrained sweeps.

The cells below sweep the main modes with a scatter plot for a 2-modulus example. The first uses the default (linear-programming) path; the second applies rescale_points to enforce a common norm; the third turns on use_rays=True so positive combinations of cone generators are used instead of interior perturbation.

As a quick illustration, we can compare the following ways of sampling moduli. The default option for sampling modes "cone" and "stretched_cone" is using a linear solver from gurobi to solve the hyperplane constraints \(H\cdot x\geq 0\).

sampling_modes = ["box","cone","stretched_cone","tip_ray","random_ray","random_rays"]
fig = plt.figure(dpi=150,figsize=(8,6))
for sampling_mode in sampling_modes:
    
    stretching = 0.
    if sampling_mode=="stretched_cone":
        stretching = 1.
        
    moduli = sampler.get_moduli(1000, 
                    sampling_mode = sampling_mode,
                    stretching  = stretching,
                    n_rays = 2,
                    perturbation = 1,
                    use_rays = False)
    
    sn.scatterplot(x=moduli[:,0],y=moduli[:,1],label=sampling_mode)
    
plt.xlabel(r"Im$(z^1)$")
plt.ylabel(r"Im$(z^2)$")
plt.title(r"Sampling using linear programming")
plt.show()

The rescale_points(pts, norm, maxval) helper projects sampled points onto an \(L^1\), \(L^2\), or \(L^\infty\) ball of radius maxval. It’s independent of the sampling mode: you can draw points with any get_moduli(sampling_mode=…) and then reshape them through whichever norm matches the physics you care about.

The three norms differ in how they treat anisotropic samples:

\(L^2\) (\(|x|_2 \le R\)) — isotropic Euclidean cap, rotationally symmetric. Useful when “overall length of the moduli vector” is the natural small parameter (e.g. distances in the Weil–Petersson metric).
\(L^\infty\) (\(\max_i |x_i| \le R\)) — axis-aligned cap. Equivalent to clipping each component independently; the ball is a hypercube. Useful when the per-direction box bounds are already your physical constraint and you want a post-hoc rescaling that preserves that geometry.
\(L^1\) (\(\sum_i |x_i| \le R\)) — cross-polytope cap. Biased toward axis-aligned samples; often too aggressive in squeezing interior cone points. The three plots below show the same get_moduli output rescaled under each norm at maxval=2. Compare the footprint each norm produces — you’ll see the characteristic circular / square / diamond shapes corresponding to the \(L^2\) / \(L^\infty\) / \(L^1\) unit balls.

Note. Sampler-level cone_cutoff always uses the \(L^2\) norm by convention. The rescale_points helper is the one-off override when you want a different norm for a specific analysis.

sampling_modes = ["box","cone","stretched_cone","tip_ray","random_ray","random_rays"]

for norm in ["l1","l2","inf"]:
    fig = plt.figure(dpi=100,figsize=(6,4))
    for sampling_mode in sampling_modes:

        stretching = 0.
        if sampling_mode=="stretched_cone":
            stretching = 1.

        moduli = sampler.get_moduli(1000, 
                        sampling_mode = sampling_mode,
                        stretching  = stretching,
                        n_rays = 2,
                        perturbation = 1,
                        use_rays = False)
        
        moduli = sampler.rescale_points(moduli,norm=norm,maxval=2)

        sn.scatterplot(x=moduli[:,0],y=moduli[:,1],label=sampling_mode)

    plt.xlabel(r"Im$(z^1)$")
    plt.ylabel(r"Im$(z^2)$")
    plt.title(r"Sampling using "+norm+" norm")
    plt.show()

We can also use rays for sampling modes "cone" and "stretched_cone":

sampling_modes = ["box","cone","stretched_cone","tip_ray","random_ray","random_rays"]
fig = plt.figure(dpi=150,figsize=(8,6))
for sampling_mode in sampling_modes:
    
    stretching = 0.
    if sampling_mode=="stretched_cone":
        stretching = 1.
        
    moduli = sampler.get_moduli(1000, 
                    sampling_mode = sampling_mode,
                    stretching  = stretching,
                    n_rays = 2,
                    perturbation = 1,
                    use_rays = True)
    
    sn.scatterplot(x=moduli[:,0],y=moduli[:,1],label=sampling_mode)
    
plt.xlabel(r"Im$(z^1)$")
plt.ylabel(r"Im$(z^2)$")
plt.title(r"Sampling using rays")
plt.show()

4.2 `get_axions` and `get_dilaton`#

get_axions(N) draws the real parts of the complex structure moduli from \(U[\texttt{axion\_lower}, \texttt{axion\_upper}]\). get_dilaton(N) draws the axio-dilaton from the SL(2,ℤ) fundamental domain — the real part from \(U[-0.5, 0.5]\) and the imaginary part from \([\texttt{s\_lower}, \texttt{s\_upper}]\), subject to \(|\tau| \geq 1\).

np.random.seed(2)

axions = np.asarray(sampler.get_axions(2000))
tau_samples = np.asarray(sampler.get_complex_tau(2000))

fig, ax = plt.subplots(1, 3, figsize=(11, 3), dpi=120)
ax[0].hist(axions[:, 0], bins=40, color="slategrey")
ax[0].set_title(r"$\mathrm{Re}(z^1)$ (axion)"); ax[0].set_xlabel("value")
ax[1].hist(tau_samples.real, bins=40, color="coral")
ax[1].set_title(r"$\mathrm{Re}(\tau)$"); ax[1].set_xlabel("value")
ax[2].hist(tau_samples.imag, bins=40, color="seagreen")
ax[2].set_title(r"$\mathrm{Im}(\tau) = s$"); ax[2].set_xlabel("value")
plt.tight_layout()
plt.show()

4.3 `get_fluxes`#

get_fluxes produces integer flux vectors of length \(2\cdot n_{\rm fluxes}\). Three sampling_mode options:

`sampling_mode`	What it does
`"box"` (default)	Each flux component \(\sim U[\texttt{flux\_lower}, \texttt{flux\_upper}]\) (rounded).
`"sphere"`	Uniform on the integer \(L^2\)-sphere of a given `radius`.
`"none"`	Return zero fluxes (useful for debugging ISD flows).

Additionally, mode="tadpole" triggers a rejection loop that keeps only fluxes whose tadpole \(|f^T \Sigma h|\) stays below the given radius — convenient when you want pre-filtered starting points for an ISD search.

np.random.seed(3)

fl_box    = sampler.get_fluxes(5, sampling_mode="box")
fl_sphere = sampler.get_fluxes(5, sampling_mode="sphere", radius=10.)

print("box samples (5 rows of length 2·n_fluxes):")
print(np.asarray(fl_box))
print("\nsphere samples (same shape, points on an L²-sphere of radius 10):")
print(np.asarray(fl_sphere))
print("\nsphere norms:")
print(np.linalg.norm(np.asarray(fl_sphere), axis=1))

4.4 `get_complex_moduli`, `get_complex_tau`#

These two helpers glue the component samplers into complex scalars:

get_complex_moduli(N) returns \(z_i = c_i + \mathrm{i}\, y_i\) with \(c_i\) from get_axions and \(y_i\) from get_moduli.
get_complex_tau(N) returns \(\tau = c_0 + \mathrm{i}\, s\) from get_dilaton.

initial_guesses uses them internally.

np.random.seed(4)
z_c    = sampler.get_complex_moduli(3)
tau_c  = sampler.get_complex_tau(3)
print("complex moduli (shape (3, h12)):\n", np.asarray(z_c))
print("\ncomplex tau (shape (3,)):\n", np.asarray(tau_c))

5. Filter API#

The sampler exposes four named filters. They can be composed in any of three ways:

Via the filter_moduli=True kwarg on initial_guesses / initial_guesses_ISD — combines filter_by_km ∧ filter_by_instantons.
Via the filter= callback on get_moduli — attaches a user-supplied filter to the cone-sampling rejection loop.
Standalone on a pre-computed (moduli, tau) batch.

5.1 `filter_by_instantons`#

Drops moduli points whose instanton correction to the prepotential exceeds a fraction inst_cutoff (default 1e-1) of the perturbative part:

\[ |F_{\rm inst}(z)| \leq \texttt{inst\_cutoff} \cdot |F_{\rm pert}(z)|. \]

Useful to stay in the regime where jaxvacua’s analytic prepotential is a good approximation.

np.random.seed(5)

z_raw, tau_raw = sampler.initial_guesses(500, include_fluxes=False)
cutoffs = [1e-3, 1e-2, 1e-1, 1., 10.]
retention = []
for c in cutoffs:
    mask = sampler.filter_by_instantons(z_raw, inst_cutoff=c)
    retention.append(float(np.asarray(mask).mean()))
    print(f"inst_cutoff={c:>6.0e}:  kept {int(mask.sum())}/500 "
          f"({100*retention[-1]:5.1f}%)")

plt.figure(figsize=(5, 3), dpi=120)
plt.semilogx(cutoffs, [100*r for r in retention], "-o")
plt.xlabel("inst_cutoff"); plt.ylabel("retention (%)")
plt.title("filter_by_instantons retention")
plt.grid(True, alpha=0.3)
plt.tight_layout(); plt.show()

5.2 `filter_by_km`#

Drops moduli points where the Kähler metric on the complex structure sector is not positive-definite at \((z, \tau)\). The metric is computed via model.kahler_metric(z, conj(z)) (vectorised over a batch). A point fails if any eigenvalue is ≤ 0, i.e. some kinetic-term sign is wrong — physically unacceptable starting point for Newton refinement.

np.random.seed(6)
z_raw, tau_raw = sampler.initial_guesses(500, include_fluxes=False)
mask_km = sampler.filter_by_km(z_raw, tau_raw)
print(f"filter_by_km: kept {int(mask_km.sum())}/500 "
      f"({100*float(mask_km.mean()):.1f}%)")

5.3 `filter_moduli` (compound)#

filter_moduli is the logical AND of the two filters above. It returns the filtered subset of moduli (not a mask). This is the callable wired into initial_guesses(filter_moduli=True):

filtered = sampler.filter_moduli(z_raw, tau_raw, inst_cutoff=1e-1)

5.4 `filter_points` (cone membership + `exclude_walls`)#

filter_points is the Kähler-cone hyperplane filter — orthogonal to the instanton / Kähler-metric filters above. It’s used inside get_moduli when stretching > 0 or any wall is excluded. The per-wall threshold is

\[\begin{split} \texttt{threshold}_k = \begin{cases} 0 & \text{if } k \in \texttt{exclude\_walls}, \\ \texttt{stretching} & \text{otherwise.} \end{cases} \end{split}\]

A point \(x\) passes iff \(H_k \cdot x \geq \texttt{threshold}_k\) for every wall \(k\). See §2.3 for the plot.

6. Geometry helpers#

A handful of lower-level methods expose the Kähler-cone geometry that the cone-sampling modes use internally. You rarely need them for routine work, but they’re handy for diagnostics.

6.1 `find_interior_points` and `update_interior_points`#

sampler._cone_points is a cached pool of ~1000 interior points constructed once at sampler initialisation. It’s the starting set for get_moduli(sampling_mode="cone") when use_rays=False (perturbation-based mode). You can rebuild it at any time with a larger pool via update_interior_points:

sampler.update_interior_points(num_pts=2000)

Or inspect a fresh interior set without caching via:

pts = sampler.find_interior_points(N=100, stretching=0.1)

The stretching parameter here sets the initial solver offset used when building the cone interior.

6.2 `sample_ray`#

sampler.sample_ray() returns a single random ray from the sampler’s cone generator set (_extremal_rays or _rays). It’s the primitive behind "random_ray" sampling. Useful if you want to do your own scan over a single ray direction.

6.3 `rescale_points`#

Rescale a batch of points so that their L¹, L², or L∞ norm is bounded by maxval. Default maxval=None → uses sampler.cone_cutoff (the dedicated scalar cap). This is the only place cone_cutoff is consumed:

pts_rescaled = sampler.rescale_points(pts, norm="l2")   # uses cone_cutoff
pts_linf     = sampler.rescale_points(pts, norm="inf", maxval=3.0)

The cell below demonstrates how applying different norms with rescale_points re-shapes a common input into different regions:

7. ISD sampling#

General principle of ISD sampling#

The general concept behind ISD sampling is that, instead of randomly choosing integer fluxes from some distribution, a subset of the fluxes can be replaced in favour of fixing values for the moduli, see in particular 2306.06160. In the process, the remaining fluxes can be fixed through the ISD condition \(\star G_3 = \text{i}G_3\) where the Hodge-\(\star\) depends explicitly on the moduli. A clear advantage is that the ISD condition is linear in the fluxes and typically can be easily solved for given inputs. (In contrast, when fixing all fluxes and solving the ISD condition for the moduli, one has to deal with highly non-trivial coupled equations in the moduli which can typically only be solved numerically.) However, one immediate disadvantage is that the fluxes sampled through the ISD condition are typically not quantised. But even after rounding the fluxes to integers, the chosen parameters approximately solve the ISD condition and hence serve as valuable starting guesses in a numerical search for flux vacua.

Below, we demonstrate the above idea in one single example before setting up a larger scan further below. Let us choose initial guesses z0 and tau0 for the moduli \(z^i_0\) and the axio-dilaton \(\tau_0\) together with a choice Hflux of NSNS flux quanta \(h\)

# Moduli starting guesses
z0 = jnp.array([0.3+3j , 0.36+3.1j])

# Axio-dilaton starting guess
tau0 = -0.3+6.7j

# Choices of H-fluxes
Hflux = jnp.array([39., -13.,  -4.,   0.,  -0.,  0.])

Then, the RR-fluxes \(f\) at the minimum are specified through the following version of the ISD condition

\[ f=(s\, M(z_0^i,\overline{z}_0^i)\Sigma + c_0)\, h\; ,\quad \tau_0=c_0 + \text{i} s \, . \]

This particular version of ISD sampling was employed in 2501.03984. Here, the ISD -matrix \(M\) is computed from model.ISD_matrix:

s = tau0.imag
c0 = tau0.real
M = model.ISD_matrix(z0,jnp.conj(z0))

Fflux = jnp.matmul(M,jnp.matmul(model.periods.sigma,Hflux))*s+c0*Hflux

Fflux.real

As expected, the fluxes obtained in this way are not quantised. For these values of the fluxes, the ISD condition or equivalently \(D_IW=0\) is satisfied:

model.DW(z0,jnp.conj(z0),tau0,jnp.conj(tau0),jnp.append(Fflux.real,Hflux))

To find an actual flux vacuum, we now round this choice to integers leading to

Fflux_rounded = jnp.around(Fflux.real,0)
fluxes0 = jnp.append(Fflux_rounded,Hflux)
fluxes0

The above steps have been collected in a single wrapper function model.ISD_sampling (see below for additional details)

fluxes0 = sampler.ISD_sampling(z0,jnp.conj(z0),tau0,jnp.conj(tau0),Hflux,mode="H",output="full",return_integer_flux=True).real
fluxes0

However, since we changed the choices of RR-fluxes, the initial guesses z0 and tau0 do not correspond to the actual points in moduli space at which the scalar potential is minimised as can be seen by computing \(D_I W\):

model.DW(z0,jnp.conj(z0),tau0,jnp.conj(tau0),fluxes0)

Nonetheless, the initial guesses z0 and tau0 are typically very close to an actual solution to \(D_I W=0\) which we can e.g. find by employing Newton’s method.

This is illustrated in the notebook 06_ISD_sampling_principle.ipynb

Different versions of ISD sampling#

Lastly, let us describe the four modes for solving the ISD condition when half of the fluxes are specified. Either we use the following form of the ISD condition

\[ f_1-\tau h_1=\overline{\mathcal{N}}(z^i,\overline{z}^i)\, (f_2-\tau h_2)\, , \]

in terms of the gauge kinetic matrix \(\mathcal{N}\) which can be computed using model.gauge_kinetic_matrix. This expression can be solved for fluxes \((f_1,h_1)\) (corresponding to mode="ISD+") or \((f_2,h_2)\) (associated to mode="ISD-") corresponding to the components of the RR-flux vector \(f=(f_1,f_2)\) and the NSNS-flux vector \(h=(h_1,h_2)\) for given input values for the moduli \(z^i\) and the axio-dilaton \(\tau\). The other two modes are obtained by rewriting the above expression in the following form

\[ f=(s \, M(z^i,\overline{z}^i)\Sigma + c_0)\, h\; ,\quad \tau=c_0 + \text{i} s \, . \]

This can then be solved for RR-flux vector \(f=(f_1,f_2)\) (mode="F") or the NSNS-flux vector \(h=(h_1,h_2)\) (mode="H").

Using the same choice of fluxes and initial guesses as above, let us first collect

Summary of ISD sampling modes:

Mode	Fixed input	Solved for	ISD condition used
`"H"`	NSNS fluxes \(h = (h_1, h_2)\)	RR fluxes \(f\)	\(f = (s\,M\Sigma + c_0)\,h\)
`"F"`	RR fluxes \(f = (f_1, f_2)\)	NSNS fluxes \(h\)	\(h = (s\,M\Sigma + c_0)^{-1}\,f\)
`"ISD+"`	\((f_2, h_2)\)	\((f_1, h_1)\)	\(f_1 - \tau h_1 = \bar{\mathcal{N}}\,(f_2 - \tau h_2)\)
`"ISD-"`	\((f_1, h_1)\)	\((f_2, h_2)\)	\(f_2 - \tau h_2 = \bar{\mathcal{N}}^{-1}(f_1 - \tau h_1)\)

Here \(M\) is the ISD matrix (model.ISD_matrix), \(\mathcal{N}\) is the gauge kinetic matrix (model.gauge_kinetic_matrix), \(s = \mathrm{Im}(\tau)\), and \(c_0 = \mathrm{Re}(\tau)\). Modes "H" and "F" tend to give the best initial guesses because the ISD matrix \(M\) is computed at the exact input moduli values.

h1 = jnp.array([39., -13.,  -4.])
h2 = jnp.array([ 0.,  0.,  0.])
f1 = jnp.array([4.,  -3.,  -2.])
f2 = jnp.array([-2.,   3.,   2.])

# Choices of H-fluxes
Hflux = jnp.append(h1,h2)

# Choices of F-fluxes
Fflux = jnp.append(f1,f2)

# Choice for ISD+
ISDflux_plus = jnp.append(f2,h2)

# Choice for ISD-
ISDflux_minus = jnp.append(f1,h1)

Then the four different versions of ISD sampling can be used as follows

fluxes0 = sampler.ISD_sampling(z0,jnp.conj(z0),tau0,jnp.conj(tau0),Hflux,mode="H",output="full",return_integer_flux=True).real
fluxes0,np.sum(np.abs(model.DW(z0,jnp.conj(z0),tau0,jnp.conj(tau0),fluxes0)))

fluxes0 = sampler.ISD_sampling(z0,jnp.conj(z0),tau0,jnp.conj(tau0),Fflux,mode="F",output="full",return_integer_flux=True).real
fluxes0,np.sum(np.abs(model.DW(z0,jnp.conj(z0),tau0,jnp.conj(tau0),fluxes0)))

fluxes0 = sampler.ISD_sampling(z0,jnp.conj(z0),tau0,jnp.conj(tau0),ISDflux_plus,mode="ISD+",output="full",return_integer_flux=True).real
fluxes0,np.sum(np.abs(model.DW(z0,jnp.conj(z0),tau0,jnp.conj(tau0),fluxes0)))

fluxes0 = sampler.ISD_sampling(z0,jnp.conj(z0),tau0,jnp.conj(tau0),ISDflux_minus,mode="ISD-",output="full",return_integer_flux=True).real
fluxes0,np.sum(np.abs(model.DW(z0,jnp.conj(z0),tau0,jnp.conj(tau0),fluxes0)))

We notice that we do not always obtain the same choice of fluxes0 for each version. This is because the different versions of ISD sampling find only the closest integer flux vector satisfying the ISD condition for the given input of moduli values.

Vectorising ISD sampling#

We also provide the optional argument to vectorise ISD sampling across many input points. To do so, we set vmap=True and in_axes=(0,0,None) where the latter implies that we are vectorising over moduli and the axio-dilaton, but not over fluxes. This gives us

np.random.seed(1)
moduli,tau = sampler.initial_guesses(100,include_fluxes=False)
fluxes0 = sampler.ISD_sampling(moduli,jnp.conj(moduli),tau,jnp.conj(tau),ISDflux_plus,mode="ISD+",output="full",return_integer_flux=True,
                     in_axes=(0,0,None),vmap=True)

fluxes0[:5]

As expected, certain components of the flux vector are always constant because they are specified by ISDflux_plus=(-2., 3., 2., 0., 0., 0.).

Alternatively, we can also evaluate the same function across many different input fluxes as follows:

np.random.seed(1)
moduli,tau,ISDfluxes = sampler.initial_guesses(100,include_fluxes=True)
fluxes0 = sampler.ISD_sampling(moduli,jnp.conj(moduli),tau,jnp.conj(tau),ISDfluxes[:,:model.n_fluxes],mode="ISD+",output="full",return_integer_flux=True,
                     in_axes=(0,0,0),vmap=True)

fluxes0[:5]

Take-aways#

High-level API:

When you want to…	Call
Ready-to-go starting points	`sampler.initial_guesses(N)`
Starting points filtered for physical viability	`sampler.initial_guesses(N, filter_moduli=True)`
Starting points with ISD-refined flux	`sampler.initial_guesses_ISD(N, mode="H")`
Zoom into a per-direction sub-region for one call	`minval_moduli=[...], maxval_moduli=[...]`

Sampler configuration (at construction):

When you want to…	Pass
Uniform box in every direction	`moduli_bounds=(lo, hi)`
Per-direction box bounds	`moduli_bounds=(lower_vec, upper_vec)`
Stay away from all Kähler-cone walls	`stretching=δ`
Stay close to the conifold wall but far from others	`stretching=δ, exclude_walls=[k]`
Hard L² cap on cone-sampled points	`cone_cutoff=R`

Fine-grained samplers:

Component	Method
Imaginary part of \(z\)	`get_moduli(N, sampling_mode=…)`
Real part of \(z\)	`get_axions(N)`
Axio-dilaton \(\tau\)	`get_dilaton(N)`
Complex scalars \(z\), \(\tau\)	`get_complex_moduli(N)`, `get_complex_tau(N)`
Integer fluxes inside a box or on a sphere	`get_fluxes(N, sampling_mode="box"\|"sphere")`

Filters:

Filter	Drops points where
`filter_by_instantons`	Instanton / perturbative prepotential ratio > `inst_cutoff`
`filter_by_km`	Kähler metric is not positive-definite
`filter_moduli`	Both of the above fail
`filter_points`	Kähler-cone hyperplane(s) are violated; respects `exclude_walls`

Cone geometry:

Helper	Purpose
`find_interior_points(N, stretching)`	Fresh interior-point set (non-cached)
`update_interior_points(num_pts)`	Rebuild the cached interior-point pool
`sample_ray()`	Single random ray from the cone’s generator set
`rescale_points(pts, norm, maxval)`	Project onto an L¹/L²/L∞ ball of given radius

Next steps:

NB06: ISD sampling principle — single-vacuum ISD workflow with Newton refinement
NB07: ISD sampling workflows — large-scale moduli-first and flux-first sampling using sample_SUSY_flux_vacua
NB08: flux bounding — stochastic search over bounded flux lattice

Sampling module

Contents

Sampling module#

How to use this notebook#

Outline#

1. Imports#

2. Model and sampler setup#

2.1 Scalar `moduli_bounds` (backward-compatible)#

2.2 Direction-aware `moduli_bounds`#

2.3 Cone knobs: `stretching`, `exclude_walls`, `cone_cutoff`#

3. `initial_guesses` — the one-shot API#

3.1 Minimal example#

3.2 Keyword overrides#

3.3 Filter-aware sampling#

4. Component samplers#

4.1 `get_moduli` and its 7 sampling modes#

4.2 `get_axions` and `get_dilaton`#

4.3 `get_fluxes`#

4.4 `get_complex_moduli`, `get_complex_tau`#

5. Filter API#

5.1 `filter_by_instantons`#

5.2 `filter_by_km`#

5.3 `filter_moduli` (compound)#

5.4 `filter_points` (cone membership + `exclude_walls`)#

6. Geometry helpers#

6.1 `find_interior_points` and `update_interior_points`#

6.2 `sample_ray`#

6.3 `rescale_points`#

7. ISD sampling#

General principle of ISD sampling#

Different versions of ISD sampling#

Vectorising ISD sampling#

Take-aways#

Further reading#

Sampling module

Contents

Sampling module#

How to use this notebook#

Outline#

1. Imports#

2. Model and sampler setup#

2.1 Scalar moduli_bounds (backward-compatible)#

2.2 Direction-aware moduli_bounds#

2.3 Cone knobs: stretching, exclude_walls, cone_cutoff#

3. initial_guesses — the one-shot API#

3.1 Minimal example#

3.2 Keyword overrides#

3.3 Filter-aware sampling#

4. Component samplers#

4.1 get_moduli and its 7 sampling modes#

4.2 get_axions and get_dilaton#

4.3 get_fluxes#

4.4 get_complex_moduli, get_complex_tau#

5. Filter API#

5.1 filter_by_instantons#

5.2 filter_by_km#

5.3 filter_moduli (compound)#

5.4 filter_points (cone membership + exclude_walls)#

6. Geometry helpers#

6.1 find_interior_points and update_interior_points#

6.2 sample_ray#

6.3 rescale_points#

7. ISD sampling#

General principle of ISD sampling#

Different versions of ISD sampling#

Vectorising ISD sampling#

Take-aways#

Further reading#

2.1 Scalar `moduli_bounds` (backward-compatible)#

2.2 Direction-aware `moduli_bounds`#

2.3 Cone knobs: `stretching`, `exclude_walls`, `cone_cutoff`#

3. `initial_guesses` — the one-shot API#

4.1 `get_moduli` and its 7 sampling modes#

4.2 `get_axions` and `get_dilaton`#

4.3 `get_fluxes`#

4.4 `get_complex_moduli`, `get_complex_tau`#

5.1 `filter_by_instantons`#

5.2 `filter_by_km`#

5.3 `filter_moduli` (compound)#

5.4 `filter_points` (cone membership + `exclude_walls`)#

6.1 `find_interior_points` and `update_interior_points`#

6.2 `sample_ray`#

6.3 `rescale_points`#