Skip to content

Configuration

Aragog reads its input from a single configuration file. Two formats are accepted: TOML (.toml, recommended for new work) and INI (.cfg).

When Aragog is invoked from PROTEUS the configuration is built programmatically by the PROTEUS wrapper and the file-based path is bypassed; the same key names and value semantics apply in both cases.

TOML configuration

A TOML configuration file specifies all model parameters grouped by section. Below is a representative coupled-mantle setup; the reference tables that follow document every key.

[solver]
start_time = 0
end_time = 1.0e9                  # yr
atol = 1e-9
rtol = 1e-9
tsurf_poststep_change = 30
cvode_output_points = 65          # dense-output grid per CVODE macro-step
max_steps = 100000                # internal CVODE steps per output interval
# tcore_change_limit = 3000       # optional; flag a per-solve core-T change above this [K] (off by default)

[boundary_conditions]
outer_boundary_condition = 4      # 4 = prescribed flux (PROTEUS coupling)
outer_boundary_value = 280.0      # W/m^2
inner_boundary_condition = 1      # 1 = core cooling
inner_boundary_value = 0.0        # ignored when inner_bc=1
emissivity = 1.0
equilibrium_temperature = 273.0
core_heat_capacity = 880.0
tfac_core_avg = 1.147
param_utbl = false
param_utbl_const = 1.0e-7
core_bc = "energy_balance"        # CMB BC mode; default matches PROTEUS production

[mesh]
outer_radius = 6.371e6            # m
inner_radius = 3.480e6            # m
number_of_nodes = 200
mixing_length_profile = "nearest_boundary"
core_density = 11000.0            # kg/m^3
surface_density = 4078.95095544   # kg/m^3 (PROTEUS production default; SPIDER -adams_williamson_rhos)
gravitational_acceleration = 9.81 # m/s^2 (scalar fallback)
adiabatic_bulk_modulus = 260e9    # Pa
adams_williamson_beta = 0.0
surface_pressure = 0.0            # Pa
mass_coordinates = true           # PROTEUS production default; SPIDER mesh layout
eos_method = 1                    # 1 = Adams-Williamson, 2 = external file
eos_file = ""                     # used when eos_method = 2

[energy]
conduction = true
convection = true
gravitational_separation = true
mixing = true
radionuclides = true
tidal = false
eddy_diffusivity_thermal = 1.0
eddy_diffusivity_chemical = 1.0
kappah_floor = 10.0                # PROTEUS production default; 0.0 = textbook MLT
convective_threshold = 0.0         # 0.0 = SPIDER switch; > 0: minimum drive gradient and floor band [J/kg/K/m]
phi_step_cap = 0.0                 # 0 = disabled; 0.05 caps per-call melt-fraction change in mushy band
eos_floor_margin = 0.0             # 0 (default) = off; > 0: J/kg/K above the entropy table minimum, 100 tested
phase_boundary_cap = "rate"        # "rate" (default) or "fixed" (1 yr step near a phase boundary)
bottom_up_grav_sep = true
phase_smoothing = "tanh"           # "tanh" (default, SPIDER-parity) or "cubic_hermite"
solver_method = "cvode"            # "cvode" | "radau" | "bdf"
use_jax_jacobian = true
tidal_array = [0.0]

[initial_condition]
initial_condition = 1              # 1=linear T, 2=user file, 3=adiabat
surface_temperature = 3600.0       # K
basal_temperature = 4000.0         # K
init_file = ""                     # used when initial_condition = 2

[phase_solid]
density = "data/lookup/density_solid.dat"
heat_capacity = "data/lookup/heat_capacity_solid.dat"
melt_fraction = 0.0
thermal_conductivity = 4.0
thermal_expansivity = 1.0e-5
viscosity = 1.0e21

[phase_liquid]
density = "data/lookup/density_melt.dat"
heat_capacity = "data/lookup/heat_capacity_melt.dat"
melt_fraction = 1.0
thermal_conductivity = 4.0
thermal_expansivity = 1.0e-5
viscosity = 1.0e2

