Mercurial > repos > public > sbplib_julia

comparison DiffOps/src/laplace.jl @ 282:ce6a2f3f732a boundary_conditions

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Make Laplace a TensorOperator and add tests. NOTE: Two of the tests for Laplace2D are currently failing.
author Vidar Stiernström <vidar.stiernstrom@it.uu.se>
date Thu, 09 Jan 2020 10:54:24 +0100
parents 1eefaefdd0c7
children 12a12a5cd973
comparison
equal deleted inserted replaced
281:1eefaefdd0c7 282:ce6a2f3f732a
1 struct Laplace{Dim,T<:Real,N,M,K} <: DiffOpCartesian{Dim} 1 struct Laplace{Dim,T<:Real,N,M,K} <: TensorOperator{T,Dim}
2 grid::EquidistantGrid{Dim,T} 2 grid::EquidistantGrid{Dim,T}
3 a::T 3 a::T # TODO: Better name?
4 op::D2{Float64,N,M,K} 4 op::D2{T,N,M,K}
5 end 5 end
6 6 export Laplace
7 function apply(L::Laplace{Dim}, v::AbstractArray{T,Dim} where T, I::CartesianIndex{Dim}) where Dim 7
8 LazyTensors.domain_size(H::Laplace{Dim}, range_size::NTuple{Dim,Integer}) where Dim = size(L.grid)
9
10 function LazyTensors.apply(L::Laplace{Dim,T}, v::AbstractArray{T,Dim}, I::NTuple{Dim,Index}) where {T,Dim}
8 error("not implemented") 11 error("not implemented")
9 end 12 end
10 13
11 # u = L*v 14 # u = L*v
12 function apply(L::Laplace{1}, v::AbstractVector, i::Int) 15 function LazyTensors.apply(L::Laplace{1,T}, v::AbstractVector{T}, I::NTuple{1,Index}) where T
13 uᵢ = L.a * SbpOperators.apply_2nd_derivative(L.op, inverse_spacing(L.grid)[1], v, i) 16 uᵢ = L.a*apply_2nd_derivative(L.op, inverse_spacing(L.grid)[1], v, I[1])
14 return uᵢ 17 return uᵢ
15 end 18 end
16 19
17 @inline function apply(L::Laplace{2}, v::AbstractArray{T,2} where T, I::Tuple{Index{R1}, Index{R2}}) where {R1, R2} 20
21 @inline function LazyTensors.apply(L::Laplace{2,T}, v::AbstractArray{T,2}, I::NTuple{2,Index}) where T
18 # 2nd x-derivative 22 # 2nd x-derivative
19 @inbounds vx = view(v, :, Int(I[2])) 23 @inbounds vx = view(v, :, Int(I[2]))
20 @inbounds uᵢ = L.a*SbpOperators.apply_2nd_derivative(L.op, inverse_spacing(L.grid)[1], vx , I[1]) 24 @inbounds uᵢ = L.a*apply_2nd_derivative(L.op, inverse_spacing(L.grid)[1], vx , I[1])
21 # 2nd y-derivative 25 # 2nd y-derivative
22 @inbounds vy = view(v, Int(I[1]), :) 26 @inbounds vy = view(v, Int(I[1]), :)
23 @inbounds uᵢ += L.a*SbpOperators.apply_2nd_derivative(L.op, inverse_spacing(L.grid)[2], vy, I[2]) 27 @inbounds uᵢ += L.a*apply_2nd_derivative(L.op, inverse_spacing(L.grid)[2], vy, I[2])
24 # NOTE: the package qualifier 'SbpOperators' can problably be removed once all "applying" objects use LazyTensors
25 return uᵢ 28 return uᵢ
26 end
27
28 # Slow but maybe convenient?
29 function apply(L::Laplace{2}, v::AbstractArray{T,2} where T, i::CartesianIndex{2})
30 I = Index{Unknown}.(Tuple(i))
31 apply(L, v, I)
32 end 29 end
33 30
34 quadrature(L::Laplace) = Quadrature(L.op, L.grid) 31 quadrature(L::Laplace) = Quadrature(L.op, L.grid)
35 inverse_quadrature(L::Laplace) = InverseQuadrature(L.op, L.grid) 32 inverse_quadrature(L::Laplace) = InverseQuadrature(L.op, L.grid)
36 boundary_value(L::Laplace, bId::CartesianBoundary) = BoundaryValue(L.op, L.grid, bId) 33 boundary_value(L::Laplace, bId::CartesianBoundary) = BoundaryValue(L.op, L.grid, bId)