Skip to content

Commit

Permalink
Merge pull request #106 from paresula/feature/LinearBasisWilton_new
Browse files Browse the repository at this point in the history 
Feature Linear Basis function
  • Loading branch information
@krcools
krcools authored Sep 21, 2023
2 parents 0d0b333 + 84aefd9 commit 6f65cc8
Show file tree
Hide file tree
Showing 7 changed files with 1,320 additions and 114 deletions.
2 changes: 1 addition & 1 deletion Project.toml
View file
Original file line number Diff line number Diff line change
Expand Up @@ -51,7 +51,7 @@ SauterSchwabQuadrature = "2.2.0"
SparseMatrixDicts = "0.2"
SpecialFunctions = "0.7, 0.8, 0.9, 0.10, 1, 2"
StaticArrays = "0.8.3, 0.9, 0.10, 0.11, 0.12, 1"
WiltonInts84 = "0.2.3"
WiltonInts84 = "0.2.5"
julia = "1.6"


3 changes: 2 additions & 1 deletion src/BEAST.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -215,7 +215,7 @@ include("maxwell/spotential.jl")
include("maxwell/maxwell.jl")
include("maxwell/sourcefield.jl")

# Support for the Helmholtz equatio
# Support for the Helmholtz equation
include("helmholtz2d/helmholtzop.jl")

include("helmholtz3d/hh3dexc.jl")
Expand All @@ -225,6 +225,7 @@ include("helmholtz3d/hh3dnear.jl")
include("helmholtz3d/hh3dfar.jl")
include("helmholtz3d/hh3d_sauterschwabqr.jl")
include("helmholtz3d/helmholtz3d.jl")
include("helmholtz3d/wiltonints.jl")

#suport for Volume Integral equation
include("volumeintegral/vie.jl")
Expand Down
84 changes: 44 additions & 40 deletions src/helmholtz3d/hh3d_sauterschwabqr.jl
View file
Original file line number Diff line number Diff line change
@@ -1,6 +1,6 @@
function pulled_back_integrand(op::HH3DSingleLayerFDBIO,
test_local_space::LagrangeRefSpace{<:Any,0},
trial_local_space::LagrangeRefSpace{<:Any,0},
test_local_space::LagrangeRefSpace,
trial_local_space::LagrangeRefSpace,
test_chart, trial_chart)