[phase_mixed]
latent_heat_of_fusion = 4.0e5
rheological_transition_melt_fraction = 0.4
rheological_transition_width = 0.15
solidus = "data/lookup/solidus.dat"
liquidus = "data/lookup/liquidus.dat"
phase = "composite"
phase_transition_width = 0.01
grain_size = 1.0e-3
matprop_smooth_width = 0.0
separation_viscosity = "melt"      # "melt" or "mixture"
cp_blend = "latent"                # "latent" or "linear"
const_properties = false

[radionuclide_K40]
name = "K40"
t0_years = 4.55e9
abundance = 1.1668e-4
concentration = 310.0
heat_production = 2.8761e-5
half_life_years = 1.248e9

Any number of [radionuclide_*] sections can be added; each section name must start with radionuclide_. Each isotope heats at heat_production * abundance * 1e-6 * concentration * 2**((t0_years - t) / half_life_years) per kg, with concentration in ppm; an isotope whose amplitude is zero contributes no heating at any time. With radionuclides = true, loading raises a ValueError when a half-life is not positive, when heat_production, abundance or concentration is negative, or when an isotope's heating is not finite at start_time.

Configuration sections

A configuration file containing a [scalings] block is rejected at load time with a ValueError. Aragog applies its internal nondimensionalisation around the integrator only (entropy reference \(S_\mathrm{ref}\) and time reference \(t_\mathrm{ref}\)); it is not user-configurable.

[solver]

Time-integration controls.

Key Unit Description
start_time yr Start of integration
end_time yr End of integration
atol -- Absolute tolerance, default \(10^{-8}\) (floored at \(10^{-8}\))
rtol -- Relative tolerance, default \(10^{-8}\), the value at which phase_boundary_cap = "rate" is verified. Above \(10^{-8}\) a CVODE run logs one warning: a loose rtol can lock CVODE at the convective switch (the convective mask and the \(\kappa_h\) floor). The counts below were measured at the default convective_threshold = 0, where the switch is at \(dS/dr = 0\); the probe did not cover \(D > 0\), where CVODE also stalled at 1e-8 (see Heat transport). They come from a probe variant of the energy_balance isentropic start of tests/test_phase_boundary_cap.py (24 nodes, 200 yr, fixed cap unless stated) that draws the initial entropy uniformly between 8000 and 8400 J/kg/K for each run and adds a relative noise of \(10^{-15}\); the probe script is not in the repository, and the runs were on macOS arm64. A run counts as locked when it has not reached 200 yr after 600 s of wall time: every run was screened at 60 s, and each run that missed that was run again alone to 600 s, by which time a locked run had used about \(1.7 \times 10^5\) RHS evaluations. With rtol = atol = 1e-6, CVODE locked in 12 of 1133 runs, in 4 of 474 on the quasi_steady core and in 2 of 286 with the rate cap, where 2 more runs finished only after 233 s and 346 s. With both at 1e-5, all 80 runs were run to 600 s and none locked (too few runs to show it safe; the slowest finished after 100.6 s and \(6.5 \times 10^4\) RHS evaluations), but the median run needed about \(1.6 \times 10^4\) RHS evaluations, against about \(1.3 \times 10^3\) at 1e-8. With both at 1e-8, CVODE locked in none of 1026, 627 and 1224 runs for three setups (fixed and rate cap on the energy_balance core, fixed cap on quasi_steady). On CI the AVX-512 numpy path made the test's own case, at its single fixed initial entropy, stall at 1e-6, judged by the CI test timeout rather than by the 600 s criterion.
tsurf_poststep_change K Maximum allowed surface-temperature change per coupling step (PROTEUS use)
cvode_output_points -- Number of points on the CVODE dense-output grid returned per macro-step (default 65, minimum 2). The per-call energy integrals take CVODE's accepted internal steps and these output points as nodes, so a coarse grid does not under-resolve them, and their cost follows the number of accepted steps, not this grid. The core-temperature change check (tcore_change_limit) reads only the output points. The grid feeds back into CVODE stepping, so the accepted step count and the final state shift weakly with it. Used only when solver_method = "cvode".
max_steps -- Maximum number of internal CVODE steps between two consecutive output times (SUNDIALS mxstep; default 100000, minimum 1). CVODE returns CV_TOO_MUCH_WORK and stops once one output interval reaches this count. With phase_boundary_cap = "rate" each CVODE segment applies the same per-interval limit. A stiff phase-change window can need more internal steps than the default budget; raise this value to let such a solve complete. Used only when solver_method = "cvode".
tcore_change_limit K Optional limit on the per-solve core-temperature change (unset by default; must be positive when set). Aragog always measures the largest change of the core temperature from the solve-entry value over the returned grid and reports it as tcore_change_max. When this limit is set and that change exceeds it, Aragog sets the tcore_change_exceeded flag on the result; the flag is also set whenever any sampled core temperature is non-finite, independent of whether a limit is set, since that signals a corrupted solve. The solve never raises. A caller can use the flag to reject a solve whose core temperature jumps as it crosses a phase boundary in a single accepted step, or whose core temperature is non-finite.

