Mercurial > repos > public > sbplib_julia

comparison src/SbpOperators/volumeops/laplace/laplace.jl @ 677:011863b3f24c feature/laplace_opset

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Make use of the function boundary_quadrature in Laplace constructor
author Vidar Stiernström <vidar.stiernstrom@it.uu.se>
date Sat, 06 Feb 2021 15:17:18 +0100
parents 2d56a53a1646
children 3cd582257072
comparison
equal deleted inserted replaced
676:bf48761c1345 677:011863b3f24c
3 Laplace(grid::EquidistantGrid, fn; order) 3 Laplace(grid::EquidistantGrid, fn; order)
4 4
5 Implements the Laplace operator, approximating ∑d²/xᵢ² , i = 1,...,`Dim` as a 5 Implements the Laplace operator, approximating ∑d²/xᵢ² , i = 1,...,`Dim` as a
6 `TensorMapping`. Additionally, `Laplace` stores the quadrature, and boundary 6 `TensorMapping`. Additionally, `Laplace` stores the quadrature, and boundary
7 operators relevant for constructing a SBP finite difference scheme as `TensorMapping`s. 7 operators relevant for constructing a SBP finite difference scheme as `TensorMapping`s.
8
9 Laplace(grid::EquidistantGrid, fn; order) creates the Laplace operator on an
10 equidistant grid, where the operators are read from a TOML. The laplace operator
11 is created using laplace(grid,...).
8 """ 12 """
9 struct Laplace{T, Dim, Rb, TMdiffop<:TensorMapping{T,Dim,Dim}, # Differential operator tensor mapping 13 struct Laplace{T, Dim, Rb, TMdiffop<:TensorMapping{T,Dim,Dim}, # Differential operator tensor mapping
10 TMqop<:TensorMapping{T,Dim,Dim}, # Quadrature operator tensor mapping 14 TMqop<:TensorMapping{T,Dim,Dim}, # Quadrature operator tensor mapping
11 TMbop<:TensorMapping{T,Rb,Dim}, # Boundary operator tensor mapping 15 TMbop<:TensorMapping{T,Rb,Dim}, # Boundary operator tensor mapping
12 TMbqop<:TensorMapping{T,Rb,Rb}, # Boundary quadrature tensor mapping 16 TMbqop<:TensorMapping{T,Rb,Rb}, # Boundary quadrature tensor mapping
31 e_closure_stencil = op.eClosure 35 e_closure_stencil = op.eClosure
32 d_closure_stencil = op.dClosure 36 d_closure_stencil = op.dClosure
33 37
34 # Volume operators 38 # Volume operators
35 Δ = laplace(grid, D_inner_stecil, D_closure_stencils) 39 Δ = laplace(grid, D_inner_stecil, D_closure_stencils)
36 H = DiagonalQuadrature(grid, H_closure_stencils) 40 H = quadrature(grid, H_closure_stencils)
37 H⁻¹ = InverseDiagonalQuadrature(grid, H_closure_stencils) 41 H⁻¹ = InverseDiagonalQuadrature(grid, H_closure_stencils)
38 42
39 # Pair operators with boundary ids 43 # Boundary operator - id pairs
40 bids = boundary_identifiers(grid) 44 bids = boundary_identifiers(grid)
41 # Boundary operators
42 e_pairs = ntuple(i -> Pair(bids[i],BoundaryRestriction(grid,e_closure_stencil,bids[i])),length(bids)) 45 e_pairs = ntuple(i -> Pair(bids[i],BoundaryRestriction(grid,e_closure_stencil,bids[i])),length(bids))
43 d_pairs = ntuple(i -> Pair(bids[i],NormalDerivative(grid,d_closure_stencil,bids[i])),length(bids)) 46 d_pairs = ntuple(i -> Pair(bids[i],NormalDerivative(grid,d_closure_stencil,bids[i])),length(bids))
44 # Boundary quadratures are constructed on the lower-dimensional grid defined 47 Hᵧ_pairs = ntuple(i -> Pair(bids[i],boundary_quadrature(grid,H_closure_stencils,bids[i])),length(bids))
45 # by the coordinite directions orthogonal to that of the boundary.
46 dims = collect(1:dimension(grid))
47 orth_grids = ntuple(i -> restrict(grid,dims[dims .!= dim(bids[i])]),length(bids))
48 Hᵧ_pairs = ntuple(i -> Pair(bids[i],DiagonalQuadrature(orth_grids[i],H_closure_stencils)),length(bids))
49 48
50 return Laplace(Δ, H, H⁻¹, Dict(e_pairs), Dict(d_pairs), Dict(Hᵧ_pairs)) 49 return Laplace(Δ, H, H⁻¹, Dict(e_pairs), Dict(d_pairs), Dict(Hᵧ_pairs))
51 end 50 end
52 51
53 LazyTensors.range_size(L::Laplace) = LazyTensors.range_size(L.D) 52 LazyTensors.range_size(L::Laplace) = LazyTensors.range_size(L.D)