Objectives. This study aims to demonstrate the feasibility of combining spatially fractionated radiotherapy (SFRT) and FLASH radiotherapy (FLASH-RT) techniques using very high energy electron (VHEE) through Monte Carlo simulations and to develop a simulation platform for this purpose. Approach. The simulation platform was built in TOPAS based on the beam parameters of the Wuhan Advanced Light Source (WALS) linear accelerator. Quadrupole magnets are used to focus VHEE beams to control the longitudinal depth of the maximum dose. Lateral scanning irradiation is achieved by combining deflecting magnets. The spread-out electron peak (SOEP) depth adjustment, and spot scanning were achieved through magnetic field adjustments. The dose distribution of VHEE in the water phantom was simulated, and the three-dimensional dose distribution, dose volume histogram (DVH), and peak-to-valley ratio were analyzed. Main results. Based on the fitted relationships, Monte Carlo simulations with 106 particles were conducted in TOPAS to verify the positions of maximum dose, yielding millimeter errors. The SOEP flatness between 12–17 cm depths exceeded 98.7%. SFRT simulated at four target points showed a peak-to-valley dose ratio (D10%/D90%) of 11.026, where D10% is the minimum dose within the highest 10% dose voxels, and D90% is the minimum dose received by 90% of the target volume. A volume fraction of 4.41% within the PTV received ⩾90% of the maximum dose. The dose rate in the PTV reached up to 54.25 Gy s−1 under 200 nA electron beam irradiation. Significance. Take the example of beam parameters of the linear accelerator at WALS, the simulation platform is validated for conducting research on combined SFRT and FLASH-RT for deep tumor radiation treatment. This also demonstrates the feasibility of combining SFRT and FLASH-RT using VHEE.
Export citation and abstractBibTeXRIS
Radiation therapy is a key component in the treatment of various cancers, with the potential to precisely target tumors. This technique relies on balancing the maximization of radiation dose to the tumor while minimizing damage to surrounding healthy tissues and critical organs (Eufemon Cereno et al 2023, Zivkovic et al 2025). In radiation therapy for deep-seated tumors, the radiation must penetrate to depths of 15–30 cm within the body (Whitmore et al 2021). The concept of using very high energy electron (VHEE) beams, which encompass energies from 50 to 250 MeV, was first introduced for radiotherapy applications in 2000 (DesRosiers et al 2000). This technique not only achieves sufficient penetration to treat lesions at typical clinical depths but also produces a lateral penumbra that is comparable to that of state-of-the-art megavoltage photon radiotherapy (DesRosiers et al 2000, Böhlen et al 2021, 2024), and they are easily adapted for electromagnetic (EM) scanning in pencil beam configurations (Panaino et al 2025). Compared to proton beams, VHEEs provide similar, and sometimes even better, target coverage (Yeboah et al 2002, Bazalova-Carter et al 2015, Palma et al 2015, 2016, Kokurewicz et al 2019). They also exhibit reduced sensitivity to tissue heterogeneities (Papiez et al 2002, Bazalova‐Carter et al 2015, Lagzda et al 2017, Kokurewicz et al 2019).
VHEE can be easily manipulated by magnetic components. The charge-to-mass ratio of electrons is 1836 times higher than that of protons (Böhlen et al 2024), resulting in approximately 34 times lower magnetic rigidity for electrons at similar energies (Desrosiers 2004, Ronga et al 2021). Compared to protons, electrons can be manipulated and accelerated in a significantly weaker magnetic field, resulting in a smaller trajectory. Quadrupole magnets focus the VHEE, enabling a concentrated dose distribution within the deep target region (Kokurewicz et al 2019, Whitmore et al 2021, An et al 2024). Using EM scanning, VHEE can be precisely directed to specific locations along both horizontal and vertical planes, facilitating intensity modulated VHEE therapy (Panaino et al 2025). Pencil beam scanning allows for the delivery of VHEE, enabling the implementation of spatially fractionated radiotherapy (SFRT, Fischer et al 2024). Pencil beam scanning for tumor radiation therapy eliminates the production of neutrons, photon scattering, and electron contamination that can occur with grid systems in collimation, thereby minimizing their impact on the treatment area (Karimi et al 2022).
Most traditional treatment methods still provide relatively uniform dose distribution within the tumor volume (Chargari et al 2020). When treating large tumors, these conventional approaches often result in unnecessary high dose exposure to surrounding healthy tissues, increasing the risk of toxic reactions (Price et al 2023). SFRT is an emerging treatment modality that employs non-standard dose delivery techniques to reduce normal tissue toxicity while improving the therapeutic ratio (Schneider et al 2022). It creates a highly uneven dose distribution within the tumor, inducing biological effects such as the bystander effect (Asur et al 2015) and immune activation (Lukas et al 2023), which enhances the tumor cell killing efficiency while minimizing damage to surrounding healthy tissues(Yuan et al 2025). Non-standard dose delivery refers to a method that combines high-dose regions (peaks) and low-dose regions (valleys) to produce an inhomogeneous dose distribution (Billena and Khan 2019, Fischer et al 2024). Several dosimetric parameters have been identified that may predict normal tissue sparing, including the peak-to-valley dose ratio (PVDR), beam spacing, and beam size (Anderson et al 2012, Smyth et al 2016, Fischer et al 2024). PVDR is the ratio of the peak dose to the valley dose, while beam spacing refers to the center-to-center distance between peaks. Recent studies have shown that the valley dose is a key parameter in defining normal tissue protection (Rivera et al 2020, Fernandez-Palomo et al 2022). Spatial fractionation of the dose is typically achieved by dividing the irradiated field into several narrow beams, separated by small gaps, usually 1–4 times the beam width. The beam can be pencil-shaped (narrowing along all lateral dimensions) or planar (narrowing along only one lateral dimension) (Schneider et al 2022).
FLASH radiotherapy (FLASH-RT) radiation therapy requires the delivery of ultra-high dose rates (⩾40 Gy s−1) within microsecond timescales, aiming to damage tumor cells while sparing surrounding healthy tissues (Böhlen et al 2022, 2023, No et al 2022, Valdés Zayas et al 2023). Although based on different mechanisms, both FLASH-RT and SFRT have been shown to induce radiobiological effects that are significantly different from those of conventional radiation therapy. Since the therapeutic potential of each technique has been individually established, the combination of FLASH-RT and SFRT may represent a promising alliance (Schneider et al 2022).
The key contributions of presented work include
Development of a simulation platform for scanning irradiation of tumors using focused VHEE beams.Design of a research platform for tumor radiation therapy utilizing the simulation platform, providing a reliable solution for applying radiation therapy to tumors with the linear electron accelerator at Wuhan Advanced Light Source (WALS).Provision of an effective simulation tool for the future use of VHEE beams in the SFRT of deep tumors.Validation of the feasibility of combining SFRT with FLASH-RT, based on the beam extraction parameters of the linear electron accelerator at WALS.In summary, this platform offers a simulation tool for studying SFRT using focused VHEE beam pencil scanning irradiation. The proposed simulation platform provides a physical design for the application of the WALS linear accelerator in tumor radiation therapy. The feasibility of using VHEE for SFRT and FLASH-RT treatment of deep-seated tumors has been demonstrated for the first time.
2.1. Electron beam parametersThe WALS linear accelerator uses a photoinjector to produce the beam (Dai et al 2023), and the LINAC is able to provide VHEE beams with an energy of 250 MeV and even higher. This enables research on the use of focused VHEE beams for radiation therapy of deep-seated tumors. The extracted electron beam parameters are shown in table 1, and the application platform is designed based on the input parameters. The temporal structure of the beam is shown in figure 1(a).
Figure 1. (a) Temporal structure of the beam from the WALS linear accelerator (b) model diagram built in TOPAS.
Download figure:
Standard image High-resolution imageTable 1. 250 MeV electron beam parameters.
ParameterValueUnitEnergy250MeVEnergy spread1218keVHorizontal emittance0.731
Vertical emittance0.736
Horizontal beam size
1.516mmVertical beam size
1.515mmHorizontal correlated divergence1.9567 × 10−2mradVertical correlated divergence1.9613 × 10−2mradBunch charge0.1–2nCMacro-pulse repetition frequency1–100HzMicro-pulses per macro-pulse1–4/Micro-pulse separation2ns
and
represent the standard deviation of the beam profile in the x and y directions, respectively. The average beam current will be 200 nC s−1 when the single bunch charge is 2 nC repeating at 100 Hz.
The position of the maximum dose is set at 17 cm, with the clinical treatment source-to-skin distance (SSD) set to 60 cm (which refers to the distance between the vacuum chamber exit and the water phantom entrance in this model), and 7 cm allocated for the installation space of the dose detector and other beam measurement components. The model set in the Monte Carlo code is shown in figure 1(b), where T1 to T6 are quadrupole magnets, B1 and B2 are deflecting magnets, and D1 to D12 are detectors. All parts of the model, except for the water phantom and air gap, are set to be in vacuum. The dimensions of all magnets are listed in supplementary table S1.
2.3. Program for focused VHEE scanning irradiation of the target areaRadiation dose calculation has traditionally relied on algorithms of varying complexity, including convolution/superposition methods and Monte Carlo simulations (Zheng et al 2025). The TopasOpt (Whelan et al 2023) code was further developed by embedding the optimization program within the NelderMeadOptimizer (Chang and Chang 1996), as shown in figure 2. The optimization parameters include the magnetic field gradients of the deflection magnets B1 and B2, the magnetic field gradients of quadrupole magnets T1 to Tn, and the distances between components L1 to Ln + 3. The Monte Carlo code TOPAS (Perl et al 2012) is then executed.
Figure 2. Flowchart of the program for focused VHEE scanning irradiation of the target area.
Download figure:
Standard image High-resolution imageCalculate the position error
:

