Mercurial > repos > public > sbplib_julia
changeset 1562:efa994405c38 feature/sbp_operators/laplace_curvilinear
Add attempt at laplace for mapped grids
| author | Jonatan Werpers <jonatan@werpers.com> |
|---|---|
| date | Thu, 25 Apr 2024 09:34:34 +0200 |
| parents | 6fdc81860b0c |
| children | 6e910408c51a |
| files | src/SbpOperators/volumeops/laplace/laplace.jl |
| diffstat | 1 files changed, 22 insertions(+), 0 deletions(-) [+] |
line wrap: on
line diff
--- a/src/SbpOperators/volumeops/laplace/laplace.jl Thu Apr 25 09:03:35 2024 +0200 +++ b/src/SbpOperators/volumeops/laplace/laplace.jl Thu Apr 25 09:34:34 2024 +0200 @@ -52,3 +52,25 @@ return Δ end laplace(g::EquidistantGrid, stencil_set) = second_derivative(g, stencil_set) + + +function laplace(g::MappedGrid, stencil_set) + J = jacobian_determinant(g) + J⁻¹ = map(inv, J) + + Jḡ = map(*, J, ggeometric_tensor_inverse(g)) + lg = logicalgrid(g) + + return mapreduce(+, CartesianIndices(first(ḡ))) do I + i,j = I[1], I[2] + Jgⁱʲ = componentview(Jḡ, I[1], I[2]) + + if i == j + J⁻¹∘second_derivative_variable(lg, Jgⁱʲ, stencil_set, i) + else + Dᵢ = first_derivative(lg, stencil_set, i) + Dⱼ = first_derivative(lg, stencil_set, j) + J⁻¹∘Dᵢ∘DiagonalTensor(Jgⁱʲ)∘Dⱼ + end + end +end