How to: Geant4-DNA chemistry

This chapter explains how to run a Geant4-DNA chemistry simulation in GATE 10, how to configure chemistry lists and track-structure EM physics, and how to use chemistry actors.

Warning

The chemistry functionality is currently in alpha stage. Interfaces, behavior, and outputs may still change, so results should be validated carefully before production use.

The current chemistry support is built around the Geant4-DNA workflow:

• activate track-structure EM in the region where DNA physics should run;

• select a Geant4 chemistry list;

• optionally customize that chemistry list;

• optionally define a chemistry world and scavenger reactions;

• add one or more chemistry-aware actors;

• read summary and counter results from actor outputs.

Minimal setup

A minimal chemistry simulation needs three ingredients:

1. a standard Geant4 physics list;

2. a Geant4-DNA chemistry list;

3. track-structure EM activated in the region of interest.

The smallest practical setup looks like this:

import opengate as gate

sim = gate.Simulation()

sim.physics_manager.physics_list_name = "G4EmStandardPhysics"
sim.chemistry_manager.chemistry_list_name = "G4EmDNAChemistry_option3"
sim.chemistry_manager.time_step_model = "IRT"

water_box = sim.add_volume("Box", "water_box")
water_box.material = "G4_WATER"
water_box.size = [10 * gate.g4_units.um] * 3

water_box.set_track_structure_em_physics("G4EmDNAPhysics_option2")

This does not yet add chemistry scoring, but it prepares the simulation so that Geant4-DNA chemistry can run.

Track-structure EM physics

Chemistry requires track-structure EM physics in the relevant region. In GATE 10 this is configured with track_structure_em_physics.

The recommended full Geant4 names are for example:

• G4EmDNAPhysics

• G4EmDNAPhysics_option2

• G4EmDNAPhysics_option4

• G4EmDNAPhysics_option6

The simplest way is to configure it directly on a volume:

target = sim.add_volume("Box", "target")
target.material = "G4_WATER"
target.size = [10 * gate.g4_units.um] * 3
target.set_track_structure_em_physics("G4EmDNAPhysics_option2")

This creates or reuses the appropriate GATE region and activates Geant4-DNA track-structure EM there.

You can also configure it explicitly through a region:

region = sim.physics_manager.add_region("dna_region")
region.associate_volume(target)
sim.physics_manager.set_track_structure_em_physics_in_region(
    "dna_region",
    "G4EmDNAPhysics_option4",
)

If you use a chemistry actor, the actor may also request track_structure_em_physics for its attached volume. This is a convenience feature, but configuring the region or volume explicitly is often clearer.

Choosing a chemistry list

The chemistry list is selected through ChemistryManager:

sim.chemistry_manager.chemistry_list_name = "G4EmDNAChemistry_option3"

At the moment, one simulation must resolve to one coherent chemistry list. This means:

• the manager can request a chemistry list;

• a chemistry actor can also request a chemistry list;

• all such requests must agree.

If multiple incompatible chemistry-list names are requested, GATE stops with a fatal error.

Current built-in chemistry list names include:

• G4EmDNAChemistry

• G4EmDNAChemistry_option1

• G4EmDNAChemistry_option2

• G4EmDNAChemistry_option3

Chemistry time-step model

The chemistry time-step model is configured through:

sim.chemistry_manager.time_step_model = "IRT"

Current allowed values are:

• SBS

• IRT

• IRT_syn

For many aqueous-radiolysis examples, IRT is a reasonable starting point.

Customizing the chemistry list

The active chemistry list is available as:

chem_list = sim.chemistry_manager.chemistry_list

You can extend the selected built-in Geant4 chemistry list with:

• additional chemical species;

• additional bimolecular reactions;

• additional dissociation channels.

Add a reaction

This is the most common customization pattern:

dm3_per_mole_s = (
    1e-3 * gate.g4_units.m3 / (gate.g4_units.mole * gate.g4_units.s)
)

chem_list.add_reaction(
    reactant_a="H",
    reactant_b="H",
    rate_constant=0.503e10 * dm3_per_mole_s,
    products=["H2"],
    reaction_type=0,
)

Reactions are appended to the Geant4-DNA reaction table after the selected built-in chemistry list has created its default content.

Add a species

You can add or modify a species by name:

chem_list.add_chemical_species(
    name="MySpecies",
    charge=0,
    diffusion_coefficient=1.0e-9 * gate.g4_units.m2 / gate.g4_units.s,
)

If the molecular configuration already exists, GATE updates its configurable properties. Otherwise it creates the necessary Geant4 molecule definition and configuration.

Add a dissociation channel

You can also add a dissociation channel:

chem_list.add_chemical_dissociation(
    parent="H2O2",
    products=["°OH", "°OH"],
    probability=1.0,
)

Current model