[boundary_conditions]

Thermal boundary conditions at the surface and CMB.

Key Unit Description
outer_boundary_condition int Surface BC mode. 1 = grey-body atmosphere, 4 = prescribed flux, 5 = prescribed temperature
outer_boundary_value W/m² or K Surface flux (modes 1, 4) or temperature (mode 5)
inner_boundary_condition int CMB BC mode. 1 = core cooling, 2 = prescribed flux, 3 = prescribed temperature
inner_boundary_value W/m² or K CMB flux (modes 1, 2) or temperature (mode 3)
emissivity -- Surface emissivity (used in mode 1)
equilibrium_temperature K Radiative equilibrium temperature (mode 1)
core_heat_capacity J/kg/K Core specific heat capacity
tfac_core_avg -- Core adiabat correction factor (default 1.147; Bower+2018 Table 2)
param_utbl bool Enable upper-thermal-boundary-layer parameterisation (default false)
param_utbl_const K⁻² UTBL constant in \(\Delta T = b\,T^3\)
core_bc str CMB BC mode. energy_balance (SPIDER bit-parity, length \(N+1\), default and PROTEUS production), quasi_steady (alpha-factor, state length \(N\)), gradient (entropy gradient as state, length \(N+2\), experimental), bower2018 (parity-only, not recommended), or core_module (the core evolution module, length \(N+2\) plus the shell cells with stratification = true, needs inner_boundary_condition = 1). Any other value is refused.
top_node str Entropy of the top basic node. linear (default) extrapolates the top 2 cells; limited reduces that half step smoothly near an edge of the EOS entropy table, so the node rises with the top cell and, with the top cell in the table, stays in the table (aragog.mesh.top_node_entropy). limited is refused on core_bc = "gradient": when the config is read, on a later assignment of either key on a BoundaryConfig, and when the mesh is built, which also covers a parsed parameter object changed after loading. Any other value is refused.

[mesh]

Spatial discretisation and pressure-density profile.

