COSICAF, a fission chamber simulation tool for academic purposes

. Nuclear instrumentation is a complex topic since it involves a wide range of physics phenomena like nuclear reactions, heavy ion interactions with matter, electrostatic, charge creation etc. Understanding and modelling fission chambers is a difficult task usually performed with Monte-Carlo and finite element simulations. Since a few years, analytical and simplified Monte Carlo models were introduced at the French Atomic Energy Commission to easily design detectors. It is proposed here to present the derivation of such model, called COSICAF, for academic purposes; this numerical model provided with this article, will help students and researchers to understand and design fission chambers. To demonstrate the interest and the limitation of proposed work in research field, the model is applied to simulate two real miniature fission chamber designs.


Introduction
Since the 1940s nuclear instrumentation plays a key role for industrial and scientific operation of nuclear reactors.Ionization chambers and more specifically fission chambers are widely used in power plant reactors for safety and monitoring purposes.In nuclear science, such devices allow to monitor experiments in a very precise way.It is possible to estimate the nuclear power in real time or to probe the neutron flux either in energy or space.To master the measurement uncertainties and limit the perturbation induced by detectors, it is essential to finely tailor detectors and electronic measurement systems to the needs.
Fission chamber is a special kind of ionization chamber.It is constituted of a gas tight body containing at least one fissile deposit, two polarised electrodes and a volume of noble gas.When a neutron passes across the detector, there is a chance of fission inside the deposit.Fission fragments can interact with the gas and charged particles are generated.By collecting the charged particles thanks to polarised electrodes, current pulses are generated at the output of the fission chamber.
Depending on the neutron flux and fission chamber design, such detector may work in three different modes: -Pulse mode: each pulse is distinguishable and pulse pile up is negligible.
* e-mail: gregoire.deizarra@cea.fr-Campbell mode: pulses may pile up or not being distinguishable any-more.It is still possible to relate the count rate to the variance (or higher moment) of the signal.-Current mode: the count rate is high enough to produce a current proportional to the count rate.Each mode has advantages, pulse mode can be useful for spectra ratio in zero power reactor since count rate is easily related to reaction rate while Campbell mode is useful to follow transient power pulses because it makes possible flux monitoring over 6 to 7 decades.Depending of the application, it is necessary to properly scale the detector in order to use a specific measurement mode and the associated electronics.The working mode is not the only parameter to take into account for reactor applications, sensitivity and size could also constraint the design.Fission chambers can operate under different regimes depending on the polarisation voltage as shown in Figure 1.At low polarisation voltage, charged particles tend to recombine each other and the collected charge is lower than what was created in the gas volume.At higher voltage, the U/I curve reaches a plateau where the collected charge correspond to the charge created along the track of the fission fragment, there is no charge multiplication: it is the saturation regime.Fission chambers usually operate in this regime.If the polarisation voltage is increased, Townsend avalanches are initiated in the gas allowing proportional charge amplification.Further increase of the voltage lead to non proportional charge amplification and to Geiger regime.Even if the regular fission chamber designs are now well mastered, simulation of miniature and innovative fission chambers is still of interest for the scientific community [1][2][3].
In this paper, COSICAF, a simplified open source fission chamber model is proposed mainly for academic purpose.The source code is available under a GPL licence with the article.It uses a Monte-Carlo method to simulate independent current pulses at the output of a detector in saturation regime.Despite some strong assumptions, this model is still valuable for a qualitative comprehension of the detector operation as well as rather good quantitative evaluations.Generated pulses have many applications; they can be used to taylor a detector for specific uses or electronic instrumentations but also to predict the signal of exotic configurations.Pulse charge distributions and pulse shapes might help to calibrate the detector in Campbell mode.For experiments, they make possible the tuning of electronic timing or the correction of detector efficiency (number of pulse detected per fission in the deposit).
In the first part of this paper, the phenomenology of a fission chamber is described and simple models are derived.In a second part, models are implemented in a Monte-Carlo simulation, including specific discussions on source computation and geometry effects, for instance.
To show the interest of such model for research purpose, two case studies are proposed with detectors described in the literature: a miniature fission chamber and a micro pocket fission cell detector.
The source code of the simulation is provided with the paper in order to be used by students.It is written in Octave/Matlab [4] to easily manipulate data and produce plots at the cost of low computational performance.

Modelling fission chambers
Numerous ionization or fission chamber models are available in the literature.Knoll [5] derived results about radiations detector physics, operation mode and measurement electronics, which allow to understand all the detection processes.More rencently, semi-analytical models were derived by Poujade and Lebrun [6] and Chabod Fig. 2. Schematic fission chamber.[7] to evaluate the saturation regime of fission chamber.A rather complete semi-analytical model was provided by Chabod et al. [8] to simulate fission chamber in current mode.Such model takes into account the effect of charge multiplication and recombination to describe the different operating regime of the detector.The effect of space charge is also taken into account.Filliatre et al. [9] proposed a Monte-Carlo model of fission chamber, called CHESTER, to simulate the induced current pulses.It relies on GARFIELD a computer suite for drift chamber simulation, and uses SRIM for interaction of fission fragments with matter in conjunction with MAGBOLTZ for the electron transport.Only the interaction of fission fragments with gas is assumed: no transport in solid materials like fissile deposit is considered.Another simulation framework, PYFC was developed by Elter [10].It uses TRIM to compute fission fragment trajectories in gas and a python program to project them randomly in a fission chamber geometry.Once again, only the interaction of fission fragment with gas is assumed.
The model proposed in this work is different from the ones previously mentioned, it involves simple semiempirical relations for most of the physical processes.Each fission fragment is followed individually in order to compute the related pulse at the output of the detector.Unlike models proposed by Filliatre and Elter, arbitrary geometries can be simulated and several materials like gases, fissile deposits, and fissile deposit's protective layers can be taken into account.The possible geometric effects are not discussed here but are part of the discussion about model implementation in a numerical simulation (Sect.3).

