---
title: 'pAvalancheMC: RPC Monte Carlo Simulation'
url: https://www.emergentmind.com/topics/pavalanchemc
type: topic
---

# pAvalancheMC: RPC Monte Carlo Simulation

Searching arXiv for the specified pAvalancheMC papers and closely related RPC simulation work.
I’ll gather the relevant arXiv records for the Garfield++/RPC pAvalancheMC work and its thesis follow-up before writing the article.
Searching arXiv now.
pAvalancheMC is a parallel Monte Carlo avalanche class introduced in Garfield++ for Resistive Plate Chamber (RPC) simulation in regimes where dynamic space-charge effects materially modify the applied field. It was developed to make 3D particle-based avalanche calculations tractable in operating domains that include avalanche growth, saturated avalanche, and streamer-like onset, while remaining integrated with standard Garfield++ components such as neBEM, MAGBOLTZ, TrackHeed or HEED, and ROOT random-number infrastructure. In its published form, pAvalancheMC combines OpenMP parallelization with a line-charge representation of the evolving charge cloud; in the later thesis treatment, this formulation is extended by image-charge modeling for layered dielectric electrodes, so that the total field can include applied, space-charge, and polarization contributions [2211.06361] [2507.22225].

## 1. Origin, scope, and detector context

pAvalancheMC was introduced to address a specific deficiency in conventional Garfield++ avalanche tools for RPC studies: large avalanches are computationally expensive, while the evolving electron and ion cloud can substantially distort the local electric field and thereby alter transport, multiplication, and signal formation. The class targets RPC operating regimes in which the number of charges becomes large enough to modify the applied field, including the transition from ordinary avalanche growth to saturated avalanche and to streamer-like behavior [2211.06361].

The detector context is the parallel-plate RPC, including single-gap and timing RPC variants with resistive or dielectric electrodes such as bakelite or glass, graphite coating, and segmented pickup strips. The thesis formulation explicitly situates the method in the ICAL/INO program, where glass RPCs operate predominantly in avalanche mode and detector optimization depends on accurate estimates of timing, efficiency, charge spectra, and rate capability under space-charge and electrode-polarization effects [2507.22225].

A central distinction of pAvalancheMC is that it remains a Monte Carlo transport-and-multiplication engine rather than a fully microscopic collision-by-collision simulation. This allows substantially faster event generation than AvalancheMicroscopic while retaining field-dependent transport, stochastic multiplication, and induced-signal calculation. The intended application domain is therefore detector-level performance simulation in parameter regimes where self-consistent field updates are indispensable.

## 2. Field model, transport physics, and avalanche growth

The basic field decomposition in the original Garfield++ implementation is
$$
E_t = E_a + E_s,
$$
where $E_a$ is the applied field and $E_s$ is the dynamically recomputed space-charge field. In the thesis extension, the total field includes an additional polarization term from image charges in the dielectric stack, so that transport can be updated under the modified local field at every time step [2211.06361] [2507.22225].

Transport and multiplication coefficients are taken from MAGBOLTZ. The quantities used are the drift velocity $v_d(E)$, longitudinal and transverse diffusion coefficients $D_L(E)$ and $D_T(E)$, Townsend ionization coefficient $\alpha(E)$, and attachment coefficient $\eta(E)$. The diffusion model is Gaussian and anisotropic under field. In one formulation,
$$
\sigma_L = D_L \sqrt{v_d \delta t}, \qquad \sigma_T = D_T \sqrt{v_d \delta t},
$$
while the thesis also gives the anisotropic Gaussian forms
$$
\phi_L(z,l) = \frac{1}{(\sqrt{2\pi l}) D_L} \exp\!\left(-\frac{(z-z_0)^2}{2 D_L^2 l}\right),
$$
and
$$
\phi_T(r,l) = \frac{1}{D_T^2 l} \exp\!\left(-\frac{(r-r_0)^2}{2 D_T^2 l}\right).
$$

Avalanche growth is governed by the effective Townsend coefficient $\alpha_{\mathrm{eff}}=\alpha-\eta$. The expected mean electron population along a drift path is written as
$$
N(x) = N_0 \exp\!\left(\int (\alpha(E(x))-\eta(E(x)))\,dx\right),
$$
and the stochastic implementation uses a modified Yule–Furry process with attachment. In the thesis description, each drift step is subdivided so that the subdivision length is approximately $(\alpha+\eta)^{-1}$, with probabilities $p=\alpha/nDiv$ and $q=\eta/nDiv$ for ionization and attachment, respectively [2211.06361] [2507.22225].

This formulation makes saturation a direct field-feedback effect rather than an imposed cap. As space charge accumulates, $E_t$ is reduced in parts of the avalanche, lowering $\alpha$ and increasing the relative importance of $\eta$. When $\alpha-\eta$ becomes negative in most voxels, net growth stalls. In Ar/CO$_2$ (80/20), the reported saturated avalanche shows a knee around $17.96$ ns at an electron count of approximately $10^8$, close to the Raether limit traditionally associated with streamer transition [2211.06361].