The current model is additive:

• start from one built-in Geant4 chemistry list;

• then append species, reactions, and dissociations in GATE.

If you customize species, reactions, or dissociations, make sure a valid base chemistry list is still selected.

Chemistry worlds and scavengers

The optional ChemistryWorld object defines the chemistry-space box used by Geant4 chemistry and can also carry background components and scavenger reactions.

You usually create it from a volume:

chemistry_world = sim.chemistry_manager.create_chemistry_world(volume=target)
chemistry_world.pH = 7

You can also create it explicitly from a box:

chemistry_world = sim.chemistry_manager.create_chemistry_world(
    translation=[0, 0, 0],
    half_size=[10 * gate.g4_units.um] * 3,
)

Once the chemistry world exists, you can define background components:

chemistry_world.add_component("H2O", 55.5)
chemistry_world.add_component("O2", 2.5e-2 * gate.g4_units.mole / gate.g4_units.liter)

and then add scavenger reactions:

dm3_per_mole_s = (
    1e-3 * gate.g4_units.m3 / (gate.g4_units.mole * gate.g4_units.s)
)

chemistry_world.add_scavenger_reaction(
    tracked_molecule="e_aq",
    scavenger="O2",
    products=["O2m"],
    rate_constant=1.74e10 * dm3_per_mole_s,
)

The scavenger species must already be present as a chemistry-world component. Scavenger products and tracked molecules must resolve to known chemistry species in the active Geant4-DNA chemistry setup.

pH is optional. At the moment, Geant4’s stock scavenger material expects an integer pH value.

This chemistry-world configuration is independent from confine_chemistry_to_volume:

• ChemistryWorld defines the chemistry-space box and optional scavenger reservoir;

• confine_chemistry_to_volume defines a GATE-side policy that kills chemistry tracks starting outside a chosen volume subtree.

Adding a chemistry actor

The current reference chemistry actor is ChemicalCountingActor. Further chemistry actors can be developed and made available in the future.

It is intended as a passive chemistry-scoring actor. Simulation-wide chemistry control stays on ChemistryManager.

It provides:

• chem6-like energy-loss and LET bookkeeping;

• chemistry-time species sampling;

• reaction counting;

• built-in and configured chemistry counter outputs.

A typical setup looks like this:

chem_actor = sim.add_actor("ChemicalCountingActor", "chem_actor")
chem_actor.attached_to = target
chem_actor.number_of_time_bins = 50

The actor must be attached to a volume. This volume defines the local region for the actor logic.

If you want the actor to request DNA EM on its attached volume itself, you can do:

chem_actor.track_structure_em_physics = "G4EmDNAPhysics_option2"

If you have already configured track-structure EM explicitly on the volume or region, you can leave this as None.

The actor also supports:

• track_only_primary

• primary_pdg_code

• energy_loss_min

• energy_loss_max

• min_kinetic_energy

• let_cutoff

• times_to_record

• number_of_time_bins

Example:

chem_actor.track_only_primary = True
chem_actor.primary_pdg_code = 11
chem_actor.number_of_time_bins = 50

Chemistry confinement

Chemistry confinement is configured globally through ChemistryManager, not per actor.

If you want chemistry tracks starting outside a given volume subtree to be killed before chemistry processing continues, configure:

sim.chemistry_manager.confine_chemistry_to_volume = target

where target can be either a volume object or a volume name.

This is a simulation-wide chemistry-control policy. It is intentionally kept outside chemistry actors so that chemistry actors can remain focused on probing and scoring.

Chemistry actor outputs

ChemicalCountingActor currently exposes five outputs:

• results

• molecule_counter

• reaction_counter

• configured_reaction_counter

• configured_species_counter

Summary output

The results output contains actor-level summary quantities such as:

• recorded_events

• chemistry_starts

• chemistry_stages

• pre_time_step_calls

• post_time_step_calls

• reaction_count

• killed_particles

• aborted_events

• accumulated_primary_energy_loss

• total_energy_deposit

• mean_restricted_let

• std_restricted_let

• species

• times_to_record

You can read it like this:

sim.run()
results = chem_actor.results.get_data()
print(results.reaction_count)
print(results.mean_restricted_let)

Counter outputs

The detailed chemistry histories live on the dedicated counter outputs:

• chem_actor.molecule_counter

• chem_actor.reaction_counter

• chem_actor.configured_reaction_counter

• chem_actor.configured_species_counter

All of these detailed counter outputs follow the same in-memory structure:

• a dictionary;

• keys are molecule or reaction labels;

• each value is a structured NumPy array with fields time and count.

Example:

molecule_data = chem_actor.molecule_counter.get_data()
reaction_data = chem_actor.reaction_counter.get_data()

for molecule_name, series in molecule_data.items():
    print(molecule_name)
    print(series["time"])
    print(series["count"])

