diff --git a/src/UniformField/excitation.jl b/src/UniformField/excitation.jl index 2370a2f..86ae92c 100644 --- a/src/UniformField/excitation.jl +++ b/src/UniformField/excitation.jl @@ -1,13 +1,13 @@ -struct UniformField{R,C,T} <: Excitation +struct UniformField{C,T,R} <: Excitation embedding::Medium{C} amplitude::T direction::SVector{3,R} # inner constructor: normalize direction - function UniformField(embedding::Medium{C}, amplitude::T, direction::SVector{3,R}) where {T,R,C} + function UniformField(embedding::Medium{C}, amplitude::T, direction::SVector{3,R}) where {C,R,T} dir_normalized = normalize(direction) - new{T,R,C}(embedding, amplitude, dir_normalized) + new{C,R,T}(embedding, amplitude, dir_normalized) end end diff --git a/src/UniformField/scattered.jl b/src/UniformField/scattered.jl index b716220..52a0f0c 100644 --- a/src/UniformField/scattered.jl +++ b/src/UniformField/scattered.jl @@ -122,8 +122,12 @@ Compute the electric field scattered by a layered dielectric sphere, for an inci The point and returned field are in Cartesian coordinates. """ function scatteredfield( - sphere::LayeredSphere, excitation::UniformField, point, quantity::ElectricField; parameter::Parameter=Parameter() -) + sphere::LayeredSphere{LN,LR,LC}, + excitation::UniformField{FC,FT,FR}, + point, + quantity::ElectricField; + parameter::Parameter=Parameter(), +) where {LN,LR,LC,FC,FT,FR} # using `Sihvola and Lindell, 1988, Transmission line analogy for calculating the effective permittivity of mixtures with spherical multilayer scatterers` E0 = excitation.amplitude @@ -134,7 +138,9 @@ function scatteredfield( n = length(a) perms = getfield.(vcat(sphere.embedding, sphere.filling), 1) - Bk = zeros(SMatrix{2,2,Float64}, n) + T = promote_type(LR, LC, FC, FT, FR) + + Bk = zeros(SMatrix{2,2,T}, n) B = Matrix(I, 2, 2) for k in range(n; stop=1, step=-1) @@ -147,8 +153,8 @@ function scatteredfield( B = Bk[k] * B end - Ck = zeros(Float64, n + 1) - Dk = zeros(Float64, n + 1) + Ck = zeros(T, n + 1) + Dk = zeros(T, n + 1) Ck[1] = 1 Dk[n + 1] = 0 @@ -183,8 +189,12 @@ Compute the scalar potential scattered by a layered dielectric sphere, for an in The point and returned field are in Cartesian coordinates. """ function scatteredfield( - sphere::LayeredSphere, excitation::UniformField, point, quantity::ScalarPotential; parameter::Parameter=Parameter() -) + sphere::LayeredSphere{LN,LR,LC}, + excitation::UniformField{FC,FT,FR}, + point, + quantity::ScalarPotential; + parameter::Parameter=Parameter(), +) where {LN,LR,LC,FC,FT,FR} # using `Sihvola and Lindell, 1988, Transmission line analogy for calculating the effective permittivity of mixtures with spherical multilayer scatterers` Φ0 = field(excitation, point, quantity) @@ -195,7 +205,9 @@ function scatteredfield( n = length(a) perms = getfield.(vcat(sphere.embedding, sphere.filling), 1) - Bk = zeros(SMatrix{2,2,Float64}, n) + T = promote_type(LR, LC, FC, FT, FR) + + Bk = zeros(SMatrix{2,2,T}, n) B = Matrix(I, 2, 2) for k in range(n; stop=1, step=-1) @@ -208,8 +220,8 @@ function scatteredfield( B = Bk[k] * B end - Ck = zeros(Float64, n + 1) - Dk = zeros(Float64, n + 1) + Ck = zeros(T, n + 1) + Dk = zeros(T, n + 1) Ck[1] = 1 Dk[n + 1] = 0 @@ -233,18 +245,23 @@ function scatteredfield( C = Ck[pos + 1] D = Dk[pos + 1] + return C * Φ0 - D * Φ0 / r^3 - Φ0 end """ - scatteredfield(sphere::LayeredSpherePEC, excitation::UniformField, point, quantity::ScalarPotential; parameter::Parameter=Parameter()) + scatteredfield(sphere::LayeredSpherePEC, excitation::UniformField{FC,FT,FR}, point, quantity::ScalarPotential; parameter::Parameter=Parameter()) Compute the scalar potential scattered by a layered dielectric sphere with PEC core, for an incident uniform field with polarization in x-direction. The point and returned field are in Cartesian coordinates. """ function scatteredfield( - sphere::LayeredSpherePEC, excitation::UniformField, point, quantity::ScalarPotential; parameter::Parameter=Parameter() -) + sphere::LayeredSpherePEC{LN,LD,LR,LC}, + excitation::UniformField{FC,FT,FR}, + point, + quantity::ScalarPotential; + parameter::Parameter=Parameter(), +) where {LN,LD,LR,LC,FC,FT,FR} # using `Sihvola and Lindell, 1988, Transmission line analogy for calculating the effective permittivity of mixtures with spherical multilayer scatterers` Φ0 = field(excitation, point, quantity) @@ -255,7 +272,9 @@ function scatteredfield( n = length(a) - 1 perms = getfield.(vcat(sphere.embedding, sphere.filling), 1) - Bk = zeros(SMatrix{2,2,Float64}, n) + T = promote_type(LR, LC, FC, FT, FR) + + Bk = zeros(SMatrix{2,2,T}, n) B = Matrix(I, 2, 2) for k in range(n; stop=1, step=-1) @@ -267,8 +286,8 @@ function scatteredfield( Bk[k] = 1 / (3 * perms[k]) * ([B11 B12; B21 B22]) B = Bk[k] * B end - Ck = zeros(Float64, n + 1) - Dk = zeros(Float64, n + 1) + Ck = zeros(T, n + 1) + Dk = zeros(T, n + 1) Ck[1] = 1 Dk[1] = (B[2, 1] + B[2, 2] * a[end]^3) / (B[1, 1] + B[1, 2] * a[end]^3) @@ -310,8 +329,12 @@ Compute the electric field scattered by a layered dielectric sphere with PEC cor The point and returned field are in Cartesian coordinates. """ function scatteredfield( - sphere::LayeredSpherePEC, excitation::UniformField, point, quantity::ElectricField; parameter::Parameter=Parameter() -) + sphere::LayeredSpherePEC{LN,LD,LR,LC}, + excitation::UniformField{FC,FT,FR}, + point, + quantity::ElectricField; + parameter::Parameter=Parameter(), +) where {LN,LD,LR,LC,FC,FT,FR} # using `Sihvola and Lindell, 1988, Transmission line analogy for calculating the effective permittivity of mixtures with spherical multilayer scatterers` E0 = excitation.amplitude @@ -321,7 +344,9 @@ function scatteredfield( n = length(a) - 1 perms = getfield.(vcat(sphere.embedding, sphere.filling), 1) - Bk = zeros(SMatrix{2,2,Float64}, n) + T = promote_type(LR, LC, FC, FT, FR) + + Bk = zeros(SMatrix{2,2,T}, n) B = Matrix(I, 2, 2) for k in range(n; stop=1, step=-1) @@ -334,8 +359,8 @@ function scatteredfield( B = Bk[k] * B end - Ck = zeros(Float64, n + 1) - Dk = zeros(Float64, n + 1) + Ck = zeros(T, n + 1) + Dk = zeros(T, n + 1) Ck[1] = 1 Dk[1] = (B[2, 1] + B[2, 2] * a[end]^3) / (B[1, 1] + B[1, 2] * a[end]^3) diff --git a/src/coordinateTransforms.jl b/src/coordinateTransforms.jl index 7dc89e0..4b2528c 100644 --- a/src/coordinateTransforms.jl +++ b/src/coordinateTransforms.jl @@ -104,11 +104,19 @@ function cart2sph(vec) y = vec[2] z = vec[3] - T = typeof(x) + T = eltype(x) r = sqrt(x^2 + y^2 + z^2) - t = (atan(sqrt(x^2 + y^2), z) + π) % π # map to range [0, π] - p = (atan(y, x) + 2 * π) % 2π # map to range [0, 2π] + if x ≈ T(0) && y ≈ T(0) + if z >= 0 + t = T(0) + else + t = T(π) + end + else + t = (atan(sqrt(x^2 + y^2), z) + π) % π # map to range [0, π] + end + p = (atan(y, x) + 2 * π) % 2π # map to range [0, 2π] return SVector{3,T}(r, t, p) end diff --git a/test/coordinateTransforms.jl b/test/coordinateTransforms.jl index 6ffce70..0e56c20 100644 --- a/test/coordinateTransforms.jl +++ b/test/coordinateTransforms.jl @@ -5,6 +5,11 @@ vec = SVector(1.0, π / 4, π / 4) xyz = SphericalScattering.sph2cart(vec) + ẑ = SVector(0.0, 0.0, 1.0) + rθϕ = SphericalScattering.cart2sph(-ẑ) + + @test rθϕ[2] ≈ π + @test xyz[1] ≈ 0.5 @test xyz[2] ≈ 0.5 @test xyz[3] ≈ 1 / √2