Commit 4d0355d9 authored by Paul Schütze's avatar Paul Schütze
Browse files

Merge branch 'rkf' into 'master'

Propagation Code: Refactoring, Stage 2

See merge request allpix-squared/allpix-squared!1043
parents cad937f5 0f17c80b
Loading
Loading
Loading
Loading
+24 −0
Original line number Diff line number Diff line
@@ -424,6 +424,30 @@ month={8},
  Type={NASA Technical Report},
  note={\url{http://hdl.handle.net/2060/19690021375}}
}
@article{fehlberg2,
    author={Fehlberg, Erwin},
    title={Klassische Runge-Kutta-Formeln f\"unfter und siebenter Ordnung mit Schrittweiten-Kontrolle},
    journal={Computing},
    volume={4},
    pages={93–106},
    year={1969},
    doi={10.1007/BF02234758},
}
@article{cashkarp,
  author = {Cash, J. R. and Karp, Alan H.},
  title = {A Variable Order Runge-Kutta Method for Initial Value Problems with Rapidly Varying Right-Hand Sides},
  year = {1990},
  publisher = {Association for Computing Machinery},
  address = {New York, NY, USA},
  volume = {16},
  number = {3},
  issn = {0098-3500},
  doi = {10.1145/79505.79507},
  journal = {ACM Trans. Math. Softw.},
  month = {sep},
  pages = {201–222},
  numpages = {22}
}

@article{pai,
      author         = "Apostolakis, J. and Giani, S. and Urban, L. and Maire, M.
+81 −0
Original line number Diff line number Diff line
@@ -77,6 +77,85 @@ auto step = runge_kutta.step();
The `getValue()` and `setValue()` methods allow to retrieve, alter and update the position, e.g. to include additional
displacements from diffusion processes.

Furthermore, the `advanceTime()` method can be used to advance the time of the Runge-Kutta solver by the given amount. This
is used for example in situations where the motion is temporarily halted by trapping, and continued when released from the
trap.

### Tableaus

The following Butcher tableaus are implemented and available.
In order to make use of adaptive step size changes, a tableau with error estimation should be chosen.


#### Third-Order Kutta Method (RK3)

This tableau implements a simple and fast third-order Kutta integration which only requires the calculation of three terms.
```math
\begin{array}
{c|ccc}
0                   \\
1/2 & 1/2             \\
1 &  -1 &   2       \\
\hline
  & 1/6 & 2/3 & 1/6 \\
\end{array}
```

#### Fourth-Order Runge-Kutta Method (RK4)

This tableau implements the classical fourth-order Runge-Kutta integration.
```math
\begin{array}
{c|cccc}
  0                 \\
1/2 & 1/2           \\
1/2 &   0 & 1/2     \\
  1 &   0 &   0 & 1 \\
\hline
    & 1/6 & 1/3 & 1/3 & 1/6
\end{array}
```

#### Fourth-Order Runge-Kutta-Fehlberg Method with Error Estimation (RK5)

This tableau implements a fourth-order RKF method with fifth-order error estimation \[[@fehlberg], [@fehlberg2]\].

```math
\begin{array}
{c|cccccc}
    0                                                                     \\
  1/4 &      1/4                                                          \\
  3/8 &      3/32 &       9/32                                            \\
12/13 & 1932/2197 & -7200/2197 &  7296/2197                               \\
    1 &   439/216 &         -8 &   3680/513 &   -845/4104                 \\
  1/2 &     -8/27 &          2 & -3544/2565 &   1859/4104 & -11/40        \\
\hline
      &    16/135 &          0 & 6656/12825 & 28561/56430 &  -9/50 & 2/55 \\
      &    25/216 &          0 &  1408/2565 &   2197/4104 &   -1/5 &    0 \\
\end{array}
```

#### Fourth-Order Runge-Kutta-Cash-Karp Method with Error Estimation (RKCK)

This tableau implements a fourth-order Cash-Karp method with fifth-order error estimation \[[@cashkarp]\].

```math
\begin{array}
{c|cccccc}
   0                                                                            \\
 1/5 &        1/5                                                               \\
3/10 &       3/40 &    9/40                                                     \\
 3/5 &       3/10 &   -9/10 &         6/5                                       \\
   1 &     -11/54 &     5/2 &      -70/27 &        35/27                        \\
 7/8 & 1631/55296 & 175/512 &   575/13824 & 44275/110592 &  253/4096            \\
\hline
     &     37/378 &       0 &     250/621 &      125/594 &         0 & 512/1771 \\
     & 2825/27648 &       0 & 18575/48384 &  13525/55296 & 277/14336 &      1/4 \\

\end{array}
```


## Field Data Parser

A field parser tool is provided, which parses files stored in the INIT or APF file formats and returns field data on a
@@ -116,3 +195,5 @@ file to be binary and parses the field as APF data.

[@eigen3]: http://eigen.tuxfamily.org
[@fehlberg]: https://ntrs.nasa.gov/search.jsp?R=19690021375
[@fehlberg2]: https://doi.org/10.1007%2FBF02234758
[@cashkarp]: https://doi.org/10.1145/79505.79507
+3 −1
Original line number Diff line number Diff line
@@ -20,7 +20,7 @@ To simulate impact ionization, the number of newly generated electron-hole pairs
The two parameters `propagate_electrons` and `propagate_holes` allow to control which type of charge carrier is propagated to their respective electrodes. Either one of the carrier types can be selected, or both can be propagated. It should be noted that this will slow down the simulation considerably since twice as many carriers have to be handled and it should only be used where sensible.
The direction of the propagation depends on the electric and magnetic fields field configured, and it should be ensured that the carrier types selected are actually transported to the implant side. For linear electric fields, a warning is issued if a possible misconfiguration is detected.

A fourth-order Runge-Kutta-Fehlberg method \[[@fehlberg]\] with fifth-order error estimation is used to integrate the particle propagation in the electric and magnetic fields. After every Runge-Kutta step, the diffusion is accounted for by applying an offset drawn from a Gaussian distribution calculated from the Einstein relation
A fourth-order Runge-Kutta-Fehlberg method \[[@fehlberg], [@fehlberg2]\] with fifth-order error estimation, RKF4(5), is used to integrate the particle propagation in the electric and magnetic fields. After every Runge-Kutta step, the diffusion is accounted for by applying an offset drawn from a Gaussian distribution calculated from the Einstein relation

$`\sigma = \sqrt{\frac{2k_b T}{e}\mu t}`$

@@ -99,3 +99,5 @@ charge_per_step = 25
```

[@fehlberg]: https://ntrs.nasa.gov/search.jsp?R=19690021375
[@fehlberg2]: https://doi.org/10.1007%2FBF02234758
+1 −2
Original line number Diff line number Diff line
@@ -16,7 +16,7 @@ The propagation consists of a combination of drift and diffusion simulation. The
This module implements charge multiplication by impact ionization. The multiplication model can be chosen using the `multiplication_model` parameter, the list of available models can be found in the user manual. By default, the model defaults to `none` and impact ionization is switched off, generating unity gain.
To simulate impact ionization, the number of newly generated electron-hole pairs is calculated for every propagation step and every charge carrier in the group, based on drawing a random number from a geometric distribution. This represents a stepwise approach to the avalanche generation process. The charge of a charge group is increased by the number of impact ionization processes per step and opposite-type charge carriers are generated at the end of the step.

A fourth-order Runge-Kutta-Fehlberg method \[[@fehlberg]\] is used to integrate the particle motion through the electric and magnetic fields. After every Runge-Kutta step, the diffusion is accounted for by applying an offset drawn from a Gaussian distribution calculated from the Einstein relation
A classic fourth-order Runge-Kutta method is used to integrate the particle motion through the electric and magnetic fields. After every Runge-Kutta step, the diffusion is accounted for by applying an offset drawn from a Gaussian distribution calculated from the Einstein relation

$`\sigma = \sqrt{\frac{2k_b T}{e}\mu t}`$

@@ -93,6 +93,5 @@ output_plots = true
timestep = 0.02ns
```

[@fehlberg]: https://ntrs.nasa.gov/search.jsp?R=19690021375
[@shockley]: https://doi.org/10.1063/1.1710367
[@ramo]: https://doi.org/10.1109/JRPROC.1939.228757
+2 −2
Original line number Diff line number Diff line
@@ -581,9 +581,9 @@ TransientPropagationModule::propagate(Event* event,
        return static_cast<int>(type) * mob * (efield + term1 + term2) / rnorm;
    };

    // Create the runge kutta solver with an RKF5 tableau
    // Create the runge kutta solver with an RK4 tableau, no error estimation required since we're not adapting step size
    auto runge_kutta = make_runge_kutta(
        tableau::RK5, (has_magnetic_field_ ? carrier_velocity_withB : carrier_velocity_noB), timestep_, position);
        tableau::RK4, (has_magnetic_field_ ? carrier_velocity_withB : carrier_velocity_noB), timestep_, position);

    // Continue propagation until the deposit is outside the sensor
    Eigen::Vector3d last_position = position;
Loading