From dba0ce83b70c8c9a960b96295eb2ba3cc0a5ad79 Mon Sep 17 00:00:00 2001 From: zpdg Date: Mon, 21 Apr 2025 18:46:40 +0200 Subject: [PATCH 1/7] Add support for simulation for dead layer-Step3: BulkSurfaceConstantLifetimeChargeTrappingModel --- docs/src/man/charge_drift.md | 24 +++++++- src/ChargeTrapping/ChargeTrapping.jl | 75 +++++++++++++++++++++++++ src/SolidStateDetector/Semiconductor.jl | 5 ++ src/SolidStateDetectors.jl | 1 + test/test_charge_drift_models.jl | 17 +++++- 5 files changed, 120 insertions(+), 2 deletions(-) diff --git a/docs/src/man/charge_drift.md b/docs/src/man/charge_drift.md index 7f6ca2c40..46bd90b19 100644 --- a/docs/src/man/charge_drift.md +++ b/docs/src/man/charge_drift.md @@ -219,7 +219,7 @@ In this case, the signal is calculated using the Schockley-Ramo theorem, i.e. by ### `BoggsChargeTrappingModel` -In SolidStateDetectors.jl, the only charge trapping model implemented so far is the one presented in [Boggs _et al._ (2023)](https://doi.org/10.1016/j.nima.2023.168756). +In SolidStateDetectors.jl, the Boggs charge trapping model implemented is the one presented in [Boggs _et al._ (2023)](https://doi.org/10.1016/j.nima.2023.168756). In this model, the charge cloud loses part of its charge at every point `path[i]` of the charge drift, depending on its drift and thermal velocity, as well as the trapping product $[n\sigma_{e/h}]^{-1}$ for electrons and holes. The charge signal is then given by the charge-decreased charge cloud reaching the contacts and the charges trapped on the way. @@ -261,6 +261,28 @@ parameters = Dict("parameters" => sim.detector = SolidStateDetector(sim.detector, BoggsChargeTrappingModel{T}(parameters)) ``` +### `BulkSurfaceConstantLifetimeTrappingModel` +This constant-lifetime-based charge trapping model is similar to the Boggs model, which is constant-mean-free-path based. +This idea is presented in [Dai _et al._ (2023)](10.1016/j.apradiso.2022.110638). +Besides, the model implemented has two sets of parameters for the bulk (sensitive region) and surface (dead layer/inactive layer) volume respectively. +It is implemented mainly for the pulse shape simulation of the surface layer. + +In this model, the charge cloud loses part of its charge at every point `path[i]` of the charge drift, depending on the lifetime. + +The charge signal is then given by the charge-decreased charge cloud reaching the contacts and the charges trapped on the way. + +It can be applied to an already read-in `SolidStateDetector` using for example: +```julia +using SolidStateDetectors, Unitful +T = Float32 +sim = Simulation{T}(SSD_examples[:InvertedCoax]) + +if_in_bulk = pos -> true # define where is the bulk region +parameters = Dict("parameters" => + Dict("bhl" => 1u"ms", "bel" => 1"ms", "shl" => 10u"ns", "sel" => 10"ns", "if_in_bulk" => if_in_bulk)) +sim.detector = SolidStateDetector(sim.detector, BulkSurfaceConstantLifetimeChargeTrappingModel{T}(parameters)) +``` + ## Group Effects diff --git a/src/ChargeTrapping/ChargeTrapping.jl b/src/ChargeTrapping/ChargeTrapping.jl index b26513692..ab5bb1c25 100644 --- a/src/ChargeTrapping/ChargeTrapping.jl +++ b/src/ChargeTrapping/ChargeTrapping.jl @@ -121,3 +121,78 @@ end + +""" + struct BulkSurfaceConstantLifetimeChargeTrappingModel{T <: SSDFloat} <: AbstractChargeTrappingModel{T} + +Charge trapping model based on the BoggsChargeTrappingModel, but use the holelifetime instead. And different in bulk/surface + +## Fields +* `bhl::T`: Lifetime for holes in bulk volume (sensitive region) (default: `bhl = 1ms`). +* `bel::T`: Lifetime for electrons in bulk volume (default: `bel = 1ms`). +* `shl::T`: Lifetime for holes in surface layer (dead layer/inactive layer) (default: `shl = 1ms`). +* `sel::T`: Lifetime for electrons in surface layer (default: `sel = 1ms`). +* `if_in_bulk::T`: Function to decide if the carrier is in the bulk volume (default: `(pos::CartesianPoint{T}) -> true`). +""" +struct BulkSurfaceConstantLifetimeChargeTrappingModel{T <: SSDFloat} <: AbstractChargeTrappingModel{T} + bhl::T + bel::T + shl::T + sel::T + if_in_bulk::Function +end + +function _calculate_signal( + ctm::BulkSurfaceConstantLifetimeChargeTrappingModel{T}, + path::AbstractVector{CartesianPoint{T}}, + pathtimestamps::AbstractVector{T}, + charge::T, + wpot::Interpolations.Extrapolation{T, 3}, + S::CoordinateSystemType + )::Vector{T} where {T <: SSDFloat} + + tmp_signal::Vector{T} = Vector{T}(undef, length(pathtimestamps)) + + running_sum::T = zero(T) + q::T = charge + bl::T = ifelse(charge > 0, ctm.bhl, ctm.bel) + sl::T = ifelse(charge > 0, ctm.shl, ctm.sel) + + @inbounds for i in eachindex(tmp_signal) + Δt::T = (i > 1) ? (pathtimestamps[i] - pathtimestamps[i-1]) : zero(T) + Δq::T = ctm.if_in_bulk(path[i]) ? q * Δt/bl : q * Δt/sl + q -= Δq + w::T = i > 1 ? get_interpolation(wpot, path[i], S) : zero(T) + running_sum += w * Δq + tmp_signal[i] = running_sum + w * q + end + + tmp_signal +end + + +BulkSurfaceConstantLifetimeChargeTrappingModel(args...; T::Type{<:SSDFloat}, kwargs...) = BulkSurfaceConstantLifetimeChargeTrappingModel{T}(args...; kwargs...) +function BulkSurfaceConstantLifetimeChargeTrappingModel{T}(config_dict::AbstractDict = Dict()) where {T <: SSDFloat} + bhl::T = ustrip(u"s", 1u"ms") + bel::T = ustrip(u"s", 1u"ms") + shl::T = ustrip(u"s", 1u"ms") + sel::T = ustrip(u"s", 1u"ms") + if_in_bulk::Function = (pos::CartesianPoint{T}) -> true + + if haskey(config_dict, "model") && !haskey(config_dict, "parameters") + throw(ConfigFileError("`BulkSurfaceConstantLifetimeChargeTrappingModel` does not have `parameters`")) + end + + parameters = haskey(config_dict, "parameters") ? config_dict["parameters"] : config_dict + + allowed_keys = ("bhl", "bel", "shl", "sel", "if_in_bulk") + k = filter(k -> !(k in allowed_keys), keys(parameters)) + !isempty(k) && @warn "The following keys will be ignored: $(k).\nAllowed keys are: $(allowed_keys)" + + if haskey(parameters, "bhl") bhl = _parse_value(T, parameters["bhl"], internal_time_unit) end + if haskey(parameters, "bel") bel = _parse_value(T, parameters["bel"], internal_time_unit) end + if haskey(parameters, "shl") shl = _parse_value(T, parameters["shl"], internal_time_unit) end + if haskey(parameters, "sel") sel = _parse_value(T, parameters["sel"], internal_time_unit) end + if haskey(parameters, "if_in_bulk") if_in_bulk = parameters["if_in_bulk"] end + BulkSurfaceConstantLifetimeChargeTrappingModel{T}(bhl, bel, shl, sel, if_in_bulk) +end \ No newline at end of file diff --git a/src/SolidStateDetector/Semiconductor.jl b/src/SolidStateDetector/Semiconductor.jl index adb536089..743668651 100644 --- a/src/SolidStateDetector/Semiconductor.jl +++ b/src/SolidStateDetector/Semiconductor.jl @@ -85,6 +85,11 @@ function Semiconductor{T}(dict::AbstractDict, input_units::NamedTuple, outer_tra haskey(dict["charge_trapping_model"], "model") && dict["charge_trapping_model"]["model"] == "Boggs" BoggsChargeTrappingModel{T}(dict["charge_trapping_model"], temperature = temperature) + # TODO: implementation method of the BulkSurfaceConstantLifetimeChargeTrappingModel using the Configuration file + elseif haskey(dict, "charge_trapping_model") && + haskey(dict["charge_trapping_model"], "model") && + dict["charge_trapping_model"]["model"] == "BulkSurfaceConstantLifetime" + BulkSurfaceConstantLifetimeChargeTrappingModel{T}(dict["charge_trapping_model"]) else NoChargeTrappingModel{T}() end diff --git a/src/SolidStateDetectors.jl b/src/SolidStateDetectors.jl index 092971631..a22832db8 100644 --- a/src/SolidStateDetectors.jl +++ b/src/SolidStateDetectors.jl @@ -69,6 +69,7 @@ export apply_initial_state! export calculate_electric_potential!, calculate_weighting_potential!, calculate_electric_field!, calculate_drift_fields! export ElectricFieldChargeDriftModel, ADLChargeDriftModel, IsotropicChargeDriftModel export NoChargeTrappingModel, BoggsChargeTrappingModel +export BulkSurfaceConstantLifetimeChargeTrappingModel export get_active_volume, is_depleted, estimate_depletion_voltage export calculate_stored_energy, calculate_mutual_capacitance, calculate_capacitance_matrix export simulate_waveforms diff --git a/test/test_charge_drift_models.jl b/test/test_charge_drift_models.jl index 58d66806b..3ae9c1bf5 100644 --- a/test/test_charge_drift_models.jl +++ b/test/test_charge_drift_models.jl @@ -31,7 +31,7 @@ nbcc = NBodyChargeCloud(pos, Edep, 40, radius = T(0.0005), number_of_shells = 2) @test isapprox( signalsum, T(2), atol = 5e-3 ) end -@timed_testset "Charge Trapping" begin +@timed_testset "Charge Trapping: BoggsChargeTrappingModel" begin sim.detector = SolidStateDetector(sim.detector, BoggsChargeTrappingModel{T}()) evt = Event(pos, Edep) timed_simulate!(evt, sim) @@ -80,6 +80,21 @@ end end end +@timed_testset "Charge Trapping: BulkSurfaceConstantLifetimeChargeTrappingModel" begin + bhl = bel = shl = sel = 0.1 # 0.1 ms -> signalsum: ~1.9778; 0.01 -> ~1.8029; 1 -> ~1.9963 + parameters = Dict("parameters" => Dict("bhl" => bhl*u"ms", "bel" => bel*u"ms", "shl" => shl*u"ms", "sel" => sel*u"ms")) + trapping_model=BulkSurfaceConstantLifetimeChargeTrappingModel{T}(parameters) + sim.detector = SolidStateDetector(sim.detector, trapping_model) + evt = Event(pos, Edep) + timed_simulate!(evt, sim) + signalsum = T(0) + for i in 1:length(evt.waveforms) + signalsum += abs(ustrip(evt.waveforms[i].signal[end])) + end + signalsum *= inv(ustrip(SolidStateDetectors._convert_internal_energy_to_external_charge(sim.detector.semiconductor.material))) + @info signalsum + @test signalsum < T(1.99) +end @timed_testset "Test completeness of charge drift models" begin for c in InteractiveUtils.subtypes(AbstractChargeDriftModel) From b8bfdee0ee4d787bf8804df78fdfeaf345413abd Mon Sep 17 00:00:00 2001 From: zpdg Date: Mon, 21 Apr 2025 23:01:25 +0200 Subject: [PATCH 2/7] Fix the problem: the short doi failed the github check --- docs/src/man/charge_drift.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/docs/src/man/charge_drift.md b/docs/src/man/charge_drift.md index 46bd90b19..90dfc522b 100644 --- a/docs/src/man/charge_drift.md +++ b/docs/src/man/charge_drift.md @@ -263,7 +263,7 @@ sim.detector = SolidStateDetector(sim.detector, BoggsChargeTrappingModel{T}(para ### `BulkSurfaceConstantLifetimeTrappingModel` This constant-lifetime-based charge trapping model is similar to the Boggs model, which is constant-mean-free-path based. -This idea is presented in [Dai _et al._ (2023)](10.1016/j.apradiso.2022.110638). +This idea is presented in [Dai _et al._ (2023)](https://doi.org/10.1016/j.apradiso.2022.110638). Besides, the model implemented has two sets of parameters for the bulk (sensitive region) and surface (dead layer/inactive layer) volume respectively. It is implemented mainly for the pulse shape simulation of the surface layer. From d158867f0dd739eb3f493d8a3caeca0465ca5438 Mon Sep 17 00:00:00 2001 From: zpdg Date: Thu, 24 Apr 2025 15:12:49 +0200 Subject: [PATCH 3/7] 1) Support for modeling with config (define inactive volume geometry within cofig) 2) use symbol 3) more test --- docs/src/man/charge_drift.md | 57 +++++++++++---- .../public_TrueCoaxial.yaml | 31 +++++--- src/ChargeTrapping/ChargeTrapping.jl | 71 ++++++++++--------- src/SolidStateDetector/Semiconductor.jl | 16 +++-- src/SolidStateDetectors.jl | 3 +- test/test_charge_drift_models.jl | 49 +++++++++---- 6 files changed, 153 insertions(+), 74 deletions(-) diff --git a/docs/src/man/charge_drift.md b/docs/src/man/charge_drift.md index 90dfc522b..a4e85db59 100644 --- a/docs/src/man/charge_drift.md +++ b/docs/src/man/charge_drift.md @@ -261,26 +261,57 @@ parameters = Dict("parameters" => sim.detector = SolidStateDetector(sim.detector, BoggsChargeTrappingModel{T}(parameters)) ``` -### `BulkSurfaceConstantLifetimeTrappingModel` -This constant-lifetime-based charge trapping model is similar to the Boggs model, which is constant-mean-free-path based. +### `ConstantLifetimeTrappingModel` +This constant-lifetime-based charge trapping model is similar to the BoggsChargeTrappingModel, which is constant-mean-free-path based. +In this model, the charge cloud loses part of its charge at every point `path[i]` of the charge drift, depending on the lifetime. +The charge signal is then given by the charge-decreased charge cloud reaching the contacts and the charges trapped on the way. + This idea is presented in [Dai _et al._ (2023)](https://doi.org/10.1016/j.apradiso.2022.110638). -Besides, the model implemented has two sets of parameters for the bulk (sensitive region) and surface (dead layer/inactive layer) volume respectively. -It is implemented mainly for the pulse shape simulation of the surface layer. -In this model, the charge cloud loses part of its charge at every point `path[i]` of the charge drift, depending on the lifetime. +Besides, the model implemented has two sets of parameters for the sensitive (bulk) and inactive (dead layer/surface layer) volume respectively. +The inactive layer effect will only be considered if the corresponding `virtual_drift_volume` is defined in the configuration file. -The charge signal is then given by the charge-decreased charge cloud reaching the contacts and the charges trapped on the way. +The `ConstantLifetimeTrappingModel` can be applied in the configuration file by adding a field `charge_trapping_model` to the `semiconductor` with `model: ConstantLifetime`, `parameters` defining the constant lifetime and `inactive_layer_geometry`(optional) defining the inactive volume: +```yaml +detectors: + - semiconductor: + material: #... + geometry: #... + charge_trapping_model: + model: ConstantLifetime + parameters: + τh: 1ms + τe: 1ms + τh_inactive: 80ns + τe_inactive: 80ns + inactive_layer_geometry: + tube: + r: + from: 9.0 + to: 10.0 + h: 10.0 + origin: + z: 5.0 +``` -It can be applied to an already read-in `SolidStateDetector` using for example: +The inactive layer could only be modeled with the help of the configuration file. +But for the sensitive volume, the parameters can always be applied to an already read-in `SolidStateDetector` using for example: ```julia using SolidStateDetectors, Unitful -T = Float32 -sim = Simulation{T}(SSD_examples[:InvertedCoax]) +T = Float64 +sim = Simulation{T}(SSD_examples[:TrueCoaxial]) +simulate!(sim) -if_in_bulk = pos -> true # define where is the bulk region -parameters = Dict("parameters" => - Dict("bhl" => 1u"ms", "bel" => 1"ms", "shl" => 10u"ns", "sel" => 10"ns", "if_in_bulk" => if_in_bulk)) -sim.detector = SolidStateDetector(sim.detector, BulkSurfaceConstantLifetimeChargeTrappingModel{T}(parameters)) +# only model the sensitive volume +τh = τe = 1e6 +parameters = Dict("parameters" => Dict("τh" => τh*u"ns", "τe" => τe*u"ns")) + +# model both the sensitive and inactive volumes based on the inactive volume geometry defined in the configuration file +τh = τe = 1e6 +τh_inactive = τe_inactive = 80 +parameters = Dict("parameters" => Dict("τh" => τh*u"ns", "τe" => τe*u"ns", "τh_inactive" => τh_inactive*u"ns", "τe_inactive" => τe_inactive*u"ns"), "inactive_layer_geometry" => sim.detector.semiconductor.charge_trapping_model.inactive_layer_geometry) + +sim.detector = SolidStateDetector(sim.detector, ConstantLifetimeChargeTrappingModel{T}(parameters)) ``` diff --git a/examples/example_config_files/public_TrueCoaxial.yaml b/examples/example_config_files/public_TrueCoaxial.yaml index 2550881cf..4976e1018 100644 --- a/examples/example_config_files/public_TrueCoaxial.yaml +++ b/examples/example_config_files/public_TrueCoaxial.yaml @@ -32,14 +32,29 @@ detectors: value: 1e8cm^-3 charge_drift_model: include: ADLChargeDriftModel/drift_velocity_config.yaml + charge_trapping_model: + model: ConstantLifetime + parameters: + τh: 1ms + τe: 1ms + τh_inactive: 80ns + τe_inactive: 80ns + inactive_layer_geometry: + tube: + r: + from: 9.0 + to: 10.0 + h: 10.0 + origin: + z: 5.0 geometry: cone: r: from: 1.0 to: 10.0 - h: 10 + h: 10.0 origin: - z: 5 + z: 5.0 contacts: - name: "P+" id: 1 @@ -49,10 +64,10 @@ detectors: tube: r: from: 0.5 - to: 1 - h: 10 + to: 1.0 + h: 10.0 origin: - z: 5 + z: 5.0 - name: "N+" id: 2 material: HPGe @@ -61,7 +76,7 @@ detectors: tube: r: from: 10.0 - to: 10.5 - h: 10. + to: 10.0 + h: 10.0 origin: - z: 5 \ No newline at end of file + z: 5.0 \ No newline at end of file diff --git a/src/ChargeTrapping/ChargeTrapping.jl b/src/ChargeTrapping/ChargeTrapping.jl index ab5bb1c25..866cfc750 100644 --- a/src/ChargeTrapping/ChargeTrapping.jl +++ b/src/ChargeTrapping/ChargeTrapping.jl @@ -123,27 +123,33 @@ end """ - struct BulkSurfaceConstantLifetimeChargeTrappingModel{T <: SSDFloat} <: AbstractChargeTrappingModel{T} + struct ConstantLifetimeChargeTrappingModel{T <: SSDFloat} <: AbstractChargeTrappingModel{T} -Charge trapping model based on the BoggsChargeTrappingModel, but use the holelifetime instead. And different in bulk/surface +This constant-lifetime-based charge trapping model is similar to the Boggs model, which is constant-mean-free-path based. + +The model implemented has two sets of parameters for the sensitive (bulk) and inactive (dead layer/surface layer) volume respectively. +The inactive layer effect will only be considered if the corresponding `virtual_drift_volume` is defined in the configuration file. ## Fields -* `bhl::T`: Lifetime for holes in bulk volume (sensitive region) (default: `bhl = 1ms`). -* `bel::T`: Lifetime for electrons in bulk volume (default: `bel = 1ms`). -* `shl::T`: Lifetime for holes in surface layer (dead layer/inactive layer) (default: `shl = 1ms`). -* `sel::T`: Lifetime for electrons in surface layer (default: `sel = 1ms`). -* `if_in_bulk::T`: Function to decide if the carrier is in the bulk volume (default: `(pos::CartesianPoint{T}) -> true`). +* `τh::T`: Lifetime for holes in bulk volume (sensitive region) (default: `τh = 1ms`). +* `τe::T`: Lifetime for electrons in bulk volume (default: `τe = 1ms`). +* `τh_inactive::T`: Lifetime for holes in surface layer (dead layer/inactive layer) (default: `τh_inactive = 1ms`). +* `τe_inactive::T`: Lifetime for electrons in surface layer (default: `τe_inactive = 1ms`). +* `inactive_layer_geometry`: The geometry of the inactive layer. When this is defined, the model will consider the inactive effect. +> Note: All τ must be much bigger than `Δt`. + +See also [Charge Trapping Models](@ref). """ -struct BulkSurfaceConstantLifetimeChargeTrappingModel{T <: SSDFloat} <: AbstractChargeTrappingModel{T} - bhl::T - bel::T - shl::T - sel::T - if_in_bulk::Function +struct ConstantLifetimeChargeTrappingModel{T <: SSDFloat} <: AbstractChargeTrappingModel{T} + τh::T + τe::T + τh_inactive::T + τe_inactive::T + inactive_layer_geometry end function _calculate_signal( - ctm::BulkSurfaceConstantLifetimeChargeTrappingModel{T}, + ctm::ConstantLifetimeChargeTrappingModel{T}, path::AbstractVector{CartesianPoint{T}}, pathtimestamps::AbstractVector{T}, charge::T, @@ -155,12 +161,14 @@ function _calculate_signal( running_sum::T = zero(T) q::T = charge - bl::T = ifelse(charge > 0, ctm.bhl, ctm.bel) - sl::T = ifelse(charge > 0, ctm.shl, ctm.sel) + + inactive_layer_exist = ctm.inactive_layer_geometry !== missing + τ::T = ifelse(charge > 0, ctm.τh, ctm.τe) + τ_inactive::T = ifelse(charge > 0, ctm.τh_inactive, ctm.τe_inactive) @inbounds for i in eachindex(tmp_signal) Δt::T = (i > 1) ? (pathtimestamps[i] - pathtimestamps[i-1]) : zero(T) - Δq::T = ctm.if_in_bulk(path[i]) ? q * Δt/bl : q * Δt/sl + Δq::T = (inactive_layer_exist && in(path[i], ctm.inactive_layer_geometry)) ? q * Δt/τ_inactive : q * Δt/τ q -= Δq w::T = i > 1 ? get_interpolation(wpot, path[i], S) : zero(T) running_sum += w * Δq @@ -171,28 +179,27 @@ function _calculate_signal( end -BulkSurfaceConstantLifetimeChargeTrappingModel(args...; T::Type{<:SSDFloat}, kwargs...) = BulkSurfaceConstantLifetimeChargeTrappingModel{T}(args...; kwargs...) -function BulkSurfaceConstantLifetimeChargeTrappingModel{T}(config_dict::AbstractDict = Dict()) where {T <: SSDFloat} - bhl::T = ustrip(u"s", 1u"ms") - bel::T = ustrip(u"s", 1u"ms") - shl::T = ustrip(u"s", 1u"ms") - sel::T = ustrip(u"s", 1u"ms") - if_in_bulk::Function = (pos::CartesianPoint{T}) -> true +ConstantLifetimeChargeTrappingModel(args...; T::Type{<:SSDFloat}, kwargs...) = ConstantLifetimeChargeTrappingModel{T}(args...; kwargs...) +function ConstantLifetimeChargeTrappingModel{T}(config_dict::AbstractDict = Dict()) where {T <: SSDFloat} + τh::T = ustrip(u"s", 1u"ms") + τe::T = ustrip(u"s", 1u"ms") + τh_inactive::T = ustrip(u"s", 80u"ns") + τe_inactive::T = ustrip(u"s", 80u"ns") if haskey(config_dict, "model") && !haskey(config_dict, "parameters") - throw(ConfigFileError("`BulkSurfaceConstantLifetimeChargeTrappingModel` does not have `parameters`")) + throw(ConfigFileError("`ConstantLifetimeChargeTrappingModel` does not have `parameters`")) end parameters = haskey(config_dict, "parameters") ? config_dict["parameters"] : config_dict + inactive_layer_geometry = haskey(config_dict, "inactive_layer_geometry") ? config_dict["inactive_layer_geometry"] : missing - allowed_keys = ("bhl", "bel", "shl", "sel", "if_in_bulk") + allowed_keys = ("τh", "τe", "τh_inactive", "τe_inactive") k = filter(k -> !(k in allowed_keys), keys(parameters)) !isempty(k) && @warn "The following keys will be ignored: $(k).\nAllowed keys are: $(allowed_keys)" - if haskey(parameters, "bhl") bhl = _parse_value(T, parameters["bhl"], internal_time_unit) end - if haskey(parameters, "bel") bel = _parse_value(T, parameters["bel"], internal_time_unit) end - if haskey(parameters, "shl") shl = _parse_value(T, parameters["shl"], internal_time_unit) end - if haskey(parameters, "sel") sel = _parse_value(T, parameters["sel"], internal_time_unit) end - if haskey(parameters, "if_in_bulk") if_in_bulk = parameters["if_in_bulk"] end - BulkSurfaceConstantLifetimeChargeTrappingModel{T}(bhl, bel, shl, sel, if_in_bulk) + if haskey(parameters, "τh") τh = _parse_value(T, parameters["τh"], internal_time_unit) end + if haskey(parameters, "τe") τe = _parse_value(T, parameters["τe"], internal_time_unit) end + if haskey(parameters, "τh_inactive") τh_inactive = _parse_value(T, parameters["τh_inactive"], internal_time_unit) end + if haskey(parameters, "τe_inactive") τe_inactive = _parse_value(T, parameters["τe_inactive"], internal_time_unit) end + ConstantLifetimeChargeTrappingModel{T}(τh, τe, τh_inactive, τe_inactive, inactive_layer_geometry) end \ No newline at end of file diff --git a/src/SolidStateDetector/Semiconductor.jl b/src/SolidStateDetector/Semiconductor.jl index 743668651..a90932e16 100644 --- a/src/SolidStateDetector/Semiconductor.jl +++ b/src/SolidStateDetector/Semiconductor.jl @@ -81,22 +81,26 @@ function Semiconductor{T}(dict::AbstractDict, input_units::NamedTuple, outer_tra T(293) end + inner_transformations = parse_CSG_transformation(T, dict, input_units) + transformations = combine_transformations(inner_transformations, outer_transformations) + geometry = Geometry(T, dict["geometry"], input_units, transformations) + charge_trapping_model = if haskey(dict, "charge_trapping_model") && haskey(dict["charge_trapping_model"], "model") && dict["charge_trapping_model"]["model"] == "Boggs" BoggsChargeTrappingModel{T}(dict["charge_trapping_model"], temperature = temperature) - # TODO: implementation method of the BulkSurfaceConstantLifetimeChargeTrappingModel using the Configuration file elseif haskey(dict, "charge_trapping_model") && haskey(dict["charge_trapping_model"], "model") && - dict["charge_trapping_model"]["model"] == "BulkSurfaceConstantLifetime" - BulkSurfaceConstantLifetimeChargeTrappingModel{T}(dict["charge_trapping_model"]) + dict["charge_trapping_model"]["model"] == "ConstantLifetime" + constant_lifetime_ctm_dict = deepcopy(dict["charge_trapping_model"]) + if haskey(constant_lifetime_ctm_dict, "inactive_layer_geometry") + constant_lifetime_ctm_dict["inactive_layer_geometry"] = Geometry(T, dict["charge_trapping_model"]["inactive_layer_geometry"], input_units, transformations) + end + ConstantLifetimeChargeTrappingModel{T}(constant_lifetime_ctm_dict) else NoChargeTrappingModel{T}() end - inner_transformations = parse_CSG_transformation(T, dict, input_units) - transformations = combine_transformations(inner_transformations, outer_transformations) - geometry = Geometry(T, dict["geometry"], input_units, transformations) return Semiconductor(temperature, material, impurity_density_model, charge_drift_model, charge_trapping_model, geometry) end diff --git a/src/SolidStateDetectors.jl b/src/SolidStateDetectors.jl index a22832db8..0e9bf5128 100644 --- a/src/SolidStateDetectors.jl +++ b/src/SolidStateDetectors.jl @@ -68,8 +68,7 @@ export ElectricPotential, PointTypes, EffectiveChargeDensity, DielectricDistribu export apply_initial_state! export calculate_electric_potential!, calculate_weighting_potential!, calculate_electric_field!, calculate_drift_fields! export ElectricFieldChargeDriftModel, ADLChargeDriftModel, IsotropicChargeDriftModel -export NoChargeTrappingModel, BoggsChargeTrappingModel -export BulkSurfaceConstantLifetimeChargeTrappingModel +export NoChargeTrappingModel, BoggsChargeTrappingModel, ConstantLifetimeChargeTrappingModel export get_active_volume, is_depleted, estimate_depletion_voltage export calculate_stored_energy, calculate_mutual_capacitance, calculate_capacitance_matrix export simulate_waveforms diff --git a/test/test_charge_drift_models.jl b/test/test_charge_drift_models.jl index 3ae9c1bf5..33f7dd39a 100644 --- a/test/test_charge_drift_models.jl +++ b/test/test_charge_drift_models.jl @@ -80,20 +80,43 @@ end end end -@timed_testset "Charge Trapping: BulkSurfaceConstantLifetimeChargeTrappingModel" begin - bhl = bel = shl = sel = 0.1 # 0.1 ms -> signalsum: ~1.9778; 0.01 -> ~1.8029; 1 -> ~1.9963 - parameters = Dict("parameters" => Dict("bhl" => bhl*u"ms", "bel" => bel*u"ms", "shl" => shl*u"ms", "sel" => sel*u"ms")) - trapping_model=BulkSurfaceConstantLifetimeChargeTrappingModel{T}(parameters) - sim.detector = SolidStateDetector(sim.detector, trapping_model) - evt = Event(pos, Edep) - timed_simulate!(evt, sim) - signalsum = T(0) - for i in 1:length(evt.waveforms) - signalsum += abs(ustrip(evt.waveforms[i].signal[end])) + +@timed_testset "Charge Trapping: ConstantLifetimeChargeTrappingModel" begin + signalsum_list = [] + lifetime_list = 10 .^ (-3:0.3:0) + for lifetime in lifetime_list + parameters = Dict("parameters" => Dict("τh" => lifetime*u"ms", "τe" => lifetime*u"ms")) + trapping_model=ConstantLifetimeChargeTrappingModel{T}(parameters) + sim.detector = SolidStateDetector(sim.detector, trapping_model) + evt = Event(pos, Edep) + timed_simulate!(evt, sim) + signalsum = T(0) + for i in 1:length(evt.waveforms) + signalsum += abs(ustrip(evt.waveforms[i].signal[end])) + end + signalsum *= inv(ustrip(SolidStateDetectors._convert_internal_energy_to_external_charge(sim.detector.semiconductor.material))) + push!(signalsum_list, signalsum) end - signalsum *= inv(ustrip(SolidStateDetectors._convert_internal_energy_to_external_charge(sim.detector.semiconductor.material))) - @info signalsum - @test signalsum < T(1.99) + + @info signalsum_list + @test all(signalsum .< T(2.0)) + @test all(diff(signalsum_list) .> 0) + + config_dict = SolidStateDetectors.parse_config_file(SSD_examples[:TrueCoaxial]) + @testset "Parse config file 1" begin + config_dict["detectors"][1]["semiconductor"]["charge_trapping_model"] = Dict( + "model" => "ConstantLifetime", + "parameters" => Dict( + "τh" => "1ms", + "τe" => "1ms", + ) + ) + simA = @test_nowarn Simulation{T}(config_dict) + @test simA.detector.semiconductor.charge_trapping_model isa ConstantLifetimeChargeTrappingModel{T} + @test simA.detector.semiconductor.charge_trapping_model.τh == T(0.001) + @test simA.detector.semiconductor.charge_trapping_model.τe == T(0.001) + end + end @timed_testset "Test completeness of charge drift models" begin From dc592cffd09b1c0f8bdd8576006de47d3e73040f Mon Sep 17 00:00:00 2001 From: zpdg Date: Thu, 24 Apr 2025 15:31:09 +0200 Subject: [PATCH 4/7] Fix a typo --- test/test_charge_drift_models.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/test_charge_drift_models.jl b/test/test_charge_drift_models.jl index 33f7dd39a..2609631f6 100644 --- a/test/test_charge_drift_models.jl +++ b/test/test_charge_drift_models.jl @@ -99,7 +99,7 @@ end end @info signalsum_list - @test all(signalsum .< T(2.0)) + @test all(signalsum_list .< T(2.0)) @test all(diff(signalsum_list) .> 0) config_dict = SolidStateDetectors.parse_config_file(SSD_examples[:TrueCoaxial]) From 35c6a7f7d9abe8147765ae3f8b472a98d96fd8a6 Mon Sep 17 00:00:00 2001 From: zpdg Date: Tue, 29 Apr 2025 20:06:38 +0200 Subject: [PATCH 5/7] Clearer naming --- src/ChargeTrapping/ChargeTrapping.jl | 17 +++++++++-------- 1 file changed, 9 insertions(+), 8 deletions(-) diff --git a/src/ChargeTrapping/ChargeTrapping.jl b/src/ChargeTrapping/ChargeTrapping.jl index 866cfc750..ee914d9a5 100644 --- a/src/ChargeTrapping/ChargeTrapping.jl +++ b/src/ChargeTrapping/ChargeTrapping.jl @@ -128,7 +128,7 @@ end This constant-lifetime-based charge trapping model is similar to the Boggs model, which is constant-mean-free-path based. The model implemented has two sets of parameters for the sensitive (bulk) and inactive (dead layer/surface layer) volume respectively. -The inactive layer effect will only be considered if the corresponding `virtual_drift_volume` is defined in the configuration file. +The inactive layer effect will only be considered if the corresponding `inactive_layer_geometry` is defined in the configuration file. ## Fields * `τh::T`: Lifetime for holes in bulk volume (sensitive region) (default: `τh = 1ms`). @@ -140,12 +140,12 @@ The inactive layer effect will only be considered if the corresponding `virtual_ See also [Charge Trapping Models](@ref). """ -struct ConstantLifetimeChargeTrappingModel{T <: SSDFloat} <: AbstractChargeTrappingModel{T} +struct ConstantLifetimeChargeTrappingModel{T <: SSDFloat, G <: Union{<:AbstractGeometry, Nothing}} <: AbstractChargeTrappingModel{T} τh::T τe::T τh_inactive::T τe_inactive::T - inactive_layer_geometry + inactive_layer_geometry::G end function _calculate_signal( @@ -162,13 +162,14 @@ function _calculate_signal( running_sum::T = zero(T) q::T = charge - inactive_layer_exist = ctm.inactive_layer_geometry !== missing + inactive_layer_exist = !isnothing(ctm.inactive_layer_geometry) τ::T = ifelse(charge > 0, ctm.τh, ctm.τe) τ_inactive::T = ifelse(charge > 0, ctm.τh_inactive, ctm.τe_inactive) @inbounds for i in eachindex(tmp_signal) Δt::T = (i > 1) ? (pathtimestamps[i] - pathtimestamps[i-1]) : zero(T) - Δq::T = (inactive_layer_exist && in(path[i], ctm.inactive_layer_geometry)) ? q * Δt/τ_inactive : q * Δt/τ + τi::T = (inactive_layer_exist && in(path[i], ctm.inactive_layer_geometry)) ? τ_inactive : τ + Δq::T = q * Δt/τi q -= Δq w::T = i > 1 ? get_interpolation(wpot, path[i], S) : zero(T) running_sum += w * Δq @@ -184,14 +185,14 @@ function ConstantLifetimeChargeTrappingModel{T}(config_dict::AbstractDict = Dict τh::T = ustrip(u"s", 1u"ms") τe::T = ustrip(u"s", 1u"ms") τh_inactive::T = ustrip(u"s", 80u"ns") - τe_inactive::T = ustrip(u"s", 80u"ns") + τe_inactive::T = ustrip(u"s", 80u"ns") if haskey(config_dict, "model") && !haskey(config_dict, "parameters") throw(ConfigFileError("`ConstantLifetimeChargeTrappingModel` does not have `parameters`")) end parameters = haskey(config_dict, "parameters") ? config_dict["parameters"] : config_dict - inactive_layer_geometry = haskey(config_dict, "inactive_layer_geometry") ? config_dict["inactive_layer_geometry"] : missing + inactive_layer_geometry = haskey(config_dict, "inactive_layer_geometry") ? config_dict["inactive_layer_geometry"] : nothing allowed_keys = ("τh", "τe", "τh_inactive", "τe_inactive") k = filter(k -> !(k in allowed_keys), keys(parameters)) @@ -201,5 +202,5 @@ function ConstantLifetimeChargeTrappingModel{T}(config_dict::AbstractDict = Dict if haskey(parameters, "τe") τe = _parse_value(T, parameters["τe"], internal_time_unit) end if haskey(parameters, "τh_inactive") τh_inactive = _parse_value(T, parameters["τh_inactive"], internal_time_unit) end if haskey(parameters, "τe_inactive") τe_inactive = _parse_value(T, parameters["τe_inactive"], internal_time_unit) end - ConstantLifetimeChargeTrappingModel{T}(τh, τe, τh_inactive, τe_inactive, inactive_layer_geometry) + ConstantLifetimeChargeTrappingModel{T, typeof(inactive_layer_geometry)}(τh, τe, τh_inactive, τe_inactive, inactive_layer_geometry) end \ No newline at end of file From 7c027130d2ffea319a8c5d365d636856815da192 Mon Sep 17 00:00:00 2001 From: zpdg Date: Tue, 29 Apr 2025 20:07:35 +0200 Subject: [PATCH 6/7] Clearer naming --- docs/src/man/charge_drift.md | 16 ++++++++-------- 1 file changed, 8 insertions(+), 8 deletions(-) diff --git a/docs/src/man/charge_drift.md b/docs/src/man/charge_drift.md index a4e85db59..b02660ba7 100644 --- a/docs/src/man/charge_drift.md +++ b/docs/src/man/charge_drift.md @@ -261,8 +261,8 @@ parameters = Dict("parameters" => sim.detector = SolidStateDetector(sim.detector, BoggsChargeTrappingModel{T}(parameters)) ``` -### `ConstantLifetimeTrappingModel` -This constant-lifetime-based charge trapping model is similar to the BoggsChargeTrappingModel, which is constant-mean-free-path based. +### `ConstantLifetimeChargeTrappingModel` +This constant-lifetime-based charge trapping model assumes electrons and holes to have a constant free lifetime throughout the crystal. In this model, the charge cloud loses part of its charge at every point `path[i]` of the charge drift, depending on the lifetime. The charge signal is then given by the charge-decreased charge cloud reaching the contacts and the charges trapped on the way. @@ -271,7 +271,7 @@ This idea is presented in [Dai _et al._ (2023)](https://doi.org/10.1016/j.apradi Besides, the model implemented has two sets of parameters for the sensitive (bulk) and inactive (dead layer/surface layer) volume respectively. The inactive layer effect will only be considered if the corresponding `virtual_drift_volume` is defined in the configuration file. -The `ConstantLifetimeTrappingModel` can be applied in the configuration file by adding a field `charge_trapping_model` to the `semiconductor` with `model: ConstantLifetime`, `parameters` defining the constant lifetime and `inactive_layer_geometry`(optional) defining the inactive volume: +The `ConstantLifetimeChargeTrappingModel` can be applied in the configuration file by adding a field `charge_trapping_model` to the `semiconductor` with `model: ConstantLifetime`, `parameters` defining the constant lifetime and `inactive_layer_geometry`(optional) defining the inactive volume: ```yaml detectors: - semiconductor: @@ -303,13 +303,13 @@ sim = Simulation{T}(SSD_examples[:TrueCoaxial]) simulate!(sim) # only model the sensitive volume -τh = τe = 1e6 -parameters = Dict("parameters" => Dict("τh" => τh*u"ns", "τe" => τe*u"ns")) +τh = τe = 1u"ms" +parameters = Dict("parameters" => Dict("τh" => τh, "τe" => τe)) # model both the sensitive and inactive volumes based on the inactive volume geometry defined in the configuration file -τh = τe = 1e6 -τh_inactive = τe_inactive = 80 -parameters = Dict("parameters" => Dict("τh" => τh*u"ns", "τe" => τe*u"ns", "τh_inactive" => τh_inactive*u"ns", "τe_inactive" => τe_inactive*u"ns"), "inactive_layer_geometry" => sim.detector.semiconductor.charge_trapping_model.inactive_layer_geometry) +τh = τe = 1u"ms" +τh_inactive = τe_inactive = 80u"ns" +parameters = Dict("parameters" => Dict("τh" => τh, "τe" => τe, "τh_inactive" => τh_inactive, "τe_inactive" => τe_inactive), "inactive_layer_geometry" => sim.detector.semiconductor.charge_trapping_model.inactive_layer_geometry) sim.detector = SolidStateDetector(sim.detector, ConstantLifetimeChargeTrappingModel{T}(parameters)) ``` From 9f276f44f6d46057fbfe38f8c4b0211ba66bf235 Mon Sep 17 00:00:00 2001 From: zpdg Date: Wed, 30 Apr 2025 09:20:52 +0200 Subject: [PATCH 7/7] Add test for inactive layer --- test/test_charge_drift_models.jl | 110 +++++++++++++++++++++++++------ 1 file changed, 89 insertions(+), 21 deletions(-) diff --git a/test/test_charge_drift_models.jl b/test/test_charge_drift_models.jl index 2609631f6..c2907e06d 100644 --- a/test/test_charge_drift_models.jl +++ b/test/test_charge_drift_models.jl @@ -82,41 +82,109 @@ end end @timed_testset "Charge Trapping: ConstantLifetimeChargeTrappingModel" begin + # test 1: parse the lifetimes and the inactive layer geometry + config_dict = SolidStateDetectors.parse_config_file(SSD_examples[:TrueCoaxial]) + simA = @test_nowarn Simulation{T}(config_dict) + @testset "Parse the TrueCoaxial config file" begin + @test simA.detector.semiconductor.charge_trapping_model isa ConstantLifetimeChargeTrappingModel{T} + @test simA.detector.semiconductor.charge_trapping_model.τh == T(1e-3) + @test simA.detector.semiconductor.charge_trapping_model.τe == T(1e-3) + @test simA.detector.semiconductor.charge_trapping_model.τh_inactive == T(8e-8) + @test simA.detector.semiconductor.charge_trapping_model.τe_inactive == T(8e-8) + @test simA.detector.semiconductor.charge_trapping_model.inactive_layer_geometry.origin == CartesianPoint{T}(0.0, 0.0, 0.005) + r0, r1 = T.((0.009, 0.01)) + @test simA.detector.semiconductor.charge_trapping_model.inactive_layer_geometry.r == tuple((r0, r1), (r0, r1)) + @test simA.detector.semiconductor.charge_trapping_model.inactive_layer_geometry.hZ == T(0.005) + end + + simA = Simulation{T}(SSD_examples[:TrueCoaxial]) + timed_simulate!(simA, convergence_limit = 1e-5, device_array_type = device_array_type, refinement_limits = [0.2, 0.1, 0.05], verbose = false) + simA_inactive_layer_geometry=deepcopy(simA.detector.semiconductor.charge_trapping_model.inactive_layer_geometry) + + # test2: only model the bulk volume, test the bulk signals while varying the lifetime: τ + pos = CartesianPoint{T}(0.005,0,0.005); Edep = 1u"eV" signalsum_list = [] - lifetime_list = 10 .^ (-3:0.3:0) - for lifetime in lifetime_list - parameters = Dict("parameters" => Dict("τh" => lifetime*u"ms", "τe" => lifetime*u"ms")) + τ_list = (10 .^ (-3:0.3:0))*u"ms" + for τ in τ_list + parameters = Dict("parameters" => Dict("τh" => τ, "τe" => τ)) trapping_model=ConstantLifetimeChargeTrappingModel{T}(parameters) - sim.detector = SolidStateDetector(sim.detector, trapping_model) + simA.detector = SolidStateDetector(simA.detector, trapping_model) evt = Event(pos, Edep) - timed_simulate!(evt, sim) + timed_simulate!(evt, simA) signalsum = T(0) for i in 1:length(evt.waveforms) signalsum += abs(ustrip(evt.waveforms[i].signal[end])) end - signalsum *= inv(ustrip(SolidStateDetectors._convert_internal_energy_to_external_charge(sim.detector.semiconductor.material))) + signalsum *= inv(ustrip(SolidStateDetectors._convert_internal_energy_to_external_charge(simA.detector.semiconductor.material))) push!(signalsum_list, signalsum) end - @info signalsum_list @test all(signalsum_list .< T(2.0)) @test all(diff(signalsum_list) .> 0) - config_dict = SolidStateDetectors.parse_config_file(SSD_examples[:TrueCoaxial]) - @testset "Parse config file 1" begin - config_dict["detectors"][1]["semiconductor"]["charge_trapping_model"] = Dict( - "model" => "ConstantLifetime", - "parameters" => Dict( - "τh" => "1ms", - "τe" => "1ms", - ) - ) - simA = @test_nowarn Simulation{T}(config_dict) - @test simA.detector.semiconductor.charge_trapping_model isa ConstantLifetimeChargeTrappingModel{T} - @test simA.detector.semiconductor.charge_trapping_model.τh == T(0.001) - @test simA.detector.semiconductor.charge_trapping_model.τe == T(0.001) + # test3: model both the bulk/inactive volumes + ## test 3.1: the bulk/inactive signals while varying the lifetimes: τ, τ_inactive + pos_bulk = CartesianPoint{T}(0.0085,0,0.005); Edep = 1u"eV" + @test !in(pos_bulk, simA_inactive_layer_geometry) + signalsum_list_bulk = [] + τ_list = (10 .^ (-2:0.2:0))*u"ms" + pos_inactive = CartesianPoint{T}(0.0095,0,0.005); Edep = 1u"eV" + @test in(pos_inactive, simA_inactive_layer_geometry) + signalsum_list_inactive = [] + τ_inactive_list = τ_list/100 + for (τ,τ_inactive) in zip(τ_list, τ_inactive_list) + parameters = Dict("parameters" => Dict("τh" => τ, "τe" => τ, "τh_inactive" => τ_inactive, "τe_inactive" => τ_inactive), + "inactive_layer_geometry" => simA_inactive_layer_geometry) + trapping_model=ConstantLifetimeChargeTrappingModel{T}(parameters) + simA.detector = SolidStateDetector(simA.detector, trapping_model) + evt_bulk = Event(pos_bulk , Edep) + timed_simulate!(evt_bulk, simA) + signalsum_bulk = T(0) + for i in 1:length(evt_bulk.waveforms) + signalsum_bulk += abs(ustrip(evt_bulk.waveforms[i].signal[end])) + end + signalsum_bulk *= inv(ustrip(SolidStateDetectors._convert_internal_energy_to_external_charge(simA.detector.semiconductor.material))) + push!(signalsum_list_bulk, signalsum_bulk) + + evt_inactive = Event(pos_inactive , Edep) + timed_simulate!(evt_inactive, simA) + signalsum_inactive = T(0) + for i in 1:length(evt_inactive.waveforms) + signalsum_inactive += abs(ustrip(evt_inactive.waveforms[i].signal[end])) + end + signalsum_inactive *= inv(ustrip(SolidStateDetectors._convert_internal_energy_to_external_charge(simA.detector.semiconductor.material))) + push!(signalsum_list_inactive, signalsum_inactive) end - + @info signalsum_list_bulk + @test all(signalsum_list_bulk .< T(2.0)) + @test all(diff(signalsum_list_bulk) .> 0) + @info signalsum_list_inactive + @test all(signalsum_list_inactive .< T(2.0)) + @test all(diff(signalsum_list_inactive) .> 0) + @test all(signalsum_list_bulk .> signalsum_list_inactive) + + ## test 3.2: the inactive layer signals while varying the depth + τ, τ_inactive = 1u"ms", 100u"ns" + parameters = Dict("parameters" => Dict("τh" => τ, "τe" => τ, "τh_inactive" => τ_inactive, "τe_inactive" => τ_inactive), + "inactive_layer_geometry" => simA_inactive_layer_geometry) + trapping_model=ConstantLifetimeChargeTrappingModel{T}(parameters) + simA.detector = SolidStateDetector(simA.detector, trapping_model) + signalsum_list_inactive = [] + for depth in (0.1:0.1:0.9)/1000 + pos_inactive = CartesianPoint{T}(0.01-depth,0,0.005); Edep = 1u"eV" + @test in(pos_inactive, simA_inactive_layer_geometry) + evt_inactive = Event(pos_inactive , Edep) + timed_simulate!(evt_inactive, simA, diffusion=false, Δt=1u"ns") + signalsum_inactive = T(0) + for i in 1:length(evt_inactive.waveforms) + signalsum_inactive += abs(ustrip(evt_inactive.waveforms[i].signal[end])) + end + signalsum_inactive *= inv(ustrip(SolidStateDetectors._convert_internal_energy_to_external_charge(simA.detector.semiconductor.material))) + push!(signalsum_list_inactive, signalsum_inactive) + end + @info signalsum_list_inactive + @test all(signalsum_list_inactive .< T(2.0)) + @test all(diff(signalsum_list_inactive) .> 0) end @timed_testset "Test completeness of charge drift models" begin