The general-purpose Geant4 Monte Carlo toolkit and its Geant4-DNA extension to investigate mechanisms underlying the FLASH effect in radiotherapy: Current status and challenges

FLASH radiotherapy is a promising approach to cancer treatment that offers several advantages over conventional radiotherapy. With this novel technique, high doses of radiation are delivered in a short period of time, inducing the so-called FLASH effect – a phenomenon characterized by healthy tissue sparing without alteration of tumor control. The mechanisms behind the FLASH effect remain unknown. One way to approach this problem is to gain insight into the initial parameters that can distinguish FLASH from conventional irradiation by simulating particle transport in aqueous media using the general-purpose Geant4 Monte Carlo toolkit and its Geant4-DNA extension. This review article discusses the current status of Geant4 and Geant4-DNA simulations to investigate mechanisms underlying the FLASH effect, as well as the challenges faced in this research field. One of the primary challenges is to accurately simulate the experimental irradiation parameters. Another challenge is the temporal extension of the simulations. This review also focuses on two hypotheses to explain the FLASH effect – namely the oxygen depletion hypothesis and the inter-track interactions hypothesis – and discusses how the Geant4 toolkit can be used to investigate them. The aim of this review is to provide an overview of Geant4 and Geant4-DNA simulations for FLASH radiotherapy and to highlight the challenges that need to be overcome in order to better study the FLASH effect.


Introduction
FLASH radiotherapy is an emerging field of cancer treatment that delivers radiation doses to tumors in a very short amount of time, typically less than a second, in comparison to minutes using conventional radiotherapy. This promising new approach has the potential to improve radiotherapy treatment by minimizing damage to healthy tissue while maintaining the destructive effect on cancer cellsa feature known as the FLASH effect [1]. Reducing damage to surrounding healthy tissue may allow higher doses of radiation to be delivered to tumors, which could result in more effective treatment and improved outcomes for patients, as well as diminished toxicity in healthy tissues. Additionally, the reduced exposure time to radiation may help manage patient motion, leading to fewer side effects and a better quality of life for the patient.
High hopes are placed on FLASH radiotherapy, as it is viewed as a potentially groundbreaking advancement in the field of radiotherapy. However, being a relatively new area, this novel technique still requires a considerable amount of research and development to gain a complete understanding of the observed effects and their potential benefits in the clinics. The mechanisms behind the FLASH effect are not yet elucidated, making it a subject of multiple ongoing research. Understanding the underlying mechanisms would help to induce and optimize it. Simulations of the radiation effects can play an important role in guiding experimental efforts, while experimental results can validate, inform and improve the simulations. This interplay between simulation and experiment is crucial in advancing the field and gaining a deeper understanding of the FLASH effect.
What working hypotheses are currently being explored in the field of simulation to explain the FLASH effect in radiotherapy? One hypothesis is oxygen depletion, which proposes that the effect is a result of the rapid reduction of oxygen in the irradiated environment [2]. According to this hypothesis, the fast delivery of high doses of radiation in a FLASH treatment leads to the temporary depletion of oxygen in the irradiated tissue. This oxygen depletion affects water radiolysis and the formation of reactive oxygen species, which lowers tissue radiosensitivity. Tumors are generally already hypoxic, which would explain why they do not experience the same sparing effect as normal tissue [3]. Another hypothesis is inter-track interactions, which suggests that the FLASH effect is caused by the interaction of radiolytic species generated by adjacent ionizing tracks [4]. The high dose rate of FLASH radiotherapy increases the density of simultaneously produced species, which favors inter-track reaction processes and could impact the production of reactive oxygen species. Due to their hypoxic nature, tumor tissues do not undergo the same chemical reaction dynamics as normal tissues, which could explain the differential response. The tumor microenvironment differs from normal tissues in several other ways [5]. For example, tumor tissues have abnormal and disorganized vasculature [6]. Apart from the oxygen level, it becomes challenging to investigate hypotheses based on biological pathways using Monte Carlo simulation toolkits. However, it should not be forgotten that proposing mechanisms that explain the sparing effect of FLASH irradiation on normal tissues is one step. Still, understanding why tumors are not affected by these mechanisms is another.
The aim of this review is to assess the status of the general-purpose Geant4 toolkit [7][8][9] and its Geant4-DNA extension [10][11][12][13] with regard to FLASH radiotherapy. This includes evaluating the existing methods used for simulation and identifying areas that require further development. We will discuss the two hypotheses previously presented and determine the necessary steps to simulate ultra-high dose rate (UHDR) irradiation with specific characteristics. It should be noted that not all UHDR irradiation qualifies as FLASH, as the FLASH effect is a biological effect and not all UHDR conditions may exhibit it [14]. For this reason, we will prefer the term UHDR in what follows whenever the discussion is general. In addition, many reported experiments do not provide detailed irradiation parameters, which makes the validation of calculations against published experimental results difficult.