## 3. Software architecture and algorithmic workflow

pAvalancheMC is built as a Garfield++ class that interoperates with the standard RPC simulation stack. The component roles described in the sources are summarized below.

| Component | Function | Notes |
|---|---|---|
| neBEM | Static and weighting fields | OpenMP-parallelized; FastVol interpolation |
| MAGBOLTZ | Transport and multiplication tables | $v_d$, $D_L$, $D_T$, $\alpha$, $\eta$ |
| TrackHeed / HEED | Primary ionization clusters | Positions and generation times |
| ROOT TRandom3 | Per-thread random streams | Mersenne Twister, period $2^{19937}-1$ |
| OpenMP | Parallel execution | Drift, gain, field update, signal accumulation |

A typical run initializes geometry and static fields, loads transport tables, configures OpenMP threads and per-thread TRandom3 instances through `SetNumberOfThreads(n)`, generates primary clusters, and acquires weighting-field maps. Primaries are sorted by generation time and started when the simulation clock satisfies $n\Delta t \ge t_i$. The time-stepping loop then advances electrons, applies diffusion, samples multiplication and attachment, updates the charge lists, and recomputes the field when the charge threshold for space-charge activation is exceeded [2211.06361].

The configuration described in the published examples includes a time step $\delta t = 20$ ps, a space-charge activation threshold around $N_e \ge 10^4$, and cylindrical segmentation parameters with typical values $\delta z \approx 0.001$ cm and $\delta \phi \approx 1^\circ$. Garfield++ interfaces named in the sources include `SetTimeStep(20e−12)`, `EnableSpaceCharge(true)` or `SetSpaceChargeOn(true)`, `SetSpaceChargeThreshold(1e4)`, `SetVoxelization(dz, dr, dphi)` or `SetGrid(zStep, rStep, phiStepDeg)`, `SetImageChargeOn(true/false)`, `SetImageCut(percent)`, `SetElectrodeProperties(...)`, and `RunAvalanche(...)` [2211.06361] [2507.22225].

OpenMP parallel regions are applied to electron drift and diffusion, multiplication and attachment sampling, space-charge field summation, and induced charge or current accumulation. The per-thread random-number strategy is explicit: TRandom3 instances are allocated per thread with distinct seeds to avoid correlated random streams.

## 4. Space-charge representation and dielectric polarization

The distinctive numerical device in pAvalancheMC is the cylindrical voxelization of the avalanche charge cloud. Electrons and ions are binned into $S$ slices in $z$, $R$ concentric radial rings, and $L$ azimuthal sectors with $L \delta \phi = 360^\circ$. Each occupied voxel is represented by a straight line-charge segment with constant linear charge density $\lambda$, and curved azimuthal segments are approximated by straight segments of length $r\delta \phi$ [2211.06361].

In the original implementation, the field of a finite line segment is evaluated numerically through
$$
E(P) = \frac{\lambda}{4\pi \epsilon} \int_{\mathrm{segment}} \frac{P-r'}{|P-r'|^3}\,dl',
$$
and voxelwise field components $E_s|_z$, $E_s|_r$, and $E_s|_\phi$ are accumulated. The paper notes that $E_s|_\phi$ is generally small, but not identically zero; one example gives a maximum $\lvert E_s|_\phi \rvert \approx 0.26$ of $\lvert E_s|_z \rvert$ at $10.46$ ns, indicating axial asymmetry [2211.06361].

The thesis extends this by introducing a segmented-line model together with image charges in layered dielectrics. The RPC stack is treated as a five-layer system: graphite, dielectric electrode, gas, dielectric electrode, graphite. Reflection and equivalence factors at interfaces are
$$
\alpha_{mn} = \frac{\epsilon_m-\epsilon_n}{\epsilon_m+\epsilon_n}, \qquad
\beta_{mn} = \frac{2\epsilon_n}{\epsilon_m+\epsilon_n},
$$
with $\alpha_{21}\to -1$ for a conductor. Image charges are generated iteratively and truncated using a percentage-change criterion,
$$
\Delta S_m = \frac{S_m-S_{m-1}}{S_{m-1}} \times 100\%,
$$
with a representative cut of $20\%$ [2507.22225].

The same source gives the equivalent RPC relaxation time
$$
\tau_g = 2 R_b (2 C_b + C_g) = \rho \epsilon_0 \left(\epsilon_r + \frac{d}{g}\right),
$$
and reports the example $d=g=2$ mm, $\epsilon_r=5$, $\rho \approx 2\times 10^{10}\ \Omega\cdot \mathrm{m}$, yielding $\tau_g \approx 1.05$ s. The static polarization field is thus incorporated through the image series, while the full time-domain polarization evolution $E(t)$ is not implemented [2507.22225].

A practical implication is that pAvalancheMC evolved from a solver for $E_a+E_s$ to a more general RPC field engine in which $E_{\mathrm{tot}}$ can include dielectric-image contributions. This suggests that later usage is not limited to gas-phase self-fields alone, but also addresses the layered electrostatics of realistic RPC stacks.

## 5. Signal induction, validation, and numerical performance

