This repository contains the Python-based numerical framework developed to solve the equations of hydrostatic equilibrium for multi-phase planetary interiors. It is explicitly designed to rigorously map the structural degeneracies between compact solid cores, volatile water mantles, and deeply suspended dilute "fuzzy" composition gradients.
This framework is built for robustness, explicitly handling severe thermodynamic discontinuities (such as liquid-vapor phase transitions) and adaptive atmospheric downsampling, making it capable of modeling both dense water-worlds (e.g., GJ 1214 b) and highly inflated super-puffs (e.g., Kepler-11e) within a single unified tool.
Primary Reference: If you use this code in your research, please cite our corresponding Astronomy & Astrophysics paper:
Wilkinson, C., Mazevet, S., Lagrange, A.-M., Charnay, B. (2026). A robust numerical framework for giant planet interior modeling: Constraining the core-envelope equivalence in the presence of dilute gradients. A&A.
fuzzycore/
├── data/ # Equation of State (EOS) tables (ANEOS, AQUA, Chabrier)
├── figures/ # Output directory for generated plots
├── notebooks/ # Jupyter notebooks for exploratory analysis
├── scripts/ # Executable scripts for parameter sweeps and modeling
├── setup.py # Package installation configuration
├── requirements.txt # Python dependencies
└── src/
└── fuzzycore/ # Core framework modules (physics.py, solver.py, eos.py, etc.)
This framework relies on standard scientific Python libraries and is designed to be installed as a global package within your environment.
# 1. Clone the repository
git clone [https://github.com/ChristianSWilkinson/fuzzycore.git](https://github.com/ChristianSWilkinson/fuzzycore.git)
cd fuzzycore
# 2. Create and activate a new virtual environment
conda create -n fuzzy_env python=3.10
conda activate fuzzy_env
# 3. Install dependencies and the package in editable mode
pip install -r requirements.txt
pip install -e .(Installing with -e . allows you to import fuzzycore from anywhere on your machine while instantly reflecting any changes you make to the source code).
The heart of fuzzycore is the solve_structure() function. By simply altering the params dictionary, the solver will dynamically switch its internal physics architecture to build whatever class of planet you request.
Every planet requires these base parameters:
M_core/M_rock: Mass of the solid interior (in kg).M_water: Mass of the water mantle (in kg). Set to0.0for a gas giant.iron_fraction: Mass fraction of iron in the core (e.g.,0.33for Earth-like).P_surf: Atmospheric pressure boundary (in bar).T_surf: Atmospheric temperature boundary (in K).z_base: The baseline metallicity of the upper atmosphere.z_profile: An array of discrete composition steps for the integrator to use.sigma_val: The width of the compositional gradient (0.01= Sharp,0.30= Fuzzy).
To build a standard sub-Neptune with a distinct, sharp boundary between a solid rock core and a fully convective gaseous envelope, set M_water = 0.0 and sigma_val = 0.01.
import numpy as np
import fuzzycore.constants as c
import fuzzycore.solver as solver
import fuzzycore.utils as utils
target_mass = 8.0 * c.M_EARTH
params = {
'M_core': 7.0 * c.M_EARTH, # 7 Earth-mass solid core
'M_water': 0.0, # No water layer
'iron_fraction': 0.33,
'P_surf': 1.0, # 1 bar surface
'T_surf': 500.0,
'z_base': 0.05, # 5% atmospheric metallicity
'z_profile': np.linspace(0.05, 1.0, 20),
'sigma_val': 0.01, # Sharp core-envelope boundary!
'debug': True
}
# Run the integrator
result = solver.solve_structure(
target_val=target_mass,
params=params,
mode='mass',
trial_id="standard_gas_giant",
csv_file="temp.csv",
write_lock=utils.DummyLock()
)To build an ocean world or a planet with a massive supercritical water mantle beneath a thin atmosphere, simply provide a non-zero M_water. The solver will automatically utilize the ab-initio water tables and nest the boundary-shooting algorithms.
params = {
'M_rock': 7.0 * c.M_EARTH, # Note: use 'M_rock' instead of 'M_core'
'M_water': 1.0 * c.M_EARTH, # 1 Earth-mass water mantle!
'iron_fraction': 0.33,
'P_surf': 1.0,
'T_surf': 500.0,
'z_base': 0.05,
'z_profile': np.linspace(0.05, 1.0, 20),
'sigma_val': 0.01,
}To simulate the massive radius inflation caused by heavy elements lofted into the gaseous envelope, shrink the solid M_core and increase sigma_val. The discrete entropy jumps will act as a deep thermal blanket, puffing out the planet.
params = {
'M_core': 2.0 * c.M_EARTH, # Tiny solid core
'M_water': 0.0,
'iron_fraction': 0.33,
'P_surf': 10.0, # Evaluate deeper at the convective boundary
'T_surf': 1200.0, # High internal heat
'z_base': 0.05,
'z_profile': np.linspace(0.05, 1.0, 20),
'sigma_val': 0.30, # Wide compositional gradient!
}We have provided three ready-to-run scripts in the scripts/ directory that reproduce the key analyses from our paper. Note: run these from the root repository directory.
This script runs a parallelized 1D parameter sweep, systematically shifting mass from a solid rock core into a suspended dilute envelope while maintaining a constant total planetary mass.
python scripts/run_kepler11e_parallel.pyThis script evaluates four distinct configurations (Sharp vs. Fuzzy cores under Hot vs. Cold internal states). It demonstrates the non-linear "amplification" effect of dilute cores.
python scripts/run_tint_amplification.pyA heavy-duty script utilizing Python's multiprocessing pool to map the deep 3x3 structural degeneracy grid across variations in total planetary mass and surface temperature.
python scripts/run_planet_grid.pyOnce a model successfully converges, it returns a dictionary of 1D arrays (R, M, P, T, Rho, S, Z). You can immediately visualize the internal architecture using our built-in PEP 8 compliant plotting suite:
import fuzzycore.plotting as plotting
# Plot internal density, temperature, and entropy profiles
plotting.plot_diagnostics(result, save_name="my_planet_structure")
# Plot the planet's trajectory across the thermodynamic phase space
plotting.plot_trajectory_on_eos(result, params, save_name="my_planet_eos")- Nested Boundary-Shooting: Decouples highly compressible atmospheres from incompressible cores to ensure absolute mass conservation across arbitrary phase boundaries.
- KD-Tree Adiabat Stepper: Uses multidimensional local linear regressions over 3D (P, T, S) phase spaces to integrate continuous fluid adiabats.
- Isothermal Phase Guards: Enforces strict thermodynamic monotonicity to bridge jagged low-pressure vapor-liquid phase transitions (e.g., 1–1000 bar) in ab-initio water tables.