Key Unit Description
outer_radius m Planet (mantle top) radius
inner_radius m CMB radius
number_of_nodes -- Number of basic-grid nodes (cell faces)
mixing_length_profile str "nearest_boundary" (distance to nearer mesh boundary) or "constant" (a fixed fraction of mantle thickness; see mixing_length_constant_fraction below)
mixing_length_constant_fraction -- Fraction of mantle thickness used as the mixing length when mixing_length_profile = "constant". Ignored otherwise. Default 0.25
core_density kg/m³ Mean core density
eos_method int 1 = analytic Adams–Williamson; 2 = external file (eos_file)
surface_density kg/m³ Surface mantle density (Adams-Williamson). Default 4078.95095544 (PROTEUS production; SPIDER -adams_williamson_rhos)
gravitational_acceleration m/s² Scalar gravity for Adams–Williamson; superseded by the per-node profile from eos_file when eos_method = 2. Default 9.81
adiabatic_bulk_modulus Pa \(K_S\) in Adams–Williamson. Default 260e9
adams_williamson_beta -- A–W exponent \(\beta\); 0.0 derives it from \(K_S\). Default 0.0
surface_pressure Pa Atmospheric overburden added to the pressure integration. Default 0.0
mass_coordinates bool If true, use a SPIDER-parity mass-coordinate grid with Newton-solved spatial radii. Default true (PROTEUS production); set false for uniform spacing in radius
eos_file str Path to a four-column file (r [m], P [Pa], rho [kg/m³], g [m/s²]) used when eos_method = 2. PROTEUS supplies a Zalmoxis-generated file via this key. A relative path resolves against the directory of the config file first, then the working directory.

[energy]

Heat-transport switches, transport parameters, and integrator selection.