Induced signals are computed through the Ramo–Shockley formalism using weighting fields from neBEM. The induced charge is written as
$$
Q_{\mathrm{ind}} = \int_0^t dt \sum_i q_i \left[\phi_w(r_{f,i}(t))-\phi_w(r_{i,i}(t))\right],
$$
and the induced current, when requested, as
$$
i(t) = \sum_i q_i\, v_i(t)\cdot E_w(r_i(t)).
$$
The thesis also gives an analytic approximation for the weighting field in a single-gap three-layer RPC geometry,
$$
E_z^{\mathrm{weight}} =
\frac{v_0 \epsilon_1 \epsilon_3}
{\epsilon_2 \epsilon_3 d_1 + \epsilon_1 \epsilon_3 d_2 + \epsilon_1 \epsilon_2 d_3}.
$$
These signal calculations are integral to timing RPC studies rather than an external post-processing stage [2211.06361] [2507.22225].

A timing RPC example with a $0.3$ mm gap at $43$ kV/cm reports induced-charge distributions at $1720$ V, $1730$ V, and $1735$ V fitted with a Polya function,
$$
f(Q)= a \left(\frac{Q b}{c}\right)^{b-1} \exp\!\left(-\frac{b}{c}Q\right).
$$
The mean charge parameter increases from approximately $0.34$ fC to $0.80$ fC to $1.3$ fC as the voltage increases, and the parameter $b$ shifts from $3.8$ to $4.9$ to $8.2$. The low-charge inefficiency peak was excluded from the fits, and the trend is reported as consistent with experimental observations attributed to Fonte and Lippmann [2211.06361].

The space-charge dynamics were validated against earlier solvers and against experimental trends. The 3D Cartesian implementation reproduces features reported in earlier 2D rotationally symmetric models, including radial blow-up and contraction and comparable gain dynamics. The thesis additionally states that the segmented-line field matches the ring model after a $2\pi$ scaling and that source-plus-image fields agree with neBEM near the electrode when image charges are included, with ratios approximately unity there and about $0.8$ to $1.25$ elsewhere where the field is not vanishing [2507.22225].

Parallel performance is a defining feature. In a timing RPC of area $30\times 30$ cm$^2$, gas gap $0.3$ mm, electrode thickness $2$ mm, and field $43$ kV/cm, the speed-up $v_p = T_1/T_N$ for $10^4$ avalanches was fitted as
$$
f(N) = p_0 - p_1 e^{-(p_2 N^{p_3})},
$$
with asymptotic speed-ups $p_0 \approx 5.46$ without space charge and $\approx 7.2$ with space charge. The saturation of speed-up is attributed to OpenMP overheads and the serial fraction of the code. neBEM parallelization was tested up to $24$ cores while preserving precision and significantly reducing field computation times [2211.06361].

## 6. Assumptions, limitations, and research outlook

The modeling assumptions are explicit. Space charge is represented by straight line segments on a cylindrical grid; the field is treated as uniform within each voxel; near-field singularities are moderated by finite-segment integrals and voxel averaging; and positive ions and negative ions are static unless ion mobility is supplied. This last simplification is justified in the sources by the large disparity between electron drift velocity and ion mobility on the simulated time scales [2211.06361].

The code does not implement full streamer physics. It can identify streamer-like onset through rapid gain divergence and high local fields, but photon transport and photoionization are absent. This is important because saturated avalanche and streamer onset are treated differently: the former corresponds to field suppression such that $\alpha-\eta<0$ over most of the avalanche volume, whereas the latter is inferred when electron number and $E_t^{\max}$ rise sharply without a clear saturation plateau. In C$_2$H$_2$F$_4$/i-C$_4$H$_{10}$/SF$_6$ at $E_a=50$ kV/cm, the reported regime is characterized by a peak followed by a sudden increase in electron number while $E_t^{\max}$ exceeds approximately $90$–$100$ kV/cm; full streamer simulation is deferred to future work [2211.06361].

The thesis adds further caveats: ion–ion interactions, detailed electrode surface roughness, graphite resistivity nonuniformity, and a complete transient polarization model are simplified or absent. Sensitivity to electrode relative permittivity $\epsilon_r$, thickness $d$, and gas gap $g$ is emphasized: the field increases rapidly with $\epsilon_r$ and saturates beyond about $20$ for fixed $d$ and $g$, decreases with increasing $d$ for fixed $\epsilon_r$ and $g$, and decreases with increasing $g$ for fixed $\epsilon_r$ and $d$ [2507.22225].

Within these bounds, pAvalancheMC functions as a detector-physics tool for RPC optimization rather than a universal plasma-discharge simulator. Its established strengths are self-consistent avalanche growth under dynamically updated fields, induced-signal prediction using realistic weighting fields, and computational tractability through OpenMP parallelism. The stated research direction is toward photon transport, more detailed time-dependent electrode charging, richer diffusion and field solvers, and broader acceleration strategies such as GPU or MPI-based parallelization [2507.22225].

Source: https://www.emergentmind.com/topics/pavalanchemc