Skip to content

Code architecture

CALLIOPE is a small package: seven Python modules totaling under 1000 raw lines, no compiled extensions, no thread or process pools. This page maps each module to its responsibility and shows how they call each other.

Package layout

src/calliope/
├── __init__.py              # version string only
├── constants.py             # molar masses, element list, plotting colours
├── oxygen_fugacity.py       # IW buffer parameterisations
├── chemistry.py             # ModifiedKeq class, six redox couples
├── solubility.py            # five solubility classes
├── solve.py                 # speciation, mass balance, outer solver
└── structure.py             # simple core/mantle mass split

Plus the test suite:

tests/
├── test_init.py             # version smoke test
├── test_core.py             # per-module unit tests
└── test_stoichiometry.py    # cross-module mass-balance tests

And the development tooling under tools/ (Jupyter notebooks for fitting JANAF tables to Schaefer-style functional forms; the JANAF source data files for H\(_2\)S, NH\(_3\), S\(_2\)).

Module responsibilities

constants.py

Exports the static data CALLIOPE needs:

  • M_earth, R_earth, R_core_earth, M_core_earth, ocean_moles: planetary reference values for unit conversions.
  • const_G, mol, R_gas: physical constants.
  • molar_mass: dictionary mapping species names to kg/mol.
  • volatile_species: ordered list of the eleven species CALLIOPE handles.
  • element_list: ordered list of the five elements (H, O, C, N, S).
  • dict_colors: colour scheme used in PROTEUS plots (so figures across the ecosystem stay consistent).

Pure data; no side effects, no logic.

oxygen_fugacity.py

Single class OxygenFugacity with two model methods (oneill default, fischer). Stateless: instantiate once, call repeatedly with (T, fO2_shift). See Oxygen fugacity for the equations.

chemistry.py

Single class ModifiedKeq with eight model methods (schaefer_C, schaefer_H, schaefer_CH4, janaf_CO, janaf_H2, janaf_SO2, janaf_H2S, janaf_NH3). Composes with OxygenFugacity to produce the modified equilibrium constant \(G_\mathrm{eq}(T, f_{\mathrm{O}_2})\). See Equilibrium chemistry.

solubility.py

Six concrete subclasses of an abstract Solubility base: SolubilityH2O, SolubilityCO2, SolubilityCO, SolubilityCH4, SolubilityN2, SolubilityS2. Each carries multiple alternative composition methods selectable through the constructor. The __call__ interface is uniform (solubility(p, *args) returns ppmw) but the *args differ per species: H\(_2\)O takes only \(p\), CO\(_2\) takes \((p, T)\), S\(_2\) takes \((p, T, \Delta\mathrm{IW})\), etc. This irregularity is by design; flat-out wrapping into a uniform signature would hide the actual physical dependences.

See Solubility laws for the equations.

solve.py

The orchestration layer. Five public functions:

  • get_partial_pressures(pin, ddict): walks the eleven-species speciation tree from the four primary pressures.
  • atmosphere_mass(pin, ddict): applies Bower et al. (2019) Eq. 2 to every species and aggregates atomic-mass tallies per element.
  • dissolved_mass(pin, ddict): applies the chosen solubility law for each soluble species and aggregates atomic-mass tallies per element.
  • get_target_from_params(ddict): translates hydrogen_earth_oceans, CH_ratio, nitrogen_ppmw, sulfur_ppmw into kg-per-element targets.
  • get_target_from_pressures(ddict): back-computes kg-per-element targets from prescribed initial atmospheric pressures.

Plus the inner-loop residual functions func (returns the four-vector residual) and obj (returns the L2 norm), and the outer driver equilibrium_atmosphere(target, ddict, ...) that ties it all together.

get_partial_pressures is the call that defines what species CALLIOPE knows about. To add a new species, you add an entry in volatile_species (constants.py), an entry in molar_mass (constants.py), a ModifiedKeq model method (chemistry.py) if it is a derived species, a Solubility subclass (solubility.py) if it has dissolved-melt physics, and the speciation step in get_partial_pressures plus the atomic-mass tally in atmosphere_mass. Tests need the analytical-vs-code consistency check in tests/test_stoichiometry.py.

structure.py

Single function calculate_mantle_mass(radius, mass, core_frac) that returns M_mantle = M_planet - M_core assuming Earth-like core density. This is a fallback for use cases where a full structure solve (Zalmoxis, SPIDER) is not warranted. The PROTEUS coupling does not go through this function; it gets M_mantle directly from Zalmoxis or SPIDER via hf_row.

Call graph

                                    equilibrium_atmosphere
                                              │
                  ┌───────────────────────────┼───────────────────────────┐
                  ▼                           ▼                           ▼
            get_initial_pressures        scipy.optimize.fsolve     scipy.optimize.minimize
                                              │                           │
                                              └─────────┬─────────────────┘
                                                        ▼
                                                       func    ◄─── ModifiedKeq.__call__
                                                        │              ▲
                  ┌─────────────────────────────────────┤              │
                  ▼                                     ▼              │
            atmosphere_mass                       dissolved_mass       │
                  │                                     │              │
                  │      ┌──────────────────────────────┤              │
                  │      ▼                              ▼              │
                  │      get_partial_pressures   Solubility{H2O,CO2,…}.__call__
                  │              │                                     │
                  ▼              ▼                                     │
            atmosphere_mean_molar_mass     ModifiedKeq.__call__  ◄─────┘
                                                  │
                                                  ▼
                                            OxygenFugacity.__call__

The graph is a DAG with a single feedback loop (the outer fsolve → func → ... → ModifiedKeq cycle). No module imports any other module's private state; all couplings are through the pin/ddict dictionaries.

What is not in CALLIOPE

By design:

  • No interior structure (structure.py is a one-line fallback, not a real model).
  • No time stepping (CALLIOPE is a steady-state equilibrium solver).
  • No I/O (no NetCDF writes, no file reads, no logging beyond the module-level log handle).
  • No threading (deterministic single-threaded; PROTEUS provides parallelism by sweeping at the outer level).
  • No JAX, no PyTorch (pure NumPy + SciPy).
  • No state (every call is functional; no class-level caches).

These omissions are deliberate: CALLIOPE is a calculation kernel, not a framework. Anything that needs framework-level features (logging, status files, retries, parallel sweeps, plotting orchestration) lives in the PROTEUS wrapper.

Performance

A single equilibrium_atmosphere() call takes:

  • ~10 ms with a good warm start (p_guess from previous iteration);
  • ~100 ms-1 s with a cold start (Monte-Carlo guess phase);
  • the residual evaluation (func) is \(\le\)1 ms; the bulk of the time is in fsolve step iteration.

For batch use cases (sensitivity sweeps, parameter studies), wrap a Python loop around equilibrium_atmosphere(); the cost of the function call is small enough that running thousands of cases serially is straightforward. For Monte-Carlo studies with \(10^5\)+ cases, parallelise at the outer Python level (e.g. multiprocessing.Pool.map); CALLIOPE itself is not threaded but is process-safe.

See also