Extending hi_class

Extending hi_class

How to add a new gravity model to hi_class v3, either a covariant (Lagrangian-based) theory or a parametrization of the α-functions. The short version: covariant models are specified by filling one struct of Horndeski G-functions, and everything else — background equations, α-functions, perturbations, stability tests — is derived from it automatically.

The v3 layout

All modified-gravity code lives in the gravity_smg/ directory, kept separate from the standard CLASS modules in source/ (which only contain thin hooks calling into it):

• gravity_smg/gravity_models_smg.c — everything that is model-specific: parsing gravity_model from the ini file, reading and validating parameters_smg, the G-functions of each covariant theory, the α-functions of each parametrization, expansion histories, initial conditions, and the "Modified gravity: …" banner printed at startup. This is where most of a new model goes.
• gravity_smg/gravity_functions_smg.c — model-independent machinery that turns the G-functions into everything else: gravity_functions_Es_from_Gs_smg (coefficients of the Friedmann constraint E0 + E1·H + E3·H³ = E2·H², solved algebraically for H), gravity_functions_Ps_and_Rs_from_Gs_smg (coefficients of the coupled H′ and φ″ equations), gravity_functions_building_blocks_from_Gs_smg (ρ_smg, p_smg, the effective Planck mass M² and the α-functions), and gravity_functions_As_from_alphas_smg / ..._Bs_from_Gs_smg / ..._Cs_from_Bs_smg (the coefficients entering the perturbation equations). You should not need to touch this file for a Horndeski model.
• gravity_smg/background_smg.c — the background driver: background_gravity_functions_smg initializes a G-functions struct to its GR defaults, calls the model dispatch to fill it, and chains the gravity_functions_* routines; also initial-condition plumbing (background_initial_conditions_smg), the background output columns (background_store_columntitles_smg, background_output_data_smg) and the post-integration stability_tests_smg.
• gravity_smg/perturbations_smg.c — the linear perturbation equations (perturbations_einstein_scalar_smg, perturbations_einstein_tensor_smg), the scalar-field initial conditions and their consistency tests (perturbations_tests_smg), and the quasi-static approximation machinery. Model-independent; new Horndeski models never touch it.
• gravity_smg/input_smg.c — ini-file plumbing: input_read_parameters_smg reads the common flags and delegates model parsing to gravity_models_gravity_properties_smg and gravity_models_expansion_properties_smg; input_default_params_smg holds the defaults for every _smg switch.
• gravity_smg/hmcode_smg.c — the _smg corrections used by the nonlinear (HMcode / fourier) module.
• gravity_smg/include/ — one header per file with the prototypes, plus hi_class.h and the G_functions_and_derivs struct in gravity_models_smg.h.

The central abstraction: struct G_functions_and_derivs

A covariant model in hi_class is a Horndeski (or GLPV) Lagrangian

L = G₂(φ,X) − G₃(φ,X) □φ + G₄(φ,X) R + … + G₅ and beyond-Horndeski F₄, F₅ terms,

and the code represents it as struct G_functions_and_derivs, defined in gravity_smg/include/gravity_models_smg.h. Its members are the values of G₂, G₃, G₄, G₅, F₄, F₅ and all the partial derivatives with respect to φ and the kinetic term X that the equations of motion need — G2, G2_X, G2_phi, G2_Xphi, … up to fourth derivatives for G₄/G₅ (e.g. G4_XXXphi). The macro DEFAULT_G_FUNCTIONS_AND_DERIVS in the same header initializes every entry to zero except G4 = 1/2, i.e. exactly General Relativity in the code's units (there is also DG4 = G4 − 1/2, the deviation from GR, which improves numerical precision when G₄ is close to ½).

A new Horndeski model is specified entirely by overwriting the nonzero entries of this struct as functions of φ, φ′ and the model parameters. That happens in one place: gravity_models_get_Gs_smg in gravity_models_smg.c, a chain of if (pba->gravity_model_smg == …) blocks, one per model, each computing X = (φ′/a)²/2 from the integrator variables and filling in its G's. Every quantity downstream — H, φ″, ρ_smg, the α-functions, the perturbation coefficients, the stability variables — is computed from the struct by generic code. You never write the equations of motion yourself.

