Mercurial > repos > public > sbplib_julia

comparison SbpOperators/src/laplace/secondderivative.jl @ 328:9cc5d1498b2d

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Refactor 1D diagonal inner product in quadrature.jl to separate file. Write tests for quadratures. Clean up laplace and secondderivative
author Vidar Stiernström <vidar.stiernstrom@it.uu.se>
date Thu, 24 Sep 2020 22:31:48 +0200
parents f2d6ec89dfc5
children
comparison
equal deleted inserted replaced
327:802edc9f252e 328:9cc5d1498b2d
18 #TODO: The 1D tensor mappings should not have to dispatch on 1D tuples if we write LazyTensor.apply for vararg right?!?! 18 #TODO: The 1D tensor mappings should not have to dispatch on 1D tuples if we write LazyTensor.apply for vararg right?!?!
19 # Currently have to index the Tuple{Index} in each method in order to call the stencil methods which is ugly. 19 # Currently have to index the Tuple{Index} in each method in order to call the stencil methods which is ugly.
20 # I thought I::Vararg{Index,R} fell back to just Index for R = 1 20 # I thought I::Vararg{Index,R} fell back to just Index for R = 1
21 21
22 # Apply for different regions Lower/Interior/Upper or Unknown region 22 # Apply for different regions Lower/Interior/Upper or Unknown region
23 @inline function LazyTensors.apply(D2::SecondDerivative{T}, v::AbstractVector{T}, I::Index{Lower}) where T 23 function LazyTensors.apply(D2::SecondDerivative{T}, v::AbstractVector{T}, I::Index{Lower}) where T
24 return @inbounds D2.h_inv*D2.h_inv*apply_stencil(D2.closureStencils[Int(I)], v, Int(I)) 24 return @inbounds D2.h_inv*D2.h_inv*apply_stencil(D2.closureStencils[Int(I)], v, Int(I))
25 end 25 end
26 26
27 @inline function LazyTensors.apply(D2::SecondDerivative{T}, v::AbstractVector{T}, I::Index{Interior}) where T 27 function LazyTensors.apply(D2::SecondDerivative{T}, v::AbstractVector{T}, I::Index{Interior}) where T
28 return @inbounds D2.h_inv*D2.h_inv*apply_stencil(D2.innerStencil, v, Int(I)) 28 return @inbounds D2.h_inv*D2.h_inv*apply_stencil(D2.innerStencil, v, Int(I))
29 end 29 end
30 30
31 @inline function LazyTensors.apply(D2::SecondDerivative{T}, v::AbstractVector{T}, I::Index{Upper}) where T 31 function LazyTensors.apply(D2::SecondDerivative{T}, v::AbstractVector{T}, I::Index{Upper}) where T
32 N = length(v) # TODO: Use domain_size here instead? N = domain_size(D2,size(v)) 32 N = length(v) # TODO: Use domain_size here instead? N = domain_size(D2,size(v))
33 return @inbounds D2.h_inv*D2.h_inv*Int(D2.parity)*apply_stencil_backwards(D2.closureStencils[N-Int(I)+1], v, Int(I)) 33 return @inbounds D2.h_inv*D2.h_inv*Int(D2.parity)*apply_stencil_backwards(D2.closureStencils[N-Int(I)+1], v, Int(I))
34 end 34 end
35 35
36 @inline function LazyTensors.apply(D2::SecondDerivative{T}, v::AbstractVector{T}, index::Index{Unknown}) where T 36 function LazyTensors.apply(D2::SecondDerivative{T}, v::AbstractVector{T}, index::Index{Unknown}) where T
37 N = length(v) # TODO: Use domain_size here instead? 37 N = length(v) # TODO: Use domain_size here instead?
38 r = getregion(Int(index), closuresize(D2), N) 38 r = getregion(Int(index), closuresize(D2), N)
39 I = Index(Int(index), r) 39 I = Index(Int(index), r)
40 return LazyTensors.apply(D2, v, I) 40 return LazyTensors.apply(D2, v, I)
41 end 41 end
42 42
43 LazyTensors.apply_transpose(D2::SecondDerivative{T}, v::AbstractVector{T}, I::Index) where {T} = LazyTensors.apply(D2, v, I)
43 44
44 @inline function LazyTensors.apply_transpose(D2::SecondDerivative, v::AbstractVector, I::Index) 45 closuresize(D2::SecondDerivative{T,N,M,K}) where {T<:Real,N,M,K} = M
45 return LazyTensors.apply(D2, v, I)
46 end
47
48
49 function closuresize(D2::SecondDerivative{T,N,M,K}) where {T<:Real,N,M,K}
50 return M
51 end