The big picture
A schematic fission chamber is presented in Figure 2. It consists of a gas tight body containing at least one fissile layer, two polarised electrodes and a volume of gas, usually argon.When the detector is positioned in a neutron flux, there is a probability to induce fissions inside the fissile deposit: some neutrons and two fragments of the fissioned nucleus are generated for each fission with a large amount of kinetic energy.The fissile deposit is tailored to make possible the exit of at least one fission fragment.Then, the fission fragment is slowed down in the gas through successive collisions with atoms: excited atoms, ions and electrons are produced all along the trajectory of the heavy ion.At last, thanks to the electric field in the detector, electrons and ions start to migrate toward the electrodes and induce a current pulse at the output of the detector.
Two main fission chamber designs are commonly used, the flat and the cylindrical configuration as depicted in Figure 3.The inter-electrode gap is usually of the order of a millimeter and the electric field is oriented as depicted on the figure except for specific cases.For higher sensitivity (in pulse count per n/cm 2 /s), some fission chambers have more than one set of electrodes but the configuration is still flat or cylindrical.
The main objective of fission chamber modelling is to estimate the signal at the output of the detector by taking into account all the physical processes described previously.
For the derivation of the model, we will consider only a fission chamber with uranium 235 as the active material and a thermal neutron flux.However, the approach is generic enough to use the same model with different fissile materials and neutron spectra with only small changes.Only the saturation regime is considered, nevertheless, the proportional regime can be simulated through the use of Townsend first ionisation coefficient and the modification of few COSICAF's functions.

Fissions process
A fission chamber is set in an energy and position dependent neutron flux φ(r, E).As a first approximation, the neutron flux is not perturbed by the detector.The number of fissions per second occurring in the fissile deposit can be estimated with: where S is the volumetric fission rate (cm −3 .s−1 ), N f is the concentration of fissionable atoms (cm −3 ), and σ(E) is the microscopic fission cross section of the considered material (in barn or 10 −24 cm 2 ).If we assume a uniform thermal flux, equation ( 1) can be simplified by defining the Westcott conventional flux: with n the total neutron density in (cm −3 ) and v 0 the neutron most probable velocity at 20 with σ 0 the U235 fission cross section at thermal energy, 558.9 b.It is clear that the sensitivity (in count/n/cm 2 /s) of a detector is directly related to the volume and to the kind of fissile deposit used.The fissile mass allowing an optimum sensitivity is constrained by the volume of the detector, by the thickness of the deposit which limits the available energy of the fission fragments to ionize the gas and by the perturbation of the neutron flux.If a fissile deposit is too thick, fission fragments might not be able to reach the gaseous medium or their energy is too low to produce a measurable signal at the detector output.Two effects may perturb the neutron flux, the first one is called flux depression [11], it corresponds to a flux lowering around the detector due to the introduction of a "neutron sink".The second one is called self-shielding or self protection: the centre of the deposit is shielded by its surrounding and the viewed flux is lowered.Jammes et al. [12] showed this last is negligible with a global perturbation of the neutron flux below 1% for a 2 mg/cm 2 U 3 O 8 deposit.
Since pulses are simulated independently and it is assumed that self shielding is negligible, reaction rate computation is uncoupled from signal generation.In reality, the reaction rate has an impact on the signal generation, but with this assumption, the model does not permit to describe the generation of space charge in the interelectrode space, increasing with the flux level.Thereby, COSICAF produces signal without dealing with neutron interaction with matter.
The fission process is complex to model, fortunately, for detector simulation, we are only interested in the fission yield and in the energy released during the fission process.The energy released during a fission is ).For 235 U, an average of 2.4 neutrons are released as well as two fission fragments (f.f.) with an average total kinetic energy of 169.5 MeV [13,14].The fission yield for 235 U is depicted in Figure 4; whatever the fissile material considered, the fission yield is similar.Two groups of fission fragments are noticeable, heavy ones with a most probable mass of 138.5 uma and light ones with a most probable mass of 95 uma.For the sake of simplicity, the distribution is reduced to two mean fission fragments only.It simplifies the simulation of fission events and decreases drastically the amount of data needed in the model.This approximation is quite common in nuclear instrumentation.Chung and Prelas [15] and Mat'ev [16] showed that this approximation only has a minor impact on energy losses in gases when a large quantity of heavy ion is generated per second.For regular in-core and ex-core fission chambers where only a fraction of the f.f.kinetic energy is lost in the gas, this assumption has a minor impact.For other cases, where all the energy is lost in the gas, only the mean pulses are simulated: on the charge spectrum, the geometric effect is noticeable but no energy spread is related to fission yield.The kinetic energy distribution of mean fission fragments is derived with: where E lff and E hff are respectively the kinetic energy of light and heavy fission fragment.m hff ,v hff are the mass and the velocity of the heavy f.f.Table 1 shows the average fission fragments for different fissile materials and their kinetic energy.
To summarise, the fission process in the detector is modelled in the following way: -The fission rate might be estimated by either using equation ( 1) and a Monte-Carlo neutron transport code or from equation ( 3).-No self shielding is considered.
-Reaction rate estimation is uncoupled from current pulses generation.-Only two mean fission fragments are considered after a fission, their characteristics are listed in Table 1.
Once heavy ions are created, their interactions with detector materials remain to be evaluated.

Heavy ion transport
When heavy ions (like f.f.) are travelling in a medium, collisions with target atoms slow them down.At a microscopic level, interactions are numerous: ions are able to capture bound or free electrons, to lose electrons, to become excited,to excite or ionize the target and, at last, to make elastic collisions.All the cross sections related to those processes are not available or precise enough to describe precisely ion trajectories [17].
Heavy ions slowing down is usually described using a mesoscopic quantity called stopping power which is defined as a friction force.According to Bohr [18], it is composed of two different parts.The first one is called the electronic stopping power, it describes all the processes occurring between the projectile ion and the electrons of the target medium.It takes into account mainly ionisation, excitation and charge transfer processes.The second part is the nuclear stopping power, which is related to elastic collisions.
The state of the art of heavy ion interaction with matter simulation is the SRIM simulation suite [19].It contains semi-empirical stopping powers, range computation, and a module, TRIM, that computes ions trajectories with a Monte-Carlo method.TRIM considers that collisions with electrons lead to continuous energy loss and no trajectory change, while some collisions with atoms are susceptible to change the energy discontinuously and the trajectory of the ion.The module uses the Magic algorithm derived by Biersack to speed-up the collision integral evaluation with a rather good precision [20,21].
To simplify ion transport, we assumed that heavy ions have straight trajectories in the deposit and in gas.Regarding dimensions of a fissile layer thickness (≈1µm), of an inter-electrode gap (≈1 mm) and heavy ion initial energy, this approximation is acceptable.Figure 5 shows the stopping power for a 235 U light f.f.As it can be seen, the nuclear stopping power is preponderant only at low energy.Depending on the level of approximation needed, various derivations of the ion transport could be considered.The most simple one consists in integrating the stopping power over the ion rectilinear path to link energy losses to the travelled distance: x where S(E) is the total stopping power in MeV/m or, to be independent of the material volumic mass, in MeV/(mg.cm 2 ).This implies a straight path of the ions i.e. a shorter travelled distance than the real one.A more accurate description consists in taking into account the mean effect of heavy ion deflection and computing the mean projected path along the direction of propagation (also called range).Due to the deflection, ion tend to lose more energy per unit distance than with equation (6).To understand the process without deep knowledge in physic, the approach of Biersack [22] is interesting.It solves the transport by considering the diffusion of heavy ion to get: where S(E) is the total stopping power, S n the nuclear stopping power, E the projectile energy and µ = M t /M p the ratio of the target and projectile particle mass.This equation is a part of PRAL [23,24] (Projected Range ALgorithm), which is the core of SRIM software.PRAL includes as well an equation defining the straggling, a measurement of ions' trajectories spread.Compared to equation ( 6), an additional term slows the ions.Results of equations ( 6) and ( 8) for light fission fragments slowing down in argon are available in Figure 6.A maximum discrepancy of 5% is noticeable when no deflection is considered for range computation.Regarding the validity of approximating the ion trajectory by a straight line, we estimated the lateral straggling in the case of 100.5 MeV light f.f. in argon, which turns out to be 4.5% of the range.In a solid, heavier material like uranium oxide, the lateral straggling is around 10% of the range.The straight trajectory hypothesis is thus rather valid.
To make problems more tractable, Nguyen and Grossman [25], among others, obtained an analytical moderation law by solving equation ( 8) with Bohr's electronic stopping power: where E(x) is the energy of the particle after it travelled a distance x, R E0 is the range of the particle with an initial energy E 0 and n a semi-empirical parameter linked to the shape of the stopping power [26].Such moderation laws are used for example in nuclear lasers to compute quantities like the deposited energy in the lasing medium [15,16].
To describe ion slowing down in COSICAF, we use moderation laws.As they are used in a numerical model, it is possible to replace the function of equation ( 9) by a polynomial fitted on f (1 − x/R E0 ) = E/E 0 data: with a, b, c, d, e the fit parameters.COSICAF includes the parameters and the implementation of equations ( 9) and ( 10) for f.f. in all the noble gases and also in solid U 3 O 8 .To expand the database, PRAL with SRIM or ASTAR stopping power libraries can be used to compute the range as a function of energy.Then, by computing the reduced travelled distance x/R E0 as a function of E/E 0 , it is possible to fit equation (10) and to gather new data for COSICAF.Considering the total stopping power brings some errors in the models since the main source of electrons is related to the electronic stopping power: if an uranium 235 light f.f. is completely slowed down in argon, 96.5 MeV is spent in electronic collisions and 4 MeV in nuclear collisions.Not taking into account the different stopping powers leads to an anomalously high amount of electron/ion pairs created at the end of the ionisation track.
To correct this effect, a function was fitted on electronic energy losses: where E e is the remaining electronic energy to spend, E 0 the initial energy of the projectile and E the total kinetic energy of the particle.With this correction, the error on the electronic energy losses is only of the order of 0.5% of the total energy.Deviations are mainly located at low energy, below 1 MeV, where nuclear stopping power becomes preponderant (see Fig. 5).In summary, the heavy ion transport is simulated in COSICAF using the mean fission fragments with the following assumptions: -Projectiles have straight trajectories.
-Moderation laws (Eqs.( 9) or (10)) are used to link the travelled distance to the energy losses.-Moderation laws parameters were obtained through fit on data computed with PRAL algorithm and SRIM stopping powers.-A correction is applied to only consider electronic energy losses i.e. ionising collisions.

Charge pair creation
By studying the slowing down of fission fragments in matter, the energy transferred to the target medium can be derived as in previous section but it gives no information on the amount of charges created during this process.In the early 20th century, the study of ionizing particles interacting with gases led to the estimation of W , a mean energy quantity that the projectile has to lose to create where E 0 is the initial particle energy and n e the total number of electrons created along the particle track.This quantity was first derived experimentally before theoretical work [27] allowed to compute it.Many ionizing particles were studied, from electrons to alpha particles and fission fragments [28].Whatever the particle, the measured values of W in a medium are similar.From a classical mechanics point of view, it is expected that electron ionising collisions are more effective than fission fragment ionising collisions since the mass ratio between bounded electrons and heavy ion is not favourable to energy transfer.In practice, this is not the case which means that electron/ion pair creation is mostly driven by secondary electrons generated along the ionization track and by the electronic shells characteristics of the target medium.The experimental evaluations of W for heavy ions give a lower value than expected because of the socalled "ionisation defect" due to energy losses by elastic collision.
In COSICAF, we chose to include W value for heavy ions found in the literature [29], results for different noble gases are available in Table 2.

Charge transport
An electric field is applied in the detector to collect the charges generated during the fission fragment slowing down.Electrons and ions are moving respectively toward the anode and the cathode.Because of the surrounding neutral atoms, charged particles make collisions and scatter, preventing them to have a ballistic transport.For a given electric field, a group velocity resulting from an equilibrium between acceleration and scattering due to collisions can be defined by: where v p is the group velocity of charged particle p, it corresponds to the mean velocity of a charged particles swarm.E is the electric field.µ p is the mobility in (m 2 .V −1 .s−1 ), it is dependent of the medium density since it describe the trade off between random scattering with surrounding atoms and acceleration in the electric field direction.It can be viewed as the property of a medium, to reach a given charged particle velocity by electric field unit.For fission chamber manufacturing, the quantity of interest for optimizing charge transport is not the density but the filling gas pressure.Since the gas filling is performed at a known and constant temperature, it is proportional to the medium density.
The mobility of electrons were measured and computed for numerous gas mixtures.First measurement of electron swarms drift velocities are mentioned in the 1920s by Townsend [30] with a drift chamber containing an electron emission device and two grids to measure the drift time.Nowadays, electrons mobility and other transport parameters can be obtained from numerical simulation.Two reference codes are used while simulating fisson chambers: MAGBOLTZ and BOLSIG [31].Both of them are solving the Boltzmann equation with a two term approximation for the electrons.In COSICAF, we use electron drift velocities measured by Haddad for argon-nitrogen mixtures and tabulated by Raju [32] for other noble gases.As for heavy ion transport, functionals were fitted on data to provide a compact expression of the drift velocity that can be easily used in the calculation.Deviation of the fit is generally below 2%.
For ion drift velocities, we used Frost's semi-empirical formula [33]: where N is the gas particle density (m −3 ) which is proportional to filling gas pressure, E/N is the reduced electric field in (Td = 10 −21 V.m 2 ), a and b are two fitting parameters.Those parameters, which are included in COSICAF, depend of the gas temperature and were derived by Golyatina and Maiorov [34].
Charges diffusion also occurs during the drift, however this process can be neglected regarding its order of magnitude in a fission chamber working in saturation regime.Charge recombination is neglected as well, hence only the simulation of saturation regime is possible.

Induced signal
The signal from a fission chamber is not due to the collection of electric charges on the electrodes but to the current induced by the motion of charges within the detector body.Prior to the 1940s, the computation of this induced current was done by evaluating the rate of change of the electrostatic flux at the electrode: However, computation methods involving image charges are tedious and difficult to apply for complex geometries.To estimate the induced current, Shockley and Ramo have derived a more efficient method at the end of the 1930s [35,36].
The Shockley-Ramo theorem states that the current induced by a point charge on an electrode is given by: where v is the charge velocity and E 0 is the weighting field (see definition below).The velocity is computed by taking into account the electric field generated by the electrode and the possible space charges while the weighting field is the field which would exist without space charge and with all the electrodes are at zero potential, except the electrode of interest which has a unit potential [37].The weighting field is defined by: where ψ 0 , the weighting potential, is dimensionless.Setting it to unit allows to see the changes of the induced charge as a fraction of the moving charge.It can be computed by solving the Poisson equation: with ψ ei the boundary condition at the electrode of interest and ψ e the boundary conditions at the other electrodes.For simple geometries, solving the equation is straightforward.However, for complex ones, numerical methods like finite difference scheme and matrix inversion or relaxation methods should be used.To estimate the signal, COSICAF implements the Shockley Ramo theorem, this is discussed in detail in the next section.

COSICAF: numerical model
Models described in the last section are not sufficient to compute fission chamber signal.Practical problems like handling geometry, source definition, particle tracking must be solved as well to simulate the signal.
To establish the foundations of a Monte-Carlo simulation, the geometric description of the simulation is of prime importance because of its impact on particle tracking.Many Monte-Carlo codess like TRIPOLI or MCNP use a Constructive Solid Geometry to define a problem.It consists in defining simple shapes from surface equations.More complex objects are constructed by combining shapes thanks to boolean operations.
Since COSICAF is dedicated to education, such a geometry description system cannot be implemented because of its complexity.
When studying existing detectors, two main configurations are noticeable: planar and cylindrical.Both are made of discs and cylinder shells.It was then decided to reduce the geometry definition to cylinder and cylinder shells parallel to the z axis.All the pieces of a geometry are then defined by an outer and an inner radius plus a cylinder height.The position of the center of the cylinder base allows to place the different parts of the geometry in the space.Neither boolean operations nor intersecting volumes are permitted to keep particle tracking simple.Nevertheless, the inclusion relation must be taken into account to put electrodes or dielectric inside the detector body.The inclusion relation is managed with an inclusion index: volumes with lower indexes contain the ones with higher indexes.Some parts of the detector like electrodes or fissile deposit are processed in a specific way by the simulation code and should be spotted during geometry definition.This is done by assigning the following keywords to each volume of the simulation: -bound: define the simulation space; -regular: define a regular volume; -stop: define a volume which stops every heavy ion; -no coll: define a volume which doesn't interact with heavy ions.Sources are defined with the help of a flag plus an heavy ion definition as shown in Figure 7.
Such description still allows to build complex detectors and simplify all the simulation steps.

Source definition
To simulate the fission process, a source of specific particles is linked to a volume.Sources are not related to a real reaction rate, they are normalized.All events are considered independent: the model is not able to capture the generation of space charge induced by high reaction rate.The emission of heavy ions is modelled as uniform in a volume.The emission positions in a cylinder shell are sampled with a uniform random number generator as follows: - where rand() is a call to the random number generator, R out , R in and h are respectively the outer radius, the inner radius and the height of the cylindrical source.x c , y c , z c is the position of the cylinder base center.Since only single particles are generated, it is not possible to emit two f.f. at a fission location.This case can be taken into account at a price of a slightly more complicated source simulation.
The initial propagation directions of heavy ions are selected randomly over all the possible directions by using the following algorithm: -ang1= 2π * rand() -ang2 = acos(1 − 2.rand()) e x = cos(ang1).sin(ang2)e y = sin(ang1).sin(ang2)e z = cos(ang2) where, e x , e y , e z are the components of the propagation direction vector.Initial energy of the particle is defined in the source; because of the semi-empirical laws used, particle energy must be lower or equal to E 0 , the maximum energy of the projectile used during computation of moderation laws.The energy definition is useful to simulate the emission of α particle during various reactions for example.Once heavy ions are generated with an initial direction, particle tracking all along their trajectories is done.

Fission fragment transport and charge creation
The heavy ion transport is done assuming a rectilinear trajectory.It consists in calculating the entry and exit points of all the volumes that the particle passes through until it stops.At these remarkable points, the kinetic energy of the particle is calculated.It is thus possible to know where it stops by examining its escape from the geometry or its kinetic energy.Once the transport is calculated, generation of electron/ion pairs is performed along the ion track.
The calculation is carried out as follows.In the initial state, the volume containing the particle, its current position r 0 and its trajectory are known since they have been previously calculated.The exit point of the current volume is searched for.Because of the inclusion relation, the exit point of the ion is determined in two steps.First the intersection between the trajectory and the cylinders plus the two planes representing the current volume is calculated.The intersection of the trajectory in the direction of propagation with a plane parallel to xy is given by: with z p and z 0 respectively the location of the plane and the altitude of the particle located in r o .k must be positive to point toward the direction of propagation.The intersection is: The intersection of the trajectory with a cylinder oriented toward z direction is given by a quadratic equation: where x x and y c is the location of the center of the cylinder, x 0 and y 0 the particle coordinate, R the radius of the cylinder.
Second, intersections of the trajectory with the inner volumes are estimated.From all those possible exit points, the correct one is the closest one to r 0 .
To summarize, the following algorithm is used to find the exit point of the current volume: -Compute the intersection points of the particle trajectory with the surfaces constituting the volume containing the particle.-Check if inner volumes are included in the current volume.-For each inner volume, compute the intersection with the trajectory.Check if intersection points are belonging to inner volume or not.If not discard intersections.-Sort the intersection points found regarding the distance to the actual particle position r 0 .-Select the closest intersection point.
Exit point is stored in memory and the next volume visited by the heavy ion is found by moving toward the current volume exit point from a 10 −9 m distance and by checking in which volume this position is by using the following algorithm: -Sort volume by decreasing inclusion index.
-For each sorted volume, check if the position is included.
If yes, then break.The entered volume is set as the current volume and r 0 is updated: the new position of the particle corresponds to the entry point of the current volume.Now, a new exit point has to be computed and so on.
By keeping track of the heavy ion energy and position at entrance and exit of each volume, we propagate the particle until its energy becomes null or until the geometry boundary is reached (Fig. 8).
Charged particles are then produced along the ion track in gaseous media only.As W is orders of magnitude lower than the heavy ion energy losses, following each electron/ion pair is not feasible regarding computation time and model limitations.Instead, a fixed number of charged meta-particles which represent hundreds of ion/electron is considered.A trajectory step is computed by dividing the trajectory length by the number of meta-particles defined by the user.On each trajectory step, the energy loss is computed using moderation laws.The correction for electronic stopping power is applied and, using W values, two charged meta-particles are generated.

Charge transport and induced signal
Once charged particles are generated, they are moved according to the electric field computed at their location and to the drift velocities.The electric field needed to estimate the particles velocity and the induced signal can be computed analytically or numerically.For planar and cylindrical geometry, the electric field is respectively given by: with e V the direction of electric field, L the inter-electrode distance, R o and R i respectively the inner and the outer radius of the electrodes and r the position.COSICAF allows to define a user function which returns E, the electric field, for such a simple case.Users have to implement equation ( 22) or ( 23) by their own.User's function is then called automatically during the charge transport simulation.
To deal with more complex geometries, solving the Poisson equation numerically is mandatory.For this purpose, a simple algorithm implemented in COSICAF is proposed in next sub-section.
At each time step, meta-particles are moved by a ∆l step defined by: where ∆t is a time step defined by the user and v d the drift velocity of the meta-particle.Each charge movement induces a current on electrodes, it is given, on the k electrode by: with E 0k the weighting field associated to electrode k and q m the charge of the meta-particle and E the norm of the e.Meta-particles do not contribute to the signal any more when they encounter a solid material and therefore their tracking is stopped to speed-up computation.

Computing the electric field in a generalised way
In order to compute the electric field numerically, the geometry has to be discretised.For the sake of simplicity, we chose to use an equispaced grid.Each volume in the problem is related to voltage boundary conditions and/or to a relative dielectric permittivity.They are converted to a discrete geometry by using a Bresenham circle algorithm and a flood and fill algorithm.The first algorithm is used to define the boundaries of the volume while the second allows to set the permittivity to a proper value inside the boundaries.Once the geometry is discretised, solving the Poisson equation can be done with a finite difference scheme or by random walk.Finite difference scheme and matrix inversion or relaxation method work well but necessitate a reasonable amount of work for an implementation in 3D.As computing electric field is not the main purpose of the simulation, it was chosen to use a Monte-Carlo resolution.Kakutani [38] has shown that it is possible to solve the Poisson equation with Dirichlet boundary conditions using random walk.Among the few existing methods [39,40], we chose the fixed random walk one [41].Let us consider the following equations defining the electric potential on a three-dimensional domain G: is the dielectric permittivity and V is the electric potential.It is linked to the electric field E by: The domain G is divided into a mesh with a constant step h.To estimate a potential at r, a particle is set at this position and starts a random walk on the mesh until it reaches a boundary, the associated voltage V p is then kept in memory.This process is repeated N times, then the solution of the Dirichlet problem in r is: Results converge as √ N , a reasonable amount of random walks is needed to ensure accurate voltage and noiseless electric field.Since the electric potential on all the mesh is required, the computation can be speeded up by updating all the potentials along the random walk.
To select the direction of walk, a random number a, with 0 < a < 1 is generated at each step from a uniform distribution.The probability to move in one (+z) of the 6 directions is given by: The interval [0-1] is subdivided according to those probabilities, hence, a is used to select the direction of the random walk.An illustration of the algorithm is given in 2 dimensions in Figure 9.
In COSICAF, solving electrostatic problem implies the definition of an additional file containing all the boundary conditions and the dielectric permittivities.Random walk solver is implemented either in Octave and in C++ to provide routines with reasonable performances.

Summary: COSICAF overview
As it was shown through this section, COSICAF gives access to advanced simulation features like arbitrary geometry definition, particle tracking, electrostatic computation.Since none of those modules are independent, students should use them for case studies or for a tutorial about how simulation code works.
The software operation is summarized in Figure 10.Simulation is reduced to scripts in which every physical process is well separated from the others, making code changes or experimentation easy for students.To allow the use of the software in practical work of nuclear instrumentation and computational physics, all the basic data and models presented in Section 2 are made available to users.Hence, simple case studies and tests can be built by students.

Code validation: miniature cylindrical fission chamber
Miniature cylindrical fission chambers are the perfect test case for COSICAF validation.The geometry is simple enough to estimate energy losses in the gas and the amount of electron/ion pairs generated.In addition, the electric field is analytical, and pulse duration can be computed independently from the simulation.For code validation, we chose to simulate a CFUR fission chamber, a case already treated with Pyfc and Chester.It consists of a cylindrical body with cathode and anode of respectively 1.25 mm and 1.00 mm radius.The deposit is located on the anode and has a 14 mm length, it usually consists in micrograms of uranium oxide.The detector is filled with a Ar + 4% N 2 gas mixture at a pressure of 5 bar.Some electronic and ionic pulses were first computed to provide a sanity check of the simulation code.The deposited energy along light f.f trajectories in gas was compared to TRIM simulations (Fig. 11, top).As it can be seen, the energy deposition evaluated by COSICAF is close to the one estimated by TRIM except for really long path were the rectilinear trajectory hypothesis does not hold any more.In saturation regime, it is shown that the total electron charge generated in the inter-electrode space is equal to the total charge induced at the electrode.To ensure the correct operation of COSICAF, we have compared for several pulses (Fig. 11, bottom) the induced electron charge (black bars) and ion charge (yellow bars) to the electron charge generated in the inter-electrode space (green bars).For each pulse, the induced charge was estimated by integrating the current pulse.The induced current is estimated in a satisfactory way by the software in regard of the relative error of 0.1% between the initial electron charge and the sum of the induced charges.Some electronic and ionic pulses computed with COSI-CAF are shown in Figure 12.Charge collection times are similar to those estimated with other simulation codes.The shape and amplitude of single pulses are also consistent.Current is maximum at t = 0 s because fission fragment transport and charge collection are uncoupled: fission fragment transport is then considered as instantaneous.However, a small number of pulses with anomalously high amplitude are noticeable.These spurious pulses correspond to extra long fission tracks.Such events would have a far smaller probability if straggling was not neglected in the model.With a careful data analysis, it  is possible to detect these events and assess their effect on the mean pulse shape.The charge distribution of electronic pulses is given in Figure 13.Most probable induced charges are comparable to the one obtained with PYFC [10] and an accurate fission yield.
For simple detectors like cylindrical fission chambers, COSICAF is able to simulate the output signal with rather good characteristics regarding the low amount of data and the simplifying hypothesis used in the model.University and INL [42].They are designed for highneutron-flux operation and are small enough to be assembled in stacks.Many designs of micro-pocket detector were proposed during the last ten years [43].We have chosen to simulate one of the last design iteration [3] because of its complexity and because of the complete description of its dimensions and materials [44].

Micro-pocket fission detectors
A scheme of the detector is available in Figure 14.consists of an alumina casing of 0.47 cm diameter where an ionisation chamber and a fissile deposit are enclosed.The ionisation chamber has a height of 0.2 cm and a 0.2 cm diameter.The alumina dielectric contains 8 holes for electric wires and a window to expose two wires per ionisation chamber.The 8 holes allow to stack up to 4 different mpfds such as the one depicted in Figure 14 at a same location, in a gas tight tube.The filling gas is argon but there is no clear mention of the pressure or gas density, from the estimation of fission fragment energy deposit in gas and from other publications [43,45], it can be assumed that the filling gas pressure is 1 atm.The fissile deposit is assumed to be uranium oxide with a thickness of 0.030 µm.To describe the geometry with only disks and cylinder, a few approximations were done: alumina is decomposed in large cylindrical pieces below and above the chamber.At chamber level, the surrounding dielectric is constituted by a series of small radius cylinders (Fig. 15).

Eletrostatic computation
Depending on mpfd designs, different polarisation voltages were tested [46].We chose to compute the electric field generated with a 200 V voltage to be in the saturation regime.This value was selected in order to get a reduced electric field similar (≈4 Td) to the one of regular fission chambers working in saturation mode.FEMM software [47], which solves electrostatic problems through a finite element method, was used to get a reference 2D electric field map in the mid plane of the mpfd.Then, the  COSICAF built-in Monte-Carlo solver was used on an approximated 3D geometry to produce the electric potential output presented in Figure 16.Potential in the mid plane of the ionization chamber computed with the two codes are in agreement within a few percent which is satisfactory regarding the approximations.All the cathodes and all the anodes are grouped on a different half of a detector stack; because of symmetry, only two configurations were simulated in order to represent all the polarisation cases that can appear in a mpfd stack.Whatever the mpfd considered, only a cathode and its corresponding anode are exposed.Case A corresponds to the case where an electrode in the window is surrounded by electrodes with the same polarisation, while case B represents the case where an electrode in the window is surrounded by a cathode and an anode.For both cases, the 3d COSICAF simulation is in agreement with the 2d reference potential within a percent.
From these results, we might think that for two configurations out of the four possible (case B and its symmetrical), electrodes located in windows will be almost blind to fission fragment tracks because the electric potential gradient and thus the electric field (Eq.( 28)) is perpendicular to inter-electrode direction.The charge induced by the transport of electrons and ions can be negligible on those electrodes because of the small weighting potential difference between location of creation and annihilation of charges.
Another probable problem with the stacked mpfd configuration is cross talk: some pulses might be detected on more than two electrodes.To probe whether this effect is important, the Shockley-Ramo Weighting potential was computed again with a reference software and with COSICAF.Results are available in Figure 17 for FEMM.Surprisingly, the weighting potential of electrodes enclosed in the dielectric seems more spread than the one located in alumina window.As a result, in the configuration B, the having the weighting potential wa (Fig. 17) will see larger pulses than the electrode having the weighting potential wb.This can be understood qualitatively by considering induced charge on an electrode: a charged particle created in the inter-electrode space generates an induced charge which can be estimated by means of the difference in weighting potential between the starting point and the end point of the charged particle trajectory (Fig. 17).In configuration B the charged particles move laterally because of electric potential gradient (Fig. 16), so the induced charge on the electrode at the weighting potential wa will be greater than the one on the electrode at wb potential.Due to this lateral movement, some pulses might even be positive and negative with an almost null charge on electrode with the wb potential.Further investigations have to be done to show if the origin of the signal is not mistaken during multi-mpdf measurements.

Signal computation
A total of 300 pulses were computed for an mpfd both with the A (Fig. 19) and B (Fig. 18) electrostatic configuration.To better understand signal formation, the two electrodes related to weighting field wa and wb were considered.Due to the complexity of these configurations, pulses have many different shapes.Charges drifting in a direction perpendicular to the weighting field tend to generate pulse with a negative and a positive current.Distance to electrode also plays a key role since fission tracks generated at the opposite of the "signal" electrode almost do not contribute to the output signal.In configuration B, it is clear that electrode with the weighting potential wa sees more signal than the electrode with weighting potential wb.Due to the design, crosstalk between electrodes is ineluctable whatever the fission track considered.Results obtained on this mpfd simulation tend to show that this design is not optimum.The variety of induced current pulses makes signal processing and neutron/gamma discrimination difficult.In addition, signal from more than one mpfd might be recorded on a single electrode, which reduces the interest of different fissile deposits in an mpfd stack.The complex geometry implies huge computation time to estimate electric field and weighting field without noticeable noise: for each configuration, 60 millions of random walks were needed.Because of the reduced pressure and the small size of the ionization chamber, the straight line trajectory hypothesis is appropriate.A more complete study has to be done to confirm preliminary trends observed with the simulation of a few pulses and with electrostatic computations.
Meanwhile, a new mpfd design [48] has been considered.Experimental results are promising and, from back of the envelope calculation, this new design seems to limit drastically the amount of crosstalk and complex electric field by choosing a coaxial configuration of the electrodes.
Despite its simple modelisation of the physics involved in ionisation chambers, COSICAF is successful to simulate complex detectors for research purpose which is for instance valuable to detect possible design flaws.

Conclusion
The new fission chamber simulation code COSICAF, mainly dedicated to education purpose, was presented in this paper.Thanks to simple semi-empirical models, students do not need to spend time on gathering useful data and are thus able to focus on nuclear instrumentation and simulation problems.All the data needed to simulate fission chamber are delivered within the paper, practical work can be derived from it and detector signal could be simulated even with pocket calculator.From a research point of view, COSICAF can also be used to quickly prototype innovative detectors.From a simple geometry definition and a script, it is possible to estimate the shape of the electric field in the detector, as well as the typical fission tracks and the signal they induce on electrodes.
The simplicity of the code comes at the price of a number of simplifications.The simulation is restricted to fission chambers operating in saturation mode.The generation of fission fragments is uncoupled from the fission rate calculation, and only two average fission fragments are considered.In addition, each ionization track is simulated independently of the others.Thus, collective effects such as space charge generation cannot be simulated.Finally, the trajectories of the heavy ions are assumed to be straight.Despite these simplifications, most of the results obtained are similar to those of other existing simulation codes.However, because of the hypothesis for heavy ion transport, spurious current pulses might appear during simulation for particular geometry.Advanced result analysis and expertise is mandatory to use COSICAF for quantitative simulation.The simulation code is provided with the paper or directly on the web [49], it could be extended to simulate other nuclear instruments like boron detectors or proportional counters.

Fig. 3 .
Fig. 3. Two example of flat (left) and cylindrical (right) fission chamber designs built by the CEA fission chamber workshop.

Fig. 5 .
Fig.5.Electronic and nuclear stopping power of a U235 light fission fragment from SRIM.For convenience, the stopping power is in MeV/(mg.cm−2 ).

Fig. 7 .
Fig. 7. Example of geometry description for a detector.The type of volume is depicted as well as the inclusion index.Second and third volumes represent electrodes, density and composition are useless since the stop attribute makes all the particles stops at volume boundaries.Fourth and fifth volumes are heavy ion sources, to avoid undetermined states during heavy ion transport, one volume is made transparent to particles.

Fig. 8 .
Fig.8.Example of a fission fragment propagation across a planar geometry.Green ellipses represent the selected exit point while grey ones are those computed but discarded.The red cylinder correspond to a U3O8 deposit from which fission fragments originated.The light blue cylinder is a volume of argon while the outer dark blue cylinder shell is alumina.

Fig. 9 .
Fig. 9. Example of a potential computation using random walk in 2 dimensions.

Fig. 11 .
Fig. 11.Top: comparison between COSICAF and TRIM energy deposition for various trajectories.The straight trajectory hypothesis fails only for long trajectories.Bottom: charges induced by electrons (black) and ions (yellow) at the electrode compared to the total amount of electron charge generated by the heavy ion in the inter-electrode gap (green).

Fig. 12 .
Fig. 12. Left: electronic pulses computed with the CFUR geometry.Right: ionic pulses computed for few fission tracks generated by heavy and light fission fragments.

5. 1
Geometrical model Micro-pocket fission detectors (mpfd) are a new kind of nuclear instrumentation developed by Kansas State

Fig. 14 .
Fig. 14.Top and lateral cross sections of a micro-pocket fission detector.Only two of the 8 electrodes are exposed to the ionization chamber.

Fig. 15 .
Fig. 15.Top view in transparency of the approximated mpfd geometry drawn by COSICAF.Alumina casing on the side of the ionisation chamber is approximated by cylinder shells and rods (in yellow) while electrodes are depicted in red.

Fig. 16 .
Fig.16.Electric potential computed with CASICAF's MC electrostatic solver on a 3D micro pocket detector geometry (Fig.15).Two cases were studied to consider all the possible configurations.The four anodes at +200 V are delimited by blue squares while cathodes at zero potential are not noticeable.The voxel geometry has a 0.05 mm discretization step.Few million random walks were used to reduce the uncertainty.

Fig. 17 .
Fig. 17.Weighting potential of two electrodes computed with FEMM on a 2D micro pocket detector.

Fig. 18 .
Fig. 18.Pulses simulated for the mpfd with the electrostatic potential B and the weighting fields wa (left) and wb (right).Pulses with the same line style correspond to the same fission track.

Fig. 19 .
Fig. 19.Pulses simulated for the mpfd with the electrostatic potential A and the weighting fields wa (left) and wb (right).

Table 1 .
Characteristics of the two mean fission fragments generated during the fission of uranium 235 and plutonium 239.

Table 2 .
Average energy required to produce an ion pair in different gas for fission fragments.