Key Type Default Description
conduction bool -- Conductive flux \(-k\,[(T/c_p)\partial S/\partial r + (\partial T/\partial P)_S\,\partial P/\partial r]\)
convection bool -- Convective flux \(\rho T\kappa_h(-\partial S/\partial r)\), driven when \(\partial S/\partial r < 0\), with the floor up to convective_threshold
gravitational_separation bool -- Stokes-regime separation of melt and solid; phase-boundary smoothed
mixing bool -- Chemical mixing flux (SPIDER bracket form involving \(\partial S_\mathrm{liq}/\partial P\) and \(\partial S_\mathrm{sol}/\partial P\))
radionuclides bool -- Internal radiogenic heating from any [radionuclide_*] sections
tidal bool -- Tidal heating from tidal_array
eddy_diffusivity_thermal float 1.0 Scalar multiplier on \(\kappa_h\). Negative values pin \(\kappa_h\) to the absolute value (SPIDER convention)
eddy_diffusivity_chemical float 1.0 Scalar multiplier on \(\kappa_c\). Negative values pin to absolute
kappah_floor m²/s 10.0 Phase-modulated lower bound on \(\kappa_h\). Default 10.0 (PROTEUS production); set 0.0 for textbook MLT
convective_threshold J kg⁻¹ K⁻¹ m⁻¹ 0 Entropy gradient \(D\): the mixing-length drive acts only where \(\partial S/\partial r < 0\) and uses at least \(D\) there, and the kappah_floor applies up to \(\partial S/\partial r = D\), so nodes with \(0 \le \partial S/\partial r < D\) get the floor alone. Finite and \(\ge 0\); the default 0.0 gives the \(\partial S/\partial r < 0\) switch of SPIDER; \(D = 10^{-6}\) gives one end state in 3 replayed coupled calls that end in one of 2 states under that switch, but \(D > 0\) can stall CVODE at a node whose gradient stays near 0, also at rtol = \(10^{-8}\). See Heat transport
phi_step_cap -- 0.0 Per-call \(\Delta\Phi_\mathrm{global}\) cap. When > 0 and the mantle straddles the rheological transition, a SUNDIALS root function fires at the step where the mass-weighted global melt fraction \(\Phi_\mathrm{global}\) has changed by cap from its value at solve() entry. 0.05 is a useful upper bound for 1 M\(_\oplus\) runs. Default 0.0 (disabled).
eos_floor_margin J kg⁻¹ K⁻¹ 0 EOS floor margin \(m\); the default 0 leaves the floor off, and a positive margin turns it on (100 is the tested value). With the floor off, a drain past the table minimum continues at the clamped table-edge temperature and the EOS lookup logs an out-of-range WARNING (on the 1st, 10th, 100th, ... occurrence of each lookup); it raises instead only for an EOS that a caller builds with strict_range=True. With the floor on, a root function (CVODE) or a terminal event (scipy) ends the call when the smaller of the top basic-node entropy, which sets the surface temperature, and the smallest cell entropy reaches the table minimum plus \(m\), or plus \(m/2\) when the call starts below the minimum plus \(m\); the call returns a shorter step. A call that starts at or below the minimum plus \(m/2\) arms a level at half its distance to the minimum and does not advance (status 1) when its entropy falls fast enough to reach that level before the call ends, and a caller that retries with a shorter span, such as the PROTEUS retry ladder, can then advance. A start at the table minimum (within \(2 \times 10^{-6} m\)) does not advance on any fall at any span, and a start at or below the minimum, also of the extrapolated top node alone, or so close to it that no level fits between them, does not advance at all. A call whose stop lands below the minimum through the root location error, which grows with the call time, is discarded as a refusal (see the solver tuning page). A missing key keeps the default 0; 0 or a negative value disables it (None too, outside TOML, which has no None), and a margin as large as the entropy span of the tables is refused. PROTEUS has no key for it, so PROTEUS runs have the floor off. See Solver method tuning
phase_boundary_cap str "rate" Step-size cap near or inside the two-phase band. "fixed" uses 1 yr. "rate", the default (also when unset), uses event-driven CVODE segments with max_step set to 0.1 of the shortest time to reach a phase boundary, clipped to [1, 100] yr. Inside the stiff zone (\(\delta = \max(3w \max_j(S_{\mathrm{liq},j}-S_{\mathrm{sol},j}), 10\ \mathrm{J\,kg^{-1}\,K^{-1}})\), with \(w = \mathtt{matprop\_smooth\_width}\)) or below the rate floor, the estimate uses the distance to the nearer boundary. The gradient core and the scipy integrators use 1 yr in both modes. On CVODE the segments of a "rate" call arm within the larger of phase_boundary_entropy_margin and \(\delta\); the step caps, "fixed", the gradient core and the scipy integrators use phase_boundary_entropy_margin. At the default rtol \(= 10^{-8}\), "rate" stays within \(1.3 \times 10^{-4}\) K of a "fixed" run at \(10^{-10}\) in the cases on the energy equation page; an rtol above \(10^{-7}\) makes the rate cap log one warning.
bottom_up_grav_sep bool true Apply SPIDER's bottom-up gating on the gravitational-separation flux
phase_smoothing str "tanh" Phase-boundary smoothing. "tanh" (default) is SPIDER's two-branch get_smoothing with width matprop_smooth_width = 0.01; "cubic_hermite" is the fallback \(16\,g\phi^2(1-g\phi)^2\) form.
solver_method str "cvode" ODE integrator. "cvode" selects SUNDIALS CVODE via scikits.odes (default); "radau" and "bdf" use scipy solve_ivp. When scikits.odes is not installed the solver falls back to Radau with a warning.
use_jax_jacobian bool true When solver_method = "cvode", install a JAX-traced analytic Jacobian instead of the finite-difference Jacobian. Requires JAX and the PROTEUS-side factory registration; falls back to FD Jacobian otherwise.
tidal_array array[float] [0.0] Per-staggered-node tidal heating in W/kg. Length must be 1 (broadcast scalar) or number_of_nodes - 1

[initial_condition]

Initial entropy profile (mapped from the temperature inputs via the EOS).

Key Unit Description
initial_condition int 1 = linear T(r) profile (default), 2 = read from init_file, 3 = adiabatic
surface_temperature K Top temperature for modes 1 and 3
basal_temperature K Bottom temperature for mode 1
init_file str Path to a two-column r, T (or r, S) file used when initial_condition = 2. A relative path resolves against the directory of the config file first, then the working directory.

When the PROTEUS wrapper drives Aragog the initial entropy is supplied directly via EntropySolver.set_initial_entropy(); the section keys are still parsed but the values are overridden by the call.

[phase_liquid] and [phase_solid]

End-member phase properties. Float values mean a constant; string values are file paths to one-column lookup tables.