Proposed key steps to study FLASH irradiation
Geant4 is a general-purpose toolkit for simulating the transport of particles through matter using Monte Carlo techniques. Various fields of research use the Geant4 toolkit, including high energy physics, nuclear physics, space science, and medical physics [9]. Some of its key capabilities include geometry and material modeling, particle tracking, and interaction simulation. Geant4 can model complex geometries and a wide range of materials such as elements, compounds, mixtures, and their properties like density [9]. It is able to track the trajectories of many particle typesphotons, electrons, protons, and ionsand can simulate the physical interactions they undergo within matter. These include electromagnetic and hadronic interactions, and take into account the production of secondary particles.
Geant4 uses condensed history processes to simulate the behavior of particles in matter [15,16]. In the condensed approach, the interactions of particles are modeled using a continuous description, which includes a continuous energy loss component and the multiple scattering theory to model the angular deviation of the particle trajectory. The advantage of this approach lies in its ability to handle particles more efficiently than the discrete approach, which models each particle-matter interaction and tracks individually all secondary particles. The condensed history approach makes it possible to simulate particles passing through large volumes of matter, such as the human body, at a reasonable computational cost.
In some scientific research applications such as radiobiology or micro/nano-dosimetry, the condensed approach lacks accuracy, as continuous models do not describe the underlying physics in detail [17].
In such cases, a discrete description must be applied to precisely simulate the complex interaction pattern between ionizing radiation and its surrounding environment. Geant4-DNA is specifically designed for this kind of problem and includes comprehensive physics and chemistry models for simulating the interactions of ionizing radiation with simplified geometrical structures of biological targets. Mostly composed of water, biological systems undergo water radiolysis when exposed to radiation [18]. Geant4-DNA models water radiolysis in three main stages, identified as the physical, physico-chemical, and chemical stages (Fig. 1). During the physical stage, the ionizing particle excites and ionizes water molecules in a so-called track-structure mode of simulationthe discrete approach which simulates all interactions of the primary particle with water and tracks each secondary in the same way. From 1 fs to 1 ps post-irradiation, the excited and ionized water molecules dissociate, autoionize, or relax. The electrons thermalize and become solvated electrons. These processes occur during the physicochemical stage and generate radiolytic species. The latter diffuse and react with each other during the chemical stage, which is initially inhomogeneous and then becomes homogeneous, typically around a microsecond [19,20].
The default geometry used in Geant4-DNA to simulate water radiolysis is an ultra-large water cube to prevent particles and radio-induced molecules from reaching the edges of the simulation volume. This volume does not correspond to an irradiated sample, characterized by an aqueous solution confined in a vessel. The radiolytic species can freely diffuse and travel far from interaction sites. Over time, the species move away from each other, reducing the likelihood of mutual reactions as the ionizing particles are simulated independently. This is why the simulation of the radiolysis process in Geant4-DNA ends after 1 µs by default, a time considered to mark the completion of the inhomogeneous chemical stage [19,20]. The user has the option of extending the simulation beyond the default duration, while knowing that the results may become increasingly questionable as the extension goes beyond a certain point.
Presently, there is no established approach for simulating the dose rate effect on water radiolysis using Geant4-DNA; ionizing particles are regarded as independent, and therefore do not have the possibility to interact with each other in a direct or inter-track manner. In this context, we propose two key steps to apprehend the dose rate, and more specifically, the simulation of UHDR experiments. Fig. 2 displays the two main components of the approach: Geant4 to model the full-scale beamline and experimental setup; Geant4-DNA to simulate water radiolysis in a micrometer-scale volume irradiated by the Geant4-generated source. The main issue is always the interplay between the simulation complexity and the computing resource requirements. Due to current computational limitations, it is not possible to perform water radiolysis simulations using Geant4-DNA for real size experiments. Therefore, we need to merge the features of Geant4 with the capabilities of Geant4-DNA.
In this framework, Geant4 serves to generate a particle source tailored to the Geant4-DNA simulations, and specific to the studied irradiation modality. The idea is to model the entire beamline together with the aqueous sample under investigation. For example, this could be a radiation beam hitting a water tank containing the studied sample. The object of interest is the irradiated sample, in particular the energy and direction distributions of the ionizing particles passing through it. These distributions characterize the ionizing particle source for Geant4-DNA, but do not yet include the irradiation time structure, which should be added in a second step. Within Geant4-DNA, we indeed propose to define both the beam time structure and the sample microenvironment, typically to inform questions such as: does the sample contain oxygen, and in which concentration? should we simulate the presence of other scavenger molecules? does the beam time structure include pulses or is simply continuous? In all situations, time remains an issue. One pulse may last a couple of microseconds, but multiple consecutive pulses spread over a millisecond scale. As presented in Fig. 1, water radiolysis simulation as classically implemented in Geant4-DNA stops around 1 µs post-irradiation. Continuing the simulation for a few milliseconds by using the same method would not be appropriate as explained previously, not to mention the computational cost this undertaking would entail. We suggest that a combination between the Monte Carlo approach of Geant4-DNA and a deterministic approach for homogeneous chemistry is the way forward to study UHDR irradiation and its effects on systems over a longer duration. The proposed steps and issues outlined above will be discussed in more detail in the following three subsections.