The counts are cumulative.

Configured counters

The configured counters let you track selected observables without relying only on the stock Geant4 built-in counter classes.

configured_reaction_counter tracks cumulative counts of selected chemistry reaction signatures versus chemistry time.

chem_actor.counters.configured_reaction_counter.tracked_reactions = [
    {
        "reactant_a": "H",
        "reactant_b": "H",
        "products": ["H2"],
    }
]
chem_actor.counters.configured_reaction_counter.active = True

configured_species_counter tracks cumulative appearances of selected species versus chemistry time.

chem_actor.counters.configured_species_counter.tracked_species = [
    "e_aq",
    "O2m",
]
chem_actor.counters.configured_species_counter.active = True

This can be especially useful in scavenger workflows, where counting the appearance of signature product species is sometimes more robust than trying to count a direct reaction callback.

Naming conventions

Chemistry names appear at several layers, and the exact Geant4-facing labels should be treated as examples rather than a stable long-term naming contract:

• Python chemistry configuration usually uses names such as e_aq, O2, or O2m;

• molecule-counter outputs often expose Geant4 runtime names such as e_aq^-1, O_2^-1, or HO_2°^0.

Configured counters accept chemistry-species names that resolve through the Geant4 molecule table, while returned output keys may follow the runtime naming used by the underlying counter. As development evolves, these runtime labels may change. In practice, inspect the keys returned by get_data() rather than assuming one naming convention everywhere.

Activating and deactivating outputs

By default, the predefined counter outputs are active. You can deactivate them through the normal actor-output interface:

chem_actor.molecule_counter.active = False
chem_actor.reaction_counter.active = False
chem_actor.configured_reaction_counter.active = False
chem_actor.configured_species_counter.active = False

The same can also be controlled through the actor-owned counter objects:

chem_actor.counters.molecule_counter.active = False
chem_actor.counters.reaction_counter.active = False
chem_actor.counters.configured_reaction_counter.active = False
chem_actor.counters.configured_species_counter.active = False

For most user scripts, using the actor output interface directly is the clearest choice.

Complete example

This combines the main pieces in one small script:

import opengate as gate

sim = gate.Simulation()
um = gate.g4_units.um
keV = gate.g4_units.keV

sim.world.material = "G4_WATER"
sim.world.size = [1 * gate.g4_units.km] * 3

sim.physics_manager.physics_list_name = "G4EmStandardPhysics"
sim.chemistry_manager.chemistry_list_name = "G4EmDNAChemistry_option3"
sim.chemistry_manager.time_step_model = "IRT"

target = sim.add_volume("Box", "chem_box")
target.material = "G4_WATER"
target.size = [10 * um, 10 * um, 10 * um]
target.set_track_structure_em_physics("G4EmDNAPhysics_option2")
sim.chemistry_manager.confine_chemistry_to_volume = target

source = sim.add_source("GenericSource", "source")
source.particle = "e-"
source.energy.mono = 2 * keV
source.position.type = "point"
source.position.translation = [0, 0, 0]
source.direction.type = "momentum"
source.direction.momentum = [0, 0, 1]
source.n = 1

chem_actor = sim.add_actor("ChemicalCountingActor", "chem_actor")
chem_actor.attached_to = target
chem_actor.number_of_time_bins = 50

sim.run()

results = chem_actor.results.get_data()
species = chem_actor.molecule_counter.get_data()
reactions = chem_actor.reaction_counter.get_data()

print(results.reaction_count)
print(species.keys())
print(reactions.keys())

Current limitations

The chemistry support on this branch is already usable, but a few constraints are worth keeping in mind:

• one simulation currently needs one coherent chemistry list;

• one simulation currently has at most one global GATE-side chemistry confinement policy;

• chemistry counter writing to disk is not implemented yet on the python side;

• chemistry counter merged data currently uses a simple successive-run merge;

• scavenger-related observables may require configured species counting rather than direct reaction counting, depending on which Geant4 callback path the process uses;

• some advanced reaction-counter scenarios may still be fragile;

• ChemicalCountingActor currently assumes at most one molecule counter for its built-in species-sampling path.

These are implementation limits, not usage mistakes.

Practical tips

• Start with G4EmStandardPhysics plus regional track_structure_em_physics rather than trying to replace the whole physics list with a DNA-only list.

• Configure DNA EM only where you need it. Track-structure physics is much more expensive than standard condensed-history EM.

• Keep the chemistry list request consistent across the simulation and all chemistry actors.

• Use sim.chemistry_manager.confine_chemistry_to_volume if you want to restrict chemistry spatially at the GATE level.

• If you only need summary chemistry information, chem_actor.results is the easiest entry point.

• If you need time-resolved species or reaction histories, use the dedicated counter outputs.