Key Unit Description
density kg/m³ Reference density (or path to lookup)
viscosity PaĀ·s Dynamic viscosity (or path)
heat_capacity J/kg/K Specific heat capacity
melt_fraction -- Reference melt fraction (0 for solid, 1 for liquid)
thermal_conductivity W/m/K Thermal conductivity
thermal_expansivity 1/K Thermal expansivity
entropy J/kg/K Optional reference entropy

In the production PROTEUS path the per-phase keys are not consumed for material properties; the EOS tables provide \(\rho\), \(c_p\), \(\alpha\), \(k\), and \(T\) as functions of \((P, S)\). The values are kept for the standalone constant-properties path (see const_properties below).

[phase_mixed]

Mixed-phase (mushy) parameters and the optional constant-properties analytical mode.

Key Unit Default Description
latent_heat_of_fusion J/kg -- Reference latent heat (legacy; production runs use \(L(P)\) from the EOS)
rheological_transition_melt_fraction -- -- Critical melt fraction for the viscosity transition
rheological_transition_width -- -- Width of the tanh viscosity blend
solidus str -- Path to a solidus reference (legacy T-based; not consumed by the entropy solver)
liquidus str -- Path to a liquidus reference (legacy T-based)
phase str -- Phase label
phase_transition_width -- -- Legacy smoothing width
grain_size m 1.0e-3 Grain size \(a\) in the melt-solid mobility \(F(\zeta)=K(\zeta)/\zeta\)
matprop_smooth_width -- 0.0 SPIDER-parity blending width across phase boundaries
separation_viscosity str "melt" Drag viscosity in \(v_\mathrm{rel} = \lvert\Delta\rho\rvert g F(\phi)/\eta\): "melt" (fixed single-phase liquid viscosity, SPIDER parity) keeps separation active below phi_rheo, which collapses the CMB temperature in coupled caps-off runs; "mixture" (rheological-transition-blended bulk viscosity) locks separation up below phi_rheo and avoids that collapse, which is why PROTEUS defaults to "mixture"
cp_blend str "latent" "latent" (SPIDER parity) or "linear" for the mushy \(C_p\) rule
const_properties bool false Enable the analytic constant-properties path (no EOS tables needed)
const_rho kg/m³ 4000 Constant density when const_properties = true
const_Cp J/kg/K 1000 Constant heat capacity
const_alpha 1/K 1e-5 Constant thermal expansivity
const_cond W/m/K 4.0 Constant thermal conductivity
const_log10visc -- 2.0 Constant log₁₀ viscosity
const_T_ref K 3500 Reference temperature for \(T(S) = T_\mathrm{ref}\exp((S - S_\mathrm{ref})/C_p)\)
const_S_ref J/kg/K 3000 Reference entropy in the same expression

Loading configuration in Python

from aragog.parser import Parameters

# From a TOML or INI file (auto-detected by suffix)
params = Parameters.from_file("path/to/config.toml")

# From a dictionary (used by PROTEUS via Config.from_dict)
from aragog.config import Config
params = Config.from_dict(config_dict)

Both calls return a Parameters dataclass; the Config facade exists primarily for the dict-driven path used by the PROTEUS wrapper.

Parameters.from_file dispatches on file suffix:

Suffix Loader Semantics
.toml stdlib tomllib TOML: quoted strings ("nearest_boundary") and # ... inline comments. Strings round-trip without surrounding quotes; comments are stripped before type coercion.
.cfg, .ini, others typed_configparser INI: bare unquoted strings (nearest_boundary). No inline comment support; values are read verbatim, quotes included.

Use the format matching the file you are writing: a .toml file must use TOML conventions, an .ini / .cfg file must use INI conventions. Mixing the two (e.g. quoted strings in a .cfg) is a common cause of Mixing length profile is unknown errors.

Legacy INI format

The legacy .cfg format is still parsed. Example files live under src/aragog/cfg/. Keys and section names match the TOML schema; only the surface syntax differs (no quotes around strings, no # inline comments).