The irradiation modality and beam structure
FLASH radiotherapy is a recent area of research, especially in the field of modeling. In the literature, studies on water radiolysis that examine the effects of UHDR conditions rely typically on a generic particle source delivering mono-energetic particles either isotropically or mono-directionally [22][23][24][25]: a fair approximation to investigate, from a fundamental point of view, the structure and impact of a radiation interaction pattern, but maybe less appropriate for extrapolating these results to experimental data collected in a specific irradiation configuration. Suppose that instead of a mono-energetic beam of particles, the sample of interest undergoes an electron irradiation characterized by a broad energy spectrum. In such a case, this aspect should be considered. Depending on their energy, the scattering of charged particles varies as they pass through water. In general, as particles lose energy, their linear energy transfer (LET) changes, which modifies the interaction pattern  [21]. Due to the simulation approach, Geant4-DNA sets a cut-off point for the chemical stage at 1 µs post-irradiation, which is regarded as the end of the inhomogeneous phase. and therefore the production of radiolytic species [26][27][28]. Within the context of FLASH effect studies, we should address this issue and, if relevant, take into account the energy spectrum and directional distribution of the particles within the irradiated sample.
Since version 11.1 of Geant4, an advanced example is available to simulate the ElectronFLASH linac (SIT, Vicenza, Italy) [29]. Also, a Geant4 model of the Oriatron eRT6 electron linac [30] was presented at the Flash Radiotherapy & Particle Therapy Conference (FRPT 2022) [31]. The eRT6 is an irradiation device that has been used in a wide range of research on FLASH radiotherapy [2,[32][33][34][35][36][37][38][39][40]. Modeling the irradiation modalities makes it possible to investigate UHDR experiments closer to the real conditions. For electrons, photons, protons, and heavy ions, Geant4 provides the means to compute the energy spectrum and directional distribution of primary and secondary charged particles within the irradiated sample. All these distributions effectively characterize the source irradiating the sample. However, this does not include a temporal structure, which is believed to be the key element at the origin of the FLASH effect. Despite numerous studies on the subject, we do not currently know the physical parameters that lead to the FLASH effect: it remains a biological effect and each UHDR modality must be tested on animal models before it can be considered a FLASH validated modality [14]. Modeling the beam structure requires careful consideration as to whether an UHDR or a FLASH irradiation is being investigated. This aspect is particularly important for consideration of the inter-track interactions hypothesis. If a deviation in the production of radiolytic species is observed due to inter-track interactions, it is essential to confirm that the modeled irradiation is credible within the context of FLASH radiotherapy. This ensures that the results accurately reflect the behavior of FLASH irradiation.
Practically, we can introduce a temporal structure into the simulation by sampling random times and associating them with the primary particles of the source. Each particle has therefore a time that corresponds to its simulation during water radiolysis in Geant4-DNA. The physical parameters of irradiation, such as the dose per pulse and repetition frequency for pulsed beams, often require the simulation of water radiolysis beyond several microseconds. We propose a few suggestions on this subject in subsection 2.2 on the temporal extension of water radiolysis simulations.

Temporal extension of water radiolysis simulations
How to approach the beam temporal structure at the millisecond/ second scale? This question is important to investigate the inter-track interactions hypothesis [4]. As a reminder, this hypothesis is based on the density at a defined time of radio-induced species. By increasing the dose rate, this density will increase which opens the possibility for species from adjacent tracks to interact. A different chemical reaction dynamic could impact the production of reactive oxygen species. To know if a FLASH irradiation induces inter-track interactions and to quantify its effect, it is crucial to properly model its temporal structure.
To be more precise, the temporal and spatial distributions of the ionizing particles are required to describe the situation accurately. Presently, Geant4-DNA models water radiolysis using a large water volume in which independent ionizing particles are simulated until 1 µs postirradiation. The user can vary this time point but to a certain extent. First, a limiting factor is the computational cost of increasing the duration after irradiation, especially for the step-by-step simulation method [41]. In addition, the dimensions of the simulation volume are such that the species move away from each other as time passes by, which progressively lowers their reaction probability. The current simulation mode of Geant4-DNA is not adapted for UHDR water radiolysis computations in the context of FLASH radiotherapy.
Tran et al. developed a parallel approach to overcome some of these limitations [25]. The model is available as a prototype in Geant4 version 11.1, until a dedicated example is provided to the community. The new model extends the water radiolysis simulation in the homogeneous chemical stage using a compartment-based approach together with the particle-based step-by-step model of Geant4-DNA. Briefly, after the physical and physico-chemical stages, the simulation enters a chemical phase where the geometry is such that the radio-induced species are confined within the simulation volume. This stage begins with the stepby-step model of Brownian diffusion to faithfully simulate the radioinduced spurs. From highly heterogeneous, the system progressively becomes homogeneous as species diffuse. To account for this evolution and to speed up the simulation, the volume is divided into voxels whose content is defined as homogeneous. The size of these voxels increases over time to reflect the level of homogeneity in the simulation volume. Once the remaining voxels merge into a single voxel, the whole system becomes homogeneous, and the Gillespie algorithm is applied to continue the simulation to the second or minute.
The approach developed by Tran et al. extends the water radiolysis simulation to the homogeneous chemical stage. The main limitation is however the restricted temporal structure of the beam that can be achieved (a pulse of a few tens of nanoseconds). The step-by-step method is resource intensive, in particular when the simulated dose increases. To reach a pulse width of 2 µs for example, we should consider the independent reaction times (IRT) method of Geant4-DNA [42].
The IRT option considerably lowers the computational cost of the heterogeneous chemical stage, with the drawback of losing the information on the spatial distribution of the species. An alternative to the aforementioned model would be to modify the Geant4-DNA version of the IRT method in order to include the possibility of simulating particle tracks arriving successively. This feature is not available in the current implementation, as it does not record the position of species during the simulation. The IRT approach is an interesting solution that should be developed as well as its combination, for example, with a deterministic model (differential equation solver) to simulate in a simple way the homogeneous chemistry following the end of the irradiation. It is worth mentioning that other codes besides Geant4-DNA are being developed for the same purpose, such as TOPAS-nBio (based on Geant4-DNA) and TRAX-CHEM [22,43].

Water radiolysis with scavengers and oxygen depletion
Properly modeling the particle beam and its temporal structure is key to studying the inter-track interactions hypothesis [4]. However, with regard to the oxygen depletion hypothesis [2], simulating the radiolysis of pure water remains insufficient. A way to include molecules like oxygen in the radiolysis process is the use of scavengers [44]. Instead of considering oxygen molecules as individual point-like objects within a continuum (water), the scavenger model simulates oxygen as a continuous medium in itself. The assumption underlying this approach is therefore a homogeneous distribution of molecules characterized by a given concentration. It is possible to simulate the oxygen molecules individually [24], but at the expense of a high computational cost. The current method used in Geant4-DNA is the scavenger model.
Since Geant4 version 11.0, an example is provided to easily handle scavengers with the IRT method: examples/extended/medical/dna/ scavenger [45]. The public version for the step-by-step method is not yet available. But both methods were validated and agreed within 10 % relative difference for two different systemsthe combination of nitrite and nitrate scavengers and the computation of oxygen depletion at different initial concentrations of oxygen in water [45]. The need to determine oxygen depletion stems from the hypothesis aiming to explain the FLASH effect by the depletion and repletion of oxygen in the area being subjected to irradiation. As a reminder, this hypothesis states that FLASH irradiation quickly depletes oxygen which makes the irradiated medium hypoxic and lowers its radiosensitivity. With the FLASH mode of irradiation, the organism does not have time to reoxygenate the affected area during irradiation, which differs from the conventional mode where irradiation is much slower and allows constant reoxygenation.
To compute oxygen depletion for conventional irradiation, it is sufficient to simulate independent particle tracks. Indeed, for an electron beam at minimum ionizing energy, equation (1): fluence = dose mass stopping power (1) makes it possible to estimate the fluence for a given dose delivered at a certain dose rate, using mass_stopping_power = 2 MeV cm 2 /g [46]. For a 5 Gy electron irradiation delivered at 0.1 Gy/s, the fluence is approximately 156 µm − 2 , resulting in an average time interval of 321 ms between two consecutive electrons impacting a 1 µm 2 area. This example shows that the interval is sufficient to neglect inter-track interactions between species from different tracks. For a proton beam, the calculation is similar. If we consider the PSI Comet Cyclotron, which can deliver 235 MeV protons at 0.1 Gy/s [47], the average time interval is approximately 641 ms. The oxygen depletion computed using Geant4-DNA and the independent particle method agrees well with experimental data obtained under conventional irradiation and with the results of the TRAX-CHEM Monte Carlo code [45,48]. To compute oxygen depletion for a FLASH irradiation, we suggest following the steps presented in Fig. 2 and explained in section 2, in order to combine the scavenger model with the one characterizing the source of irradiation both energetically and temporally. Many ongoing developments of Geant4 are directed towards this purpose. However, we do not expect a high level of oxygen depletion in FLASH compared to conventional irradiation. The literature tends to show that this depletion is even lower in UHDR mode [24,49,50], which does not support the oxygen depletion hypothesis to explain the FLASH effect.

Challenges and current/future developments
To validate Monte Carlo codes such as Geant4-DNA for the study of the FLASH effect, it is essential to compare simulation results on the same timescale as the UHDR experiments. Ongoing investigations show that H 2 O 2 yields around the microsecond scale, and computed using an UHDR electron beam as well as the Geant4-DNA step-by-step method, do not differ from conventional irradiation. On a timescale around the millisecond, no decrease in H 2 O 2 is also observed using the model developed by Tran et al. [25]. This is in contradiction with measurements made with the eRT6 beam, which show a decrease in secondary H 2 O 2 yield under UHDR irradiation compared to conventional irradiation [2,47]. The same behavior was observed under proton irradiation with the ARRONAX facility [51]. However, ongoing measurements of primary yields with the eRT6 beam do not seem to repeat this decrease in H 2 O 2 [52]. This suggests that the differential effect between UHDR and conventional dose rate modes may occur later in the homogeneous chemical phase. This should be confirmed not only for the H 2 O 2 yield but also for other radiolytic products such as • OH and e − aq . To our knowledge, such measurements are still missing under FLASH conditions.
The ability to simulate beyond the μs-ms timescale should be introduced into Geant4-DNA as an example. The objective of the model developed by Tran et al. was to extend the water radiolysis simulation to the homogeneous chemical phase. However, this model does not allow the simulation to be extended beyond the millisecond or the second, depending on the complexity of the reactions involved. Ramos-Méndez et al. also showed the ability of the IRT method to extend the computations to the millisecond to simulate the Fricke dosimeter system [22]. To develop models to longer timescales, differential equation solvers, which are frequently applied to study homogeneous chemistry, can also be used. This also allows the study of much more complex chemistry involved in radiobiology, as for example in the study by Labarbe et al. [53]. In conclusion, we see that a comprehensive model-based study of the FLASH effect needs to combine valid codes on different timescales, as is the case in the recent study of Espinosa-Rodriguez et al. [54].
A limitation of Geant4-DNA that applies not only to FLASH effect studies, but to any type of water radiolysis or radiobiological modeling, is the energy range of the physics models that may not cover the beam energy considered in some experiments. Developments to extend the electron models up to 10 MeV and the proton models up to 300 MeV in Geant4-DNA are ongoing [55].
Although this review article has focused on the ability of Geant4 and Geant4-DNA to model water radiolysis with oxygen and under FLASH conditions, it is worth discussing the ability of the toolkit to compute more complex effects involved in radiobiology. Various works on DNA damage computation with Geant4-DNA have been published but, to our knowledge, damage computation under FLASH conditions has not yet been performed [56][57][58][59]. In this context, it would be necessary to properly model the presence of scavengers and not to approximate it using a kill distance parameter as is often the case in this type of simulation. Furthermore, DNA damage simulations should be extended beyond a few nanoseconds. Other complex chemical processes, such as lipid peroxidation, should also be taken into account in future simulations. Recent experimental data obtained with the eRT6 beam showed a clear differential effect on lipid peroxidation that should be investigated by simulations [60]. Geant4-DNA is a useful simulation tool for studying the effects of radiation on aqueous solutions and simplified biological structures. However, its ability is limited in complex cases such as tumor tissues. This highlights the multidisciplinary nature of FLASH radiotherapy research. Understanding the differential response of normal and tumor tissues requires biological investigations on which simulations can build a simplified test model.

Conclusion
FLASH radiotherapy is an emerging field of research with great potential for improving patient outcomes. Despite its promising benefits over traditional radiotherapy methods, much remains to be understood about this new modality and the effects of ultra-high dose rate irradiation. Why do certain irradiations not induce the FLASH effect, although delivered at an ultra-high dose rate? What are the physical parameters of the beam to be taken into account? Modeling the radiation effects on aqueous solutions, whether by Monte Carlo techniques, a deterministic approach or a combination of both has the advantage of providing the possibility to test the impact of different irradiation configurations. The development of the Geant4 toolkit and its Geant4-DNA extension is ongoing but only a limited number of studies are currently published. This highlights the need to strengthen research and development in the area of FLASH radiotherapy in order to explore its potential as a clinical tool for routine cancer treatment. With the growing interest and collaboration among researchers but also the increased investment in this field, it is likely that we will see significant progress in the near future, ultimately advancing our understanding of FLASH radiotherapy.

Declaration of Competing Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.