where
represents the gradient of the magnetic field of the deflecting magnets, the gradient of the magnetic field of the quadrupole magnets, their dimensions, and the distances between components. (
,
,
) is the coordinate of the point with the maximum dose within the target area, (
,
,
) is the coordinate of the target position.
Calculate the full width at half maximum (FWHM) along the longitudinal direction, which in this study refers to the z-axis
:

where
is the z-coordinate corresponding to half of the maximum dose value on the extracted dose-depth curve, located after the point of maximum dose.
is the z-coordinate corresponding to half of the maximum dose value on the extracted dose-depth curve, located before the point of maximum dose.
Calculate the difference in the FWHM of the x–y profile at the location of the maximum dose
:

where
is the FWHM in the x−direction of the x−y profile at the location of the maximum dose, and
is the FWHM in the y-direction of the x–y profile at the location of the maximum dose.
Calculate the beam envelope error

where
is the actual beam envelope obtained through tracking, and
is the set limit for the beam envelope. This error is calculated when
.
The objective function of this optimization program
:

where
is the weight for the position error,
is the weight for the longitudinal FWHM,
is the weight for the difference in the FWHM of the x–y profile at the location of the maximum dose, and
is the weight for the beam envelope error.
Firstly, we calculate the total dose accumulated from multiple beams in the smallest grids along the central axis. In this study, it refers to the total dose within a 0.297 cm × 0.297 cm × 30 cm region along the target axis,

