Skip to content

Commit 647f109

Browse files
feat(ehem): add EnforceLayer type and functor implementation
skip ci
1 parent ecabe68 commit 647f109

3 files changed

Lines changed: 172 additions & 0 deletions

File tree

src/engine/ehem/enforcelayer.jl

Lines changed: 136 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,136 @@
1+
"""
2+
$(TYPEDEF)
3+
4+
An EHEM formulation that creates a homogeneous earth model by enforcing the properties of a single, specified layer from a multi-layer model.
5+
6+
# Attributes
7+
$(TYPEDFIELDS)
8+
"""
9+
struct EnforceLayer <: AbstractEHEMFormulation
10+
"Index of the earth layer to enforce. `-1` selects the bottommost layer."
11+
layer::Int
12+
13+
@doc """
14+
$(TYPEDSIGNATURES)
15+
16+
Constructs an `EnforceLayer` instance.
17+
18+
# Arguments
19+
- `layer::Int`: The index of the layer to enforce.
20+
- `-1` (default): Enforces the properties of the bottommost earth layer.
21+
- `2`: Enforces the properties of the topmost earth layer (the one directly below the air).
22+
- `> 2`: Enforces the properties of a specific layer by its index.
23+
"""
24+
function EnforceLayer(; layer::Int = -1)
25+
@assert (layer == -1 || layer >= 2) "Invalid layer index. Must be -1 (bottommost) or >= 2."
26+
new(layer)
27+
end
28+
end
29+
30+
function get_description(f::EnforceLayer)
31+
if f.layer == -1
32+
return "Assume bottom layer"
33+
elseif f.layer == 2
34+
return "Assume top earth layer"
35+
else
36+
return "Assume layer $(f.layer)"
37+
end
38+
end
39+
40+
"""
41+
$(TYPEDSIGNATURES)
42+
43+
Functor implementation for `EnforceLayer`.
44+
45+
Builds a 2-layer (air + one enforced earth layer) data pack as three matrices
46+
ρ, ε, μ of size (2 × n_freq), already converted to `T`.
47+
48+
# Returns
49+
- `(ρ, ε, μ) :: (Matrix{T}, Matrix{T}, Matrix{T})`
50+
with row 1 = air, row 2 = enforced earth layer.
51+
"""
52+
function (f::EnforceLayer)(
53+
model::EarthModel,
54+
freq::AbstractVector{<:REALSCALAR},
55+
::Type{T},
56+
) where {T <: REALSCALAR}
57+
58+
nL = length(model.layers)
59+
nF = length(freq)
60+
61+
layer_idx = f.layer == -1 ? nL : f.layer
62+
(2 <= layer_idx <= nL) || error(
63+
"Invalid layer index: $layer_idx. Model has $nL layers (including air). " *
64+
"Valid earth layer indices are 2:$nL.",
65+
)
66+
67+
Lair = model.layers[1]
68+
Lsel = model.layers[layer_idx]
69+
70+
ρ = Matrix{T}(undef, 2, nF)
71+
ε = similar(ρ)
72+
μ = similar(ρ)
73+
74+
@inbounds for j in 1:nF
75+
ρ[1, j] = T(Lair.rho_g[j])
76+
ε[1, j] = T(Lair.eps_g[j])
77+
μ[1, j] = T(Lair.mu_g[j])
78+
79+
ρ[2, j] = T(Lsel.rho_g[j])
80+
ε[2, j] = T(Lsel.eps_g[j])
81+
μ[2, j] = T(Lsel.mu_g[j])
82+
end
83+
84+
return ρ, ε, μ
85+
end
86+
87+
# """
88+
# $(TYPEDSIGNATURES)
89+
90+
# Functor implementation for `EnforceLayer`.
91+
92+
# Takes a multi-layer `EarthModel` and returns a new two-layer model (air + one effective earth layer) based on the properties of the layer specified in the `EnforceLayer` instance.
93+
94+
# # Returns
95+
# - A `Vector{EarthLayer}` containing two layers: the original air layer and the selected earth layer.
96+
# """
97+
# function (f::EnforceLayer)(
98+
# model::EarthModel,
99+
# freq::Vector{<:REALSCALAR},
100+
# T::DataType,
101+
# )
102+
# num_layers = length(model.layers)
103+
104+
# # Determine the index of the layer to select
105+
# layer_idx = f.layer == -1 ? num_layers : f.layer
106+
107+
# # Validate the chosen index
108+
# if !(2 <= layer_idx <= num_layers)
109+
# Base.error(
110+
# "Invalid layer index: $layer_idx. The model only has $num_layers layers (including air). Valid earth layer indices are from 2 to $num_layers.",
111+
# )
112+
# end
113+
114+
# # The air layer is always the first layer in the original model
115+
# air_layer = model.layers[1]
116+
117+
# # The enforced earth layer is the one at the selected index
118+
# enforced_layer = model.layers[layer_idx]
119+
120+
# # Create a NamedTuple for the air layer with type-promoted property vectors
121+
# air_data = (
122+
# rho_g = T.(air_layer.rho_g),
123+
# eps_g = T.(air_layer.eps_g),
124+
# mu_g = T.(air_layer.mu_g),
125+
# )
126+
127+
# # Create a NamedTuple for the enforced earth layer
128+
# earth_data = (
129+
# rho_g = T.(enforced_layer.rho_g),
130+
# eps_g = T.(enforced_layer.eps_g),
131+
# mu_g = T.(enforced_layer.mu_g),
132+
# )
133+
134+
# # Return a new vector containing only these two layers
135+
# return [air_data, earth_data]
136+
# end
Lines changed: 18 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,18 @@
1+
struct Lossless <: InsulationAdmittanceFormulation end
2+
get_description(::Lossless) = "Lossless insulation (ideal dielectric)"
3+
4+
@inline function (f::Lossless)(r_in::T, r_ex::T, epsr_i::T, freq::T
5+
) where {T <: REALSCALAR}
6+
7+
if r_ex == r_in
8+
# TODO: Implement consistent handling of admittance for bare conductors
9+
return zero(Complex{T})
10+
end
11+
12+
# Constants
13+
# ω = 2π * freq
14+
eps_i = T(ε₀) * epsr_i
15+
16+
17+
return Complex{T}(log(r_ex / r_in) / (2π * eps_i))
18+
end
Lines changed: 18 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,18 @@
1+
2+
struct Lossless <: InsulationImpedanceFormulation end
3+
get_description(::Lossless) = "Lossless insulation (ideal dielectric)"
4+
5+
@inline function (f::Lossless)(r_in::T, r_ex::T, mur_i::T, freq::T
6+
) where {T <: REALSCALAR}
7+
8+
if r_ex == r_in
9+
# TODO: Implement consistent handling of admittance for bare conductors
10+
return zero(Complex{T})
11+
end
12+
13+
# Constants
14+
mu_i = T(μ₀) * mur_i
15+
ω = 2π * freq
16+
17+
return Complex{T}(im * ω * mu_i * log(r_ex / r_in) / 2π)
18+
end

0 commit comments

Comments
 (0)