Derivatives you leave untouched stay at their GR value, so a minimally coupled quintessence model sets only four entries and a conformally coupled cubic galileon about eight. Consistency matters: if you set G4, also set DG4 and every φ/X derivative of your G₄ that is nonzero.

Checklist: adding a covariant model

The complete list of places a new model touches, taken from the most recent worked example — commit 31382d65 in the development repository, which adds the conformally coupled cubic galileon (cccg_exp). That commit changes exactly six files:

gravity_models/cccg_exp.ini
gravity_models/cccg_exp_canonical.ini
gravity_smg/gravity_models_smg.c
gravity_smg/input_smg.c
hi_class.ini
include/background.h

Working through them in a sensible order:

1. Register the model enum — include/background.h. Add your model's name to enum gravity_model. If your model has variants (branches, sub-cases fixing some parameters), add them to enum gravity_submodel; the parsed value is stored in pba->gravity_submodel_smg (its default, unspecified, is set in input_default_params_smg in input_smg.c).

2. Parse the parameters — gravity_models_gravity_properties_smg in gravity_models_smg.c. This function is a chain of if (strcmp(string1,"model_name") == 0) blocks. Yours must:

• set pba->gravity_model_smg to your enum value and flag2 = _TRUE_ (the flag that marks the model name as recognized);
• set pba->field_evolution_smg = _TRUE_ — this is what tells the code to integrate φ and solve the Friedmann equation self-consistently, so covariant models need no expansion_model;
• set pba->parameters_size_smg and read the parameter vector with class_read_list_of_doubles("parameters_smg", pba->parameters_smg, pba->parameters_size_smg) (the macro also checks the user supplied the right number);
• set the shooting defaults: pba->tuning_index_smg (which entry of parameters_smg the code varies to hit the requested Ω_smg today) and pba->tuning_dxdy_guess_smg (an estimate of d(parameter)/dΩ_smg to start the shooting), guarded by has_tuning_index_smg so a user-supplied choice wins. The shooting itself lives in source/input.c and needs nothing model-specific beyond these two numbers.
• optionally handle submodels: read gravity_submodel with parser_read_string, translate the user-facing parameters of each branch into the internal parameters_smg ordering, and set pba->gravity_submodel_smg. The cccg_exp block is a full example, including a secondary shooting parameter (tuning_index_2_smg) that tunes the Planck mass today.

Also extend the class_test at the bottom of the function — the error string listing all valid gravity_model values — so a typo in another model's name produces an up-to-date message.

3. The G-functions — gravity_models_get_Gs_smg in gravity_models_smg.c. Add an else if (pba->gravity_model_smg == your_model) block that unpacks pba->parameters_smg and fills the G_functions_and_derivs struct, as described above. For cccg_exp this block is ~25 lines: G₂ = c₂X − 3V₀H₀², G3_X, and G₄ = e^{βφ}/2 with DG4, G4_phi, G4_phiphi.

4. Initial conditions — gravity_models_initial_conditions_smg in gravity_models_smg.c. A switch (pba->gravity_model_smg) over models; add a case that sets pvecback_integration[pba->index_bi_phi_smg] and pvecback_integration[pba->index_bi_phi_prime_smg] at the initial scale factor a (deep in radiation domination; rho_rad is passed in for attractor-type conditions). This is called from background_initial_conditions_smg in background_smg.c, which adds generic sanity checks — you normally only touch the gravity_models_ side. Prefer physically meaningful IC parameters: cccg_exp takes Ω_smg(z = 10¹⁰) and M²(a_ini) from the user and converts them to φ_ini, φ′_ini on the early-time attractor.

5. Startup banner — gravity_models_print_stdout_smg in gravity_models_smg.c. Add a case printing the model name and parameters when background_verbose > 0.

6. Expansion history — usually nothing. Because field_evolution_smg = _TRUE_, the expansion is computed from your Lagrangian and the expansion_model machinery (gravity_models_expansion_properties_smg, gravity_models_get_back_par_smg) is bypassed entirely. Only parametrized models use it.

7. Documentation — hi_class.ini. Add one line to the enumerated list of covariant theories under the gravity_model explanation.

8. Example ini — gravity_models/your_model.ini. Ship a runnable example: a comment block stating the Lagrangian and the meaning/ordering of parameters_smg, the gravity_model (+ gravity_submodel) and parameters_smg lines, and — crucially — any stability workarounds the model needs, with comments explaining why (e.g. cccg_exp.ini sets a_min_stability_test_smg = 1e-6 and zero scalar-field perturbation ICs because at early times its kinetic determinant is at machine-noise level). Remember Omega_Lambda = 0, Omega_fld = 0, Omega_smg = -1 (the −1 requests shooting on Ω_smg).

Adding a parametrized model

Parametrizations skip the Lagrangian: they keep pba->field_evolution_smg = _FALSE_, choose an expansion_model for the background, and specify the α-functions directly as functions of time. Two places in gravity_models_smg.c:

• Parsing, again in gravity_models_gravity_properties_smg: like propto_omega, set the model enum, field_evolution_smg = _FALSE_, M2_evolution_smg = _TRUE_ (so the running of M² is integrated), parameters_2_size_smg, and read into pba->parameters_2_smg (parametrizations use the parameters_2_smg vector; parameters_smg belongs to the expansion history).
• The α-functions, in gravity_models_get_alphas_par_smg: a chain of if (pba->gravity_model_smg == …) blocks, each filling pvecback[pba->index_bg_kineticity_smg], index_bg_braiding_smg, index_bg_tensor_excess_smg, index_bg_M2_running_smg and the Planck-mass entries at the given scale factor. The propto_omega block (each α proportional to Ω_smg(a)) is the template. Quantities like Ω_smg = ρ_smg/ρ_tot are already available in pvecback.

If you also need a new background evolution, add it to enum expansion_model in include/background.h, its parsing to gravity_models_expansion_properties_smg, and its ρ_smg, p_smg to gravity_models_get_back_par_smg. See Models: Parametrized alphas and Models: Expansion history for the user-facing side.

Testing a new model

• Watch the background. Run with background_verbose = 4 or higher: you get the shooting iterations, the initial φ, φ′, and (at high verbosity) the M²-tuning steps. gravity_models_print_stdout_smg confirms the parameters were parsed as you intended.
• Dump the diagnostics. Set write_background = yes and output_background_smg = 3: the background table then includes the α-functions, the stability variables c_s², D, M², c_t², the field variables φ, φ′, φ″ with the Friedmann-constraint residual, and the perturbation coefficients λ_i. Column meanings are in Output files.
• Respect the stability tests. After integration, stability_tests_smg (in background_smg.c) aborts on ghost (D < 0, M² < 0) or gradient (c_s² < 0, c_t² < 0) instabilities, and perturbations_tests_smg checks the scalar-field initial conditions. If a healthy model trips them through early-time rounding noise (common when Ω_smg starts at ~10⁻¹⁵), the sanctioned knobs are a_min_stability_test_smg, the *_safe_smg tolerances and the choice of pert_initial_conditions_smg — document whichever you use in the example ini. See Stability and precision.
• Compare limits. Every new model should reproduce a known one in some corner of its parameter space: couplings → 0 should give ΛCDM/quintessence, and special parameter choices often reduce to a shipped model (e.g. cccg_exp with β = 0 contains the cubic galileon). Compare C_ℓ, P(k) and H(z) against the limit model at the 10⁻⁴ level or better, and against an independent implementation if one exists.

Worked example

Commit 31382d65 in the development repository adds a complete model (cccg_exp, the conformally coupled cubic galileon of arXiv:2003.06396) in a single self-contained diff — enum, parsing with two submodels and a secondary tuning, G-functions, attractor initial conditions, documentation and two example inis — and can be read as a template for all of the above. The user-facing documentation of the shipped covariant models is in Models: Covariant theories.