where
is the total dose,
is a coefficient factor (varied),
is the dose profile corresponding to the ith electron beam peak position used in the superposition, N is the number of beams.
Next, the dose homogeneity within this region was assessed. The average dose
was computed, followed by the standard deviation
of the dose values about this mean value:


These values were then used to define a flatness
:

The Differential Evolution method (Qiang and Mitchell 2014) was employed to find the optimal set of weights
that minimized the flatness objective function. The input data for this optimization (Data1.csv to Datan.csv) were generated by a program simulating VHEE pencil beam scanning (An et al 2024).
A PTV with dimensions of
is defined within the longitudinal depth range of 12 cm–17 cm in the water phantom. The remaining regions of the water phantom are designated as Non-PTV. Four target points were placed as shown in figure 3, with their center positions located at (−1.5 cm, −1.5 cm), (−1.5 cm, 1.5 cm), (1.5 cm, −1.5 cm), and (1.5 cm, 1.5 cm), respectively. The treatment plan is based on typical geometries and dose measurement parameters used for SFRT with photons (Wu et al 2020, Zeng et al 2025).
Figure 3. Layout of four target points within a 4 cm × 4 cm region.
Download figure:
Standard image High-resolution image 2.6. Dose distribution analysisThe WALS linear accelerator delivers a pulsed beam. Calculations based on the parameters in table 1 indicate that the average electron charge delivery rate can reach up to 800 nC s−1, while the pulsed beam current could reach as high as 2 nC/2 ns, i.e. 1 C s−1. For this study, simulations were conducted using the average number of electrons per second. Within the adjustable range, a bunch charge of 2 nC and a macro-pulse repetition frequency of 100 Hz were selected, resulting in an average electron beam current of 200 nC s−1. SFRT is performed using 200 nC s−1 electron beams on the PTV. The dose distribution of the PTV is calculated using Monte Carlo simulations. The feasibility of combining SFRT and FLASH-RT is validated through the analysis of the percent depth dose (PDD) curves, dose profile distributions, peak-to-valley ratios along the target axis, and average dose rates. In the case of VHEE beams, a grid of 101 × 101 × 101 voxels is set to reduce the computation time for dose distribution calculations using the ‘scorer’(An et al 2024). The relationships between the deflecting magnet B1 and the y-direction maximum dose position (max_y), as well as between B2 and the x-direction maximum dose position (max_x), have been validated using a grid of 501 × 501 × 101 voxels. The relationship between T1 and the depth of maximum dose (Dmax_depth) has been validated using a grid of 101 × 101 × 501 voxels.
3.1. Model optimization resultsThe positions of the components and the magnetic field gradients of the quadrupole magnets, as obtained from the program optimization, are presented in table 2.
Table 2. The parameters of model.
ParameterValueUnitDistanceL185cmL216cmL340.20cmL468.55cmL555.18cmL658.45cmL766.81cmL853cmL930cmMagnetic field gradientT2−0.0640T cm−1T30.1722T cm−1T4−0.1705T cm−1T50.1459T cm−1T6−0.2069T cm−13.2. Relationship between the magnetic field and the position of the maximum doseA minimal set of variables is selected to control the position of the maximum dose, based on the beamline layout described above, in order to reduce errors caused by adjusting multiple magnetic fields during radiation therapy. The magnetic field gradient of the last quadrupole magnet is adjusted to control the longitudinal depth position of the maximum dose (Whitmore et al 2021, An et al 2024). The magnetic field of the deflecting magnet B2 is varied to collect data on the x-direction maximum dose position max_x, while B1 is varied to collect data on the y-direction maximum dose position max_y. The collected data are linearly fitted, as shown in figure 4.
Figure 4. Linear fits between magnetic field parameters and dose position metrics. (a) B2 vs x-direction maximum dose position (max_x): y = 72.73x–0.37, R2 = 1. (b) B1 vs y-direction maximum dose position (max_y): y = 61.58x–0.35, R2 = 1. (c) T1 vs depth of maximum dose (Dmax_depth): y = −49 556.77x + 1976.92, R2 = 0.93.
Download figure:
Standard image High-resolution image 3.3. SFRT combined with FLASH-RTThe required deflection magnetic fields for the target point at (1.5 cm, 1.5 cm) were calculated, yielding B1 as 0.2493 T and B2 as 0.2113 T. Depths of 12 cm, 14.5 cm, and 17 cm were selected for validation using
particles, and the results are shown in table 3.
Table 3. Validation results with
particles at depths of 12 cm, 14.5 cm, and 17 cm.
The coefficients of each curve required to achieve the SOEP from 12 cm to 17 cm at the target point located at (1.5 cm, 1.5 cm), as well as the resulting flatness value, are presented in table 4.
Table 4. The coefficients and flatness result for SOEP at the target point (1.5 cm, 1.5 cm).
ParameterValueUnit
0.1643113/
0.6766686/
98.89%During the SFRT using four target points, the magnetic field directions of the quadrupole and deflecting magnets remain parallel to the x and y axes, respectively. This scanning irradiation creates regions of relatively high dose outside the target volume, as shown in figure 5. Significant differences in the dose distributions along the x and y directions were observed through the analysis of supplementary figures 1 and S2. Upon entering the water phantom, the dose distribution in the y direction is larger in size. When the target points are spaced 3 cm apart, the overlap in the y direction results in regions of relatively high dose outside the PTV.
Figure 5. Irradiation results with quadrupole and deflecting magnet fields aligned parallel.
Download figure:
Comments (0)