From 1c4e36e39a10928cc7cd0a80cbac05cb41228a97 Mon Sep 17 00:00:00 2001 From: retkiewi Date: Fri, 14 Jun 2024 15:41:02 +0200 Subject: [PATCH] Squashed commit of the following: commit 377d0f58330b4b2343d5cdf8e4cd977a334095cb Merge: 3448172 e4bc7ec Author: retkiewi <58555777+retkiewi@users.noreply.github.com> Date: Sun Jun 9 15:50:22 2024 +0200 Merge pull request #55 from jbrage/feature/43-reorganize-hadrons Feature/43 reorganize electrons commit e4bc7eca5032f54d246a3bc977bc479ed35dd648 Author: retkiewi Date: Sat Jun 8 19:57:37 2024 +0200 improve runners, beam simulation and use more numpy commit 532ac8e87e22e84e94a0cd53ae094fa129bc4cb5 Author: retkiewi Date: Sat May 25 20:45:01 2024 +0200 fix nit commit 13ed258b9a9c64956eb949cafe6291e1b66e9597 Author: retkiewi Date: Sat May 25 20:39:26 2024 +0200 move calculate method to generic solver commit 08bc5889bc6c1308286b71396a3402ab45c76c5b Author: retkiewi Date: Sat May 25 20:31:17 2024 +0200 only return f in run simulation for continous commit f2157e3a513a1a72c5a1dc1c5d18fcc098bb264a Author: retkiewi Date: Sat May 25 20:22:43 2024 +0200 unify calculate method for pulsed and continous beams commit 9a1247c26410f2cc7ba05641754ef1d389cf6f70 Author: retkiewi Date: Sat May 25 20:10:25 2024 +0200 return f as an array of results after each time step commit 2323a81b2080c7dba5eb10c83bff83d50c7f0f8e Author: retkiewi Date: Sat May 25 19:27:42 2024 +0200 fix getting solver name in runner commit 5ed86193e8339fab7d5b22bce2a575f55df52e33 Author: retkiewi Date: Sat May 25 19:14:38 2024 +0200 update simulation runners commit a2c3e3beaf80e2747194fcb1eacbe46f4f417550 Author: retkiewi Date: Sat May 25 18:32:33 2024 +0200 use dataclasses for solvers commit 55aa3aad9ebf1c5954366c32bc382c4ffb3249b0 Author: retkiewi Date: Sun May 19 16:32:29 2024 +0200 add simulation runners, fix bugs commit 3240b58f98e097112a28722ddfcc6897a87002da Author: retkiewi Date: Sun May 19 15:29:12 2024 +0200 optimize constant calculations commit 6c77bb74cf7e157fa76c1f25da617d80d702f7a1 Author: retkiewi Date: Sun May 19 15:27:00 2024 +0200 remove parameter dict in generic solver commit 70a3ef5a33d96b6bd759da0008ac3a42d0752071 Author: retkiewi Date: Sun May 19 15:09:51 2024 +0200 fix generic_electron_solver formatting commit d29a577872f8de3df5743d3fe1a0810c0d5f72cd Author: retkiewi Date: Sat May 18 19:28:05 2024 +0200 extract post beam electron density calculation to a separate method commit 13053506588754a1504b87a7f2aa110a7ad1bef2 Author: retkiewi Date: Sat May 18 19:18:20 2024 +0200 extract common init commit 841d56abd079f2a45d361b6f906cc6762b8f52ba Author: retkiewi Date: Sat May 18 18:53:48 2024 +0200 small cleanup of pulsed_e_beam.py commit 68e9ff20520c39591b1cb89cd4d6461e9670b28a Author: retkiewi Date: Sat May 18 18:50:46 2024 +0200 convert electrons/continous_e_beam.py to python commit 6b4df88585cd70994bb4f313bdbdd28b38ea22db Author: retkiewi Date: Sat May 18 18:42:14 2024 +0200 large cleanup in electrons/continous_e_beam.py commit 541d5cc6de17fb8e216e977bb3a3ff4f3f8009d7 Author: retkiewi Date: Sat May 18 18:41:46 2024 +0200 small cleanup in initial_recombination.py commit a404c0462347126cb3eca335ff73d8e12e468709 Author: retkiewi Date: Sat May 18 18:03:00 2024 +0200 convert electrons/continous_e_beam.py to python commit 3448172f91d402388636fa3f61322ec36d01b575 Author: dependabot[bot] <49699333+dependabot[bot]@users.noreply.github.com> Date: Mon Apr 1 15:23:06 2024 +0000 Bump actions/setup-python from 5.0.0 to 5.1.0 Bumps [actions/setup-python](https://github.com/actions/setup-python) from 5.0.0 to 5.1.0. - [Release notes](https://github.com/actions/setup-python/releases) - [Commits](https://github.com/actions/setup-python/compare/v5.0.0...v5.1.0) --- updated-dependencies: - dependency-name: actions/setup-python dependency-type: direct:production update-type: version-update:semver-minor ... Signed-off-by: dependabot[bot] commit d6c0035c017b4462f0ef00cf4881d84b36116e00 Author: dependabot[bot] <49699333+dependabot[bot]@users.noreply.github.com> Date: Mon Dec 11 15:54:38 2023 +0000 Bump actions/setup-python from 4.7.1 to 5.0.0 Bumps [actions/setup-python](https://github.com/actions/setup-python) from 4.7.1 to 5.0.0. - [Release notes](https://github.com/actions/setup-python/releases) - [Commits](https://github.com/actions/setup-python/compare/v4.7.1...v5.0.0) --- updated-dependencies: - dependency-name: actions/setup-python dependency-type: direct:production update-type: version-update:semver-major ... Signed-off-by: dependabot[bot] commit a246282ccab1726664296842d465ca0174ebb8ee Author: dependabot[bot] <49699333+dependabot[bot]@users.noreply.github.com> Date: Mon Oct 2 20:37:55 2023 +0000 Bump actions/checkout from 3 to 4 Bumps [actions/checkout](https://github.com/actions/checkout) from 3 to 4. - [Release notes](https://github.com/actions/checkout/releases) - [Changelog](https://github.com/actions/checkout/blob/main/CHANGELOG.md) - [Commits](https://github.com/actions/checkout/compare/v3...v4) --- updated-dependencies: - dependency-name: actions/checkout dependency-type: direct:production update-type: version-update:semver-major ... Signed-off-by: dependabot[bot] commit 8ef4655e0c0ec30b3c8a289937ba0deaa18027ac Author: dependabot[bot] <49699333+dependabot[bot]@users.noreply.github.com> Date: Mon Oct 2 15:30:14 2023 +0000 Bump actions/setup-python from 4.5.0 to 4.7.1 Bumps [actions/setup-python](https://github.com/actions/setup-python) from 4.5.0 to 4.7.1. - [Release notes](https://github.com/actions/setup-python/releases) - [Commits](https://github.com/actions/setup-python/compare/v4.5.0...v4.7.1) --- updated-dependencies: - dependency-name: actions/setup-python dependency-type: direct:production update-type: version-update:semver-minor ... Signed-off-by: dependabot[bot] commit 0dce892be5e450bd0a62dddc117f0f535489ce29 Author: dependabot[bot] <49699333+dependabot[bot]@users.noreply.github.com> Date: Mon Jan 16 15:46:26 2023 +0000 Bump actions/setup-python from 4.4.0 to 4.5.0 Bumps [actions/setup-python](https://github.com/actions/setup-python) from 4.4.0 to 4.5.0. - [Release notes](https://github.com/actions/setup-python/releases) - [Commits](https://github.com/actions/setup-python/compare/v4.4.0...v4.5.0) --- updated-dependencies: - dependency-name: actions/setup-python dependency-type: direct:production update-type: version-update:semver-minor ... Signed-off-by: dependabot[bot] commit b5ca6d82721db1f012b49f65786fe566df3ace29 Author: dependabot[bot] <49699333+dependabot[bot]@users.noreply.github.com> Date: Mon Dec 26 15:06:08 2022 +0000 Bump actions/setup-python from 4.3.1 to 4.4.0 Bumps [actions/setup-python](https://github.com/actions/setup-python) from 4.3.1 to 4.4.0. - [Release notes](https://github.com/actions/setup-python/releases) - [Commits](https://github.com/actions/setup-python/compare/v4.3.1...v4.4.0) --- updated-dependencies: - dependency-name: actions/setup-python dependency-type: direct:production update-type: version-update:semver-minor ... Signed-off-by: dependabot[bot] commit 24ecd26c13127f9bfa58b45fe6c8b4bf9547b7e4 Author: dependabot[bot] <49699333+dependabot[bot]@users.noreply.github.com> Date: Mon Dec 12 15:08:26 2022 +0000 Bump actions/setup-python from 4.3.0 to 4.3.1 Bumps [actions/setup-python](https://github.com/actions/setup-python) from 4.3.0 to 4.3.1. - [Release notes](https://github.com/actions/setup-python/releases) - [Commits](https://github.com/actions/setup-python/compare/v4.3.0...v4.3.1) --- updated-dependencies: - dependency-name: actions/setup-python dependency-type: direct:production update-type: version-update:semver-patch ... Signed-off-by: dependabot[bot] commit d2f22bafe055106c855c9ae35e42d179537afa1e Merge: ee90579 70d9bce Author: retkiewi <58555777+retkiewi@users.noreply.github.com> Date: Thu Dec 8 20:45:07 2022 +0100 Merge pull request #45 from jbrage/feature/43-reorganize-hadrons Reorganize hadrons package --- .github/workflows/main.yml | 4 +- electrons/cython/run_simulation.py | 69 ++++++ electrons/python/continuous_e_beam.py | 38 +++ electrons/python/generic_electron_solver.py | 238 +++++++++++++++++++ electrons/python/pulsed_e_beam.py | 34 +++ electrons/python/run_simulation.py | 63 +++++ hadrons/numba/initial_recombination_numba.py | 19 +- hadrons/python/initial_recombination.py | 5 +- 8 files changed, 455 insertions(+), 15 deletions(-) create mode 100644 electrons/cython/run_simulation.py create mode 100644 electrons/python/continuous_e_beam.py create mode 100644 electrons/python/generic_electron_solver.py create mode 100644 electrons/python/pulsed_e_beam.py create mode 100644 electrons/python/run_simulation.py diff --git a/.github/workflows/main.yml b/.github/workflows/main.yml index 9a9f2ed..66c028b 100644 --- a/.github/workflows/main.yml +++ b/.github/workflows/main.yml @@ -21,10 +21,10 @@ jobs: python-version: ['3.8', '3.9', '3.10'] platform: [ubuntu-latest] # TODO add Windows as well steps: - - uses: actions/checkout@v3 + - uses: actions/checkout@v4 - name: Set up Python - uses: actions/setup-python@v4.3.0 + uses: actions/setup-python@v5.1.0 with: python-version: ${{ matrix.python-version }} diff --git a/electrons/cython/run_simulation.py b/electrons/cython/run_simulation.py new file mode 100644 index 0000000..c94e313 --- /dev/null +++ b/electrons/cython/run_simulation.py @@ -0,0 +1,69 @@ +from electrons.cython.pulsed_e_beam import pulsed_beam_PDEsolver +from electrons.cython.continuous_e_beam import continuous_beam_PDEsolver + +import argparse + +SOLVER_MAP = {"continous": continuous_beam_PDEsolver, "pulsed": pulsed_beam_PDEsolver} + + +def run_simulation( + solver_name="continous", + voltage_V=300, + electrode_gap=0.1, + electron_density_per_cm3=1e9, + verbose=True, +): + parameters = { + "voltage_V": voltage_V, + "d_cm": electrode_gap, + "elec_per_cm3": electron_density_per_cm3, + "show_plot": False, + "print_parameters": False, + } + + if solver_name not in SOLVER_MAP.keys(): + print(f'Invalid solver type "{solver_name}", defaulting to Continous solver.') + solver_name = "continous" + + if verbose: + print(f"Running the simulation using the {solver_name} solver.") + print(f"Voltage: {voltage_V} [V]") + print(f"Electrode gap: {electrode_gap} [cm]") + print(f"Electron density per cm3: {electron_density_per_cm3}") + + solver = SOLVER_MAP[solver_name] + + # return the collection efficiency + result = solver(parameters) + + if solver_name == "continous": + f = result[1]["f"] + else: + f = result + + return f + + +if __name__ == "__main__": + + parser = argparse.ArgumentParser() + + parser.add_argument( + "solver_name", + type=str, + default="continous", + help="The type of the solver to use", + ) + + parser.add_argument( + "--verbose", + "-v", + type=bool, + default=True, + ) + + args = parser.parse_args() + + result = run_simulation(**vars(args)) + + print("calculated f is", result) diff --git a/electrons/python/continuous_e_beam.py b/electrons/python/continuous_e_beam.py new file mode 100644 index 0000000..d2331d5 --- /dev/null +++ b/electrons/python/continuous_e_beam.py @@ -0,0 +1,38 @@ +from __future__ import division + +from generic_electron_solver import GenericElectronSolver + + +class ContinousBeamPDEsolver(GenericElectronSolver): + def get_initialised_charge_carriers_after_beam(self): + electron_density_per_cm3_s = self.electron_density_per_cm3 * self.dt + + initialised_charge_carriers = ( + electron_density_per_cm3_s + * (self.no_xy - 2 * self.delta_border) ** 2 + * self.no_z + ) + + return initialised_charge_carriers + + def update_electron_density_arrays_after_beam( + self, + positive_array, + negative_array, + ): + electron_density_per_cm3_s = self.electron_density_per_cm3 * self.dt + + positive_array[ + self.delta_border : (self.no_xy - self.delta_border), + self.delta_border : (self.no_xy - self.delta_border), + self.no_z_electrode : (self.no_z + self.no_z_electrode), + ] += electron_density_per_cm3_s + + negative_array[ + self.delta_border : (self.no_xy - self.delta_border), + self.delta_border : (self.no_xy - self.delta_border), + self.no_z_electrode : (self.no_z + self.no_z_electrode), + ] += electron_density_per_cm3_s + + def should_simulate_beam_for_time_step(self, time_step): + return True diff --git a/electrons/python/generic_electron_solver.py b/electrons/python/generic_electron_solver.py new file mode 100644 index 0000000..63f37c2 --- /dev/null +++ b/electrons/python/generic_electron_solver.py @@ -0,0 +1,238 @@ +from __future__ import division +from abc import ABC, abstractmethod +from dataclasses import dataclass +from numpy.typing import NDArray +from typing import Tuple +import numpy as np + + +def von_neumann_expression( + dt, ion_diff, grid_spacing_cm, ion_mobility, Efield_V_cm +) -> None: + """ + Finds a time step dt which fulfils the von Neumann criterion, i.e. ensures the numericl error does not increase but + decreases and eventually damps out + """ + von_neumann_expression = False + + while not von_neumann_expression: + dt /= 1.01 + # as defined in the Deghan (2004) paper + # sx = ion_diff*dt/(self.grid_spacing_cm**2) + # sy = ion_diff*dt/(self.grid_spacing_cm**2) + sx = sy = 0.0 # uniform charge => no gradient driven diffusion in the xy plane + + sz = ion_diff * dt / (grid_spacing_cm**2) + cx = cy = 0.0 + cz = ion_mobility * Efield_V_cm * dt / grid_spacing_cm + # check von Neumann's criterion + criterion_1 = 2 * (sx + sy + sz) + cx**2 + cy**2 + cz**2 <= 1 + + criterion_2 = cx**2 * cy**2 * cz**2 <= 8 * sx * sy * sz + + von_neumann_expression = criterion_1 and criterion_2 + + return dt, sx, sy, sz, cx, cy, cz + + +@dataclass +class GenericElectronSolver(ABC): + # Simulation parameters + electron_density_per_cm3: float # fluence-rate [/cm^2/s] + voltage_V: float # [V/cm] magnitude of the electric field + electrode_gap: float # [cm] # electrode gap + grid_spacing_cm: float + ion_mobility: float = 1.73 # cm s^-1 V^-1, averaged for positive and negative ions + ion_diff: float = 3.7e-2 # cm^2/s, averaged for positive and negative ions + alpha: float = 1.60e-6 # cm^3/s, recombination constant + r_cm: float = 0.002 # radius of the sampled cylinder + buffer_radius: int = 5 + no_z_electrode: int = ( + 5 # length of the electrode-buffer to ensure no ions drift through the array in one time step + ) + dt: float = 1.0 + + # Derived properties + @property + def no_xy(self) -> int: + """Number of voxels om the xy-directions""" + + no_xy = int(2 * self.r_cm / self.grid_spacing_cm) + no_xy += 2 * self.buffer_radius + return no_xy + + @property + def no_z(self) -> int: + """Number of voxels om the z-direction""" + + return int(self.electrode_gap / self.grid_spacing_cm) + + @property + def no_z_with_buffer(self) -> int: + return 2 * self.no_z_electrode + self.no_z + + @property + def Efield_V_cm(self) -> float: + return self.voltage_V / self.electrode_gap + + @property + def separation_time_steps(self) -> int: + """ + Number of step required to drag the two charge carrier distributions apart + along with the number of tracks to be uniformly distributed over the domain + """ + return int( + self.electrode_gap / (2.0 * self.ion_mobility * self.Efield_V_cm * self.dt) + ) + + @property + def computation_time_steps(self) -> int: + return self.separation_time_steps * 3 + + @property + def delta_border(self) -> int: + return 2 + + def __post_init__( + self, + ): + # depending on the cluster/computer, the upper limit may be changed + if (self.no_xy * self.no_xy * self.no_z) > 1e8: + raise ValueError( + "Too many elements in the array: %i" + % (self.no_xy * self.no_xy * self.no_z_with_buffer) + ) + + self.dt, self.sx, self.sy, self.sz, self.cx, self.cy, self.cz = ( + von_neumann_expression( + self.dt, + self.ion_diff, + self.grid_spacing_cm, + self.ion_mobility, + self.Efield_V_cm, + ) + ) + + @abstractmethod + def should_simulate_beam_for_time_step(self, time_step: int) -> bool: + pass + + @abstractmethod + def update_electron_density_arrays_after_beam( + self, positive_array: NDArray, negative_array: NDArray + ) -> None: + pass + + @abstractmethod + def get_initialised_charge_carriers_after_beam() -> float: + pass + + def simulate_beam(self, positive_array: NDArray, negative_array: NDArray): + self.update_electron_density_arrays_after_beam(positive_array, negative_array) + + return self.get_initialised_charge_carriers_after_beam() + + def calculate(self): + positive_array = np.zeros((self.no_xy, self.no_xy, self.no_z_with_buffer)) + negative_array = np.zeros((self.no_xy, self.no_xy, self.no_z_with_buffer)) + positive_array_temp = np.zeros((self.no_xy, self.no_xy, self.no_z_with_buffer)) + negative_array_temp = np.zeros((self.no_xy, self.no_xy, self.no_z_with_buffer)) + no_recombined_charge_carriers = 0.0 + no_initialised_charge_carriers = 0.0 + + """ + Create an array which includes the number of track to be inserted for each time step + The tracks are distributed uniformly in time + """ + + positive_temp_entry = negative_temp_entry = recomb_temp = 0.0 + + f_steps_list = np.zeros(self.computation_time_steps) + + szcz_pos = self.sz + self.cz * (self.cz + 1.0) / 2.0 + szcz_neg = self.sz + self.cz * (self.cz - 1.0) / 2.0 + + sycy_pos = self.sy + self.cy * (self.cy + 1.0) / 2.0 + sycy_neg = self.sy + self.cy * (self.cy - 1.0) / 2.0 + + sxcx_pos = self.sx + self.cx * (self.cx + 1.0) / 2.0 + sxcx_neg = self.sx + self.cx * (self.cx - 1.0) / 2.0 + + cxyzsyz = ( + 1.0 + - self.cx * self.cx + - self.cy * self.cy + - self.cz * self.cz + - 2.0 * (self.sx + self.sy + self.sz) + ) + + """ + Start the simulation by evolving the distribution one step at a time + """ + + for time_step in range(self.computation_time_steps): + + """ + Refill the array with the electron density each time step + """ + if self.should_simulate_beam_for_time_step(time_step): + no_initialised_step = self.simulate_beam( + positive_array, + negative_array, + ) + else: + no_initialised_step = 0.0 + + no_initialised_charge_carriers += no_initialised_step + + # calculate the new densities and store them in temporary arrays + for i in range(1, self.no_xy - 1): + for j in range(1, self.no_xy - 1): + for k in range(1, self.no_z_with_buffer - 1): + # using the Lax-Wendroff scheme + positive_temp_entry = szcz_pos * positive_array[i, j, k - 1] + positive_temp_entry += szcz_neg * positive_array[i, j, k + 1] + + positive_temp_entry += sycy_pos * positive_array[i, j - 1, k] + positive_temp_entry += sycy_neg * positive_array[i, j + 1, k] + + positive_temp_entry += sxcx_pos * positive_array[i - 1, j, k] + positive_temp_entry += sxcx_neg * positive_array[i + 1, j, k] + + positive_temp_entry += cxyzsyz * positive_array[i, j, k] + + # same for the negative charge carriers + negative_temp_entry = szcz_pos * negative_array[i, j, k + 1] + negative_temp_entry += szcz_neg * negative_array[i, j, k - 1] + + negative_temp_entry += sycy_pos * negative_array[i, j + 1, k] + negative_temp_entry += sycy_neg * negative_array[i, j - 1, k] + + negative_temp_entry += sxcx_pos * negative_array[i + 1, j, k] + negative_temp_entry += sxcx_neg * negative_array[i - 1, j, k] + + negative_temp_entry += cxyzsyz * negative_array[i, j, k] + + # the recombination part + recomb_temp = ( + self.alpha + * positive_array[i, j, k] + * negative_array[i, j, k] + * self.dt + ) + positive_array_temp[i, j, k] = positive_temp_entry - recomb_temp + negative_array_temp[i, j, k] = negative_temp_entry - recomb_temp + + if k > self.no_z_electrode and k < ( + self.no_z + self.no_z_electrode + ): + # sum over the recombination between the virtual electrodes + no_recombined_charge_carriers += recomb_temp + + f_steps_list[time_step] = ( + no_initialised_charge_carriers - no_recombined_charge_carriers + ) / no_initialised_charge_carriers + np.copyto(positive_array, positive_array_temp) + np.copyto(negative_array, negative_array_temp) + + return f_steps_list diff --git a/electrons/python/pulsed_e_beam.py b/electrons/python/pulsed_e_beam.py new file mode 100644 index 0000000..8a295f5 --- /dev/null +++ b/electrons/python/pulsed_e_beam.py @@ -0,0 +1,34 @@ +from __future__ import division +from generic_electron_solver import GenericElectronSolver + + +class PulsedBeamPDEsolver(GenericElectronSolver): + + def get_initialised_charge_carriers_after_beam(self): + initialised_charge_carriers = ( + self.electron_density_per_cm3 + * (self.no_xy - 2 * self.delta_border) ** 2 + * self.no_z + ) + + return initialised_charge_carriers + + def update_electron_density_arrays_after_beam( + self, + positive_array, + negative_array, + ): + positive_array[ + self.delta_border : (self.no_xy - self.delta_border), + self.delta_border : (self.no_xy - self.delta_border), + self.no_z_electrode : (self.no_z + self.no_z_electrode), + ] += self.electron_density_per_cm3 + + negative_array[ + self.delta_border : (self.no_xy - self.delta_border), + self.delta_border : (self.no_xy - self.delta_border), + self.no_z_electrode : (self.no_z + self.no_z_electrode), + ] += self.electron_density_per_cm3 + + def should_simulate_beam_for_time_step(self, time_step): + return time_step == 0 diff --git a/electrons/python/run_simulation.py b/electrons/python/run_simulation.py new file mode 100644 index 0000000..8bc6e26 --- /dev/null +++ b/electrons/python/run_simulation.py @@ -0,0 +1,63 @@ +from continuous_e_beam import ContinousBeamPDEsolver +from pulsed_e_beam import PulsedBeamPDEsolver +import sys +import argparse + +SOLVER_MAP = {"continous": ContinousBeamPDEsolver, "pulsed": PulsedBeamPDEsolver} + + +def run_simulation( + solver_name="continous", + voltage_V=300, + electrode_gap=0.1, + electron_density_per_cm3=1e9, + verbose=True, +): + Solver = ( + SOLVER_MAP[solver_name] + if solver_name in SOLVER_MAP.keys() + else ContinousBeamPDEsolver + ) + + if solver_name not in SOLVER_MAP.keys(): + print(f'Invalid solver type "{solver_name}", defaulting to Continous solver.') + solver_name = "continous" + + if verbose: + print(f"Running the simulation using the {solver_name} solver.") + print(f"Voltage: {voltage_V} [V]") + print(f"Electrode gap: {electrode_gap} [cm]") + print(f"Electron density per cm3: {electron_density_per_cm3}") + + solver = Solver( + electron_density_per_cm3=electron_density_per_cm3, + voltage_V=voltage_V, + electrode_gap=electrode_gap, + grid_spacing_cm=5e-4, + ) + + return solver.calculate() + + +if __name__ == "__main__": + parser = argparse.ArgumentParser() + + parser.add_argument( + "solver_name", + type=str, + default="continous", + help="The type of the solver to use", + ) + + parser.add_argument( + "--verbose", + "-v", + type=bool, + default=True, + ) + + args = parser.parse_args() + + result = run_simulation(**vars(args)) + + print(f"calculated f is {result}") diff --git a/hadrons/numba/initial_recombination_numba.py b/hadrons/numba/initial_recombination_numba.py index ed56c47..47f73e7 100644 --- a/hadrons/numba/initial_recombination_numba.py +++ b/hadrons/numba/initial_recombination_numba.py @@ -89,18 +89,18 @@ def main_loop( negative_array_temp = np.zeros((no_x, no_x, no_z_with_buffer)) no_recombined_charge_carriers = 0.0 - for _ in range(computation_time_steps): - # calculate the new densities and store them in temporary arrays - # dunno what would be a good name for those coefficients - szcz_pos = sz + cz * (cz + 1.0) / 2.0 - szcz_neg = sz + cz * (cz - 1.0) / 2.0 + # calculate the new densities and store them in temporary arrays + # dunno what would be a good name for those coefficients + szcz_pos = sz + cz * (cz + 1.0) / 2.0 + szcz_neg = sz + cz * (cz - 1.0) / 2.0 - sycy_pos = sy + cy * (cy + 1.0) / 2.0 - sycy_neg = sy + cy * (cy - 1.0) / 2.0 + sycy_pos = sy + cy * (cy + 1.0) / 2.0 + sycy_neg = sy + cy * (cy - 1.0) / 2.0 - sxcx_pos = sx + cx * (cx + 1.0) / 2.0 - sxcx_neg = sx + cx * (cx - 1.0) / 2.0 + sxcx_pos = sx + cx * (cx + 1.0) / 2.0 + sxcx_neg = sx + cx * (cx - 1.0) / 2.0 + for _ in range(computation_time_steps): for i in range(1, no_x - 1): for j in range(1, no_x - 1): for k in range(1, no_z_with_buffer - 1): @@ -139,6 +139,7 @@ def main_loop( positive_array_temp[i, j, k] = positive_temp_entry - recomb_temp negative_array_temp[i, j, k] = negative_temp_entry - recomb_temp + no_recombined_charge_carriers += recomb_temp for i in range(1, no_x - 1): diff --git a/hadrons/python/initial_recombination.py b/hadrons/python/initial_recombination.py index 59874ff..f132b3c 100644 --- a/hadrons/python/initial_recombination.py +++ b/hadrons/python/initial_recombination.py @@ -1,5 +1,4 @@ import numpy as np -import numpy.random as rnd import time from math import exp, sqrt, pi, log, cos, sin from ..geiss_utils import Geiss_r_max, Geiss_RRD_cm @@ -127,9 +126,7 @@ def solve(self): if self.RDD_model == "Gauss": def RDD_function(r_cm): - return self.Gaussian_factor * exp( - -(r_cm**2) / self.track_radius_cm**2 - ) + return self.Gaussian_factor * exp(-(r_cm**2) / self.track_radius_cm**2) elif self.RDD_model == "Geiss":