(u,v) -> begin
Expand All @@ -20,14 +20,14 @@ function pulled_back_integrand(op::HH3DSingleLayerFDBIO,

αjG = α*G*j

SA[f[1].value * αjG * g[1].value]
SMatrix{length(f),length(g)}((f[j].value * αjG * g[i].value for i in 1:length(g) for j in 1:length(f) )...)
end
end


function pulled_back_integrand(op::HH3DHyperSingularFDBIO,
test_local_space::LagrangeRefSpace{<:Any,1},
trial_local_space::LagrangeRefSpace{<:Any,1},
test_local_space::LagrangeRefSpace,
trial_local_space::LagrangeRefSpace,
test_chart, trial_chart)

(u,v) -> begin
Expand All @@ -52,25 +52,16 @@ function pulled_back_integrand(op::HH3DHyperSingularFDBIO,
αjG = ny*α*G*j
βjG = β*G*j

A = SA[ αjG*g[1].value, αjG*g[2].value, αjG*g[3].value]
B = SA[ βjG*g[1].curl, βjG*g[2].curl, βjG*g[3].curl ]

SMatrix{3,3}((
dot(nx*f[1].value, A[1]) + dot(f[1].curl, B[1]),
dot(nx*f[2].value, A[1]) + dot(f[2].curl, B[1]),
dot(nx*f[3].value, A[1]) + dot(f[3].curl, B[1]),
dot(nx*f[1].value, A[2]) + dot(f[1].curl, B[2]),
dot(nx*f[2].value, A[2]) + dot(f[2].curl, B[2]),
dot(nx*f[3].value, A[2]) + dot(f[3].curl, B[2]),
dot(nx*f[1].value, A[3]) + dot(f[1].curl, B[3]),
dot(nx*f[2].value, A[3]) + dot(f[2].curl, B[3]),
dot(nx*f[3].value, A[3]) + dot(f[3].curl, B[3])))
A = SA[(αjG*g[i].value for i in 1:length(g))...]
B = SA[(βjG*g[i].curl for i in 1:length(g))...]

SMatrix{length(f),length(g)}((((dot(nx*f[j].value,A[i])+dot(f[j].curl,B[i])) for i in 1:length(g) for j in 1:length(f))...))
end
end

function pulled_back_integrand(op::HH3DDoubleLayerFDBIO,
test_local_space::LagrangeRefSpace{<:Any,0},
trial_local_space::LagrangeRefSpace{<:Any,1},
test_local_space::LagrangeRefSpace,
trial_local_space::LagrangeRefSpace,
test_chart, trial_chart)

(u,v) -> begin
Expand All @@ -94,16 +85,13 @@ function pulled_back_integrand(op::HH3DDoubleLayerFDBIO,
∇G = -(γ + inv_R) * G * inv_R * r
αnyj∇G = dot(ny,-α*∇G*j)

SMatrix{1,3}(
f[1].value * αnyj∇G * g[1].value,
f[1].value * αnyj∇G * g[2].value,
f[1].value * αnyj∇G * g[3].value)
SMatrix{length(f),length(g)}((f[j].value * αnyj∇G * g[i].value for i in 1:length(g) for j in 1:length(f))...)
end
end

function pulled_back_integrand(op::HH3DDoubleLayerTransposedFDBIO,
test_local_space::LagrangeRefSpace{<:Any,1},
trial_local_space::LagrangeRefSpace{<:Any,0},
test_local_space::LagrangeRefSpace,
trial_local_space::LagrangeRefSpace,
test_chart, trial_chart)

(u,v) -> begin
Expand All @@ -127,14 +115,11 @@ function pulled_back_integrand(op::HH3DDoubleLayerTransposedFDBIO,
∇G = -(γ + inv_R) * G * inv_R * r
αnxj∇G = dot(nx,α*∇G*j)

SMatrix{3,1}(
f[1].value * αnxj∇G * g[1].value,
f[2].value * αnxj∇G * g[1].value,
f[3].value * αnxj∇G * g[1].value)
SMatrix{length(f),length(g)}((f[j].value * αnxj∇G * g[i].value for i in 1:length(g) for j in 1:length(f))...)
end
end

function momintegrals!(op::HH3DSingleLayerFDBIO,
function momintegrals!(op::Helmholtz3DOp,
test_local_space::LagrangeRefSpace{<:Any,0}, trial_local_space::LagrangeRefSpace{<:Any,0},
test_triangular_element, trial_triangular_element, out, strat::SauterSchwabStrategy)

Expand All @@ -152,19 +137,21 @@ function momintegrals!(op::HH3DSingleLayerFDBIO,
trial_triangular_element.vertices[J[2]],
trial_triangular_element.vertices[J[3]])

# igd = MWSL3DIntegrand(test_triangular_element, trial_triangular_element,
# op, test_local_space, trial_local_space)
test_sign = Combinatorics.levicivita(I)
trial_sign = Combinatorics.levicivita(J)
σ = momintegrals_sign(op, test_sign, trial_sign)

igd = pulled_back_integrand(op, test_local_space, trial_local_space,
test_triangular_element, trial_triangular_element)
G = SauterSchwabQuadrature.sauterschwab_parameterized(igd, strat)
out[1,1] += G[1,1]
out[1,1] += G[1,1] * σ

nothing
end


function momintegrals!(op::HH3DHyperSingularFDBIO,
test_local_space::LagrangeRefSpace{<:Any,1}, trial_local_space::LagrangeRefSpace{<:Any,1},
function momintegrals!(op::Helmholtz3DOp,
test_local_space::LagrangeRefSpace, trial_local_space::LagrangeRefSpace,
test_triangular_element, trial_triangular_element, out, strat::SauterSchwabStrategy)

I, J, K, L = SauterSchwabQuadrature.reorder(
Expand All @@ -183,13 +170,14 @@ function momintegrals!(op::HH3DHyperSingularFDBIO,

test_sign = Combinatorics.levicivita(I)
trial_sign = Combinatorics.levicivita(J)
σ = test_sign * trial_sign

σ = momintegrals_sign(op, test_sign, trial_sign)

igd = pulled_back_integrand(op, test_local_space, trial_local_space,
test_triangular_element, trial_triangular_element)
G = SauterSchwabQuadrature.sauterschwab_parameterized(igd, strat)
for j ∈ 1:3, i ∈ 1:3
out[i,j] += G[K[i],L[j]] * σ
out[i,j] += G[K[i],L[j]] * σ
end

nothing
Expand All @@ -213,9 +201,10 @@ function momintegrals!(op::Helmholtz3DOp,
trial_triangular_element.vertices[J[2]],
trial_triangular_element.vertices[J[3]])

test_sign = Combinatorics.levicivita(I)
trial_sign = Combinatorics.levicivita(J)

σ = trial_sign
σ = momintegrals_sign(op, test_sign, trial_sign)

igd = pulled_back_integrand(op, test_local_space, trial_local_space,
test_triangular_element, trial_triangular_element)
Expand Down Expand Up @@ -247,7 +236,9 @@ function momintegrals!(op::Helmholtz3DOp,
trial_triangular_element.vertices[J[3]])

test_sign = Combinatorics.levicivita(I)
σ = test_sign
trial_sign = Combinatorics.levicivita(J)

σ = momintegrals_sign(op, test_sign, trial_sign)

igd = pulled_back_integrand(op, test_local_space, trial_local_space,
test_triangular_element, trial_triangular_element)
Expand All @@ -259,3 +250,16 @@ function momintegrals!(op::Helmholtz3DOp,

nothing
end

function momintegrals_sign(op::HH3DSingleLayerFDBIO, test_sign, trial_sign)
return 1
end
function momintegrals_sign(op::HH3DDoubleLayerFDBIO, test_sign, trial_sign)
return trial_sign
end
function momintegrals_sign(op::HH3DDoubleLayerTransposedFDBIO, test_sign, trial_sign)
return test_sign
end
function momintegrals_sign(op::HH3DHyperSingularFDBIO, test_sign, trial_sign)
return test_sign * trial_sign
end
Loading

0 comments on commit 6f65cc8

Please sign in to comment.