Mercurial > repos > public > sbplib_julia

comparison DiffOps/src/laplace.jl @ 268:f67ce2eb6019 boundary_conditions

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Remove property spacing in EquidistantGrid, due to redundancy.
author Vidar Stiernström <vidar.stiernstrom@it.uu.se>
date Thu, 05 Dec 2019 09:44:34 +0100
parents 9ad447176ba1
children ccef055233a2
comparison
equal deleted inserted replaced
267:634453a4e1d8 268:f67ce2eb6019
8 error("not implemented") 8 error("not implemented")
9 end 9 end
10 10
11 # u = L*v 11 # u = L*v
12 function apply(L::Laplace{1}, v::AbstractVector, i::Int) 12 function apply(L::Laplace{1}, v::AbstractVector, i::Int)
13 uᵢ = L.a * SbpOperators.apply(L.op, L.grid.spacing[1], v, i) 13 uᵢ = L.a * SbpOperators.apply(L.op, inverse_spacing(L.grid)[1], v, i)
14 return uᵢ 14 return uᵢ
15 end 15 end
16 16
17 @inline function apply(L::Laplace{2}, v::AbstractArray{T,2} where T, I::Tuple{Index{R1}, Index{R2}}) where {R1, R2} 17 @inline function apply(L::Laplace{2}, v::AbstractArray{T,2} where T, I::Tuple{Index{R1}, Index{R2}}) where {R1, R2}
18 # 2nd x-derivative 18 # 2nd x-derivative
19 @inbounds vx = view(v, :, Int(I[2])) 19 @inbounds vx = view(v, :, Int(I[2]))
20 @inbounds uᵢ = L.a*SbpOperators.apply(L.op, L.grid.inverse_spacing[1], vx , I[1]) 20 @inbounds uᵢ = L.a*SbpOperators.apply(L.op, inverse_spacing(L.grid)[1], vx , I[1])
21 # 2nd y-derivative 21 # 2nd y-derivative
22 @inbounds vy = view(v, Int(I[1]), :) 22 @inbounds vy = view(v, Int(I[1]), :)
23 @inbounds uᵢ += L.a*SbpOperators.apply(L.op, L.grid.inverse_spacing[2], vy, I[2]) 23 @inbounds uᵢ += L.a*SbpOperators.apply(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 24 # NOTE: the package qualifier 'SbpOperators' can problably be removed once all "applying" objects use LazyTensors
25 return uᵢ 25 return uᵢ
26 end 26 end
27 27
28 # Slow but maybe convenient? 28 # Slow but maybe convenient?
54 # TODO: Dispatch on Tuple{Index{R1},Index{R2}}? 54 # TODO: Dispatch on Tuple{Index{R1},Index{R2}}?
55 @inline function LazyTensors.apply(H::Quadrature{2}, v::AbstractArray{T,2} where T, I::NTuple{2,Integer}) 55 @inline function LazyTensors.apply(H::Quadrature{2}, v::AbstractArray{T,2} where T, I::NTuple{2,Integer})
56 I = CartesianIndex(I); 56 I = CartesianIndex(I);
57 N = size(H.grid) 57 N = size(H.grid)
58 # Quadrature in x direction 58 # Quadrature in x direction
59 @inbounds q = apply_quadrature(H.op, H.grid.spacing[1], v[I] , I[1], N[1]) 59 @inbounds q = apply_quadrature(H.op, spacing(H.grid)[1], v[I] , I[1], N[1])
60 # Quadrature in y-direction 60 # Quadrature in y-direction
61 @inbounds q = apply_quadrature(H.op, H.grid.spacing[2], q, I[2], N[2]) 61 @inbounds q = apply_quadrature(H.op, spacing(H.grid)[2], q, I[2], N[2])
62 return q 62 return q
63 end 63 end
64 64
65 LazyTensors.apply_transpose(H::Quadrature{2}, v::AbstractArray{T,2} where T, I::NTuple{2,Integer}) = LazyTensors.apply(H,v,I) 65 LazyTensors.apply_transpose(H::Quadrature{2}, v::AbstractArray{T,2} where T, I::NTuple{2,Integer}) = LazyTensors.apply(H,v,I)
66 66
81 # TODO: Dispatch on Tuple{Index{R1},Index{R2}}? 81 # TODO: Dispatch on Tuple{Index{R1},Index{R2}}?
82 @inline function LazyTensors.apply(H_inv::InverseQuadrature{2}, v::AbstractArray{T,2} where T, I::NTuple{2,Integer}) 82 @inline function LazyTensors.apply(H_inv::InverseQuadrature{2}, v::AbstractArray{T,2} where T, I::NTuple{2,Integer})
83 I = CartesianIndex(I); 83 I = CartesianIndex(I);
84 N = size(H_inv.grid) 84 N = size(H_inv.grid)
85 # Inverse quadrature in x direction 85 # Inverse quadrature in x direction
86 @inbounds q_inv = apply_inverse_quadrature(H_inv.op, H_inv.grid.inverse_spacing[1], v[I] , I[1], N[1]) 86 @inbounds q_inv = apply_inverse_quadrature(H_inv.op, inverse_spacing(H_inv.grid)[1], v[I] , I[1], N[1])
87 # Inverse quadrature in y-direction 87 # Inverse quadrature in y-direction
88 @inbounds q_inv = apply_inverse_quadrature(H_inv.op, H_inv.grid.inverse_spacing[2], q_inv, I[2], N[2]) 88 @inbounds q_inv = apply_inverse_quadrature(H_inv.op, inverse_spacing(H_inv.grid)[2], q_inv, I[2], N[2])
89 return q_inv 89 return q_inv
90 end 90 end
91 91
92 LazyTensors.apply_transpose(H_inv::InverseQuadrature{2}, v::AbstractArray{T,2} where T, I::NTuple{2,Integer}) = LazyTensors.apply(H_inv,v,I) 92 LazyTensors.apply_transpose(H_inv::InverseQuadrature{2}, v::AbstractArray{T,2} where T, I::NTuple{2,Integer}) = LazyTensors.apply(H_inv,v,I)
93 93
142 # TODO: Make this independent of dimension 142 # TODO: Make this independent of dimension
143 function LazyTensors.apply(d::NormalDerivative, v::AbstractArray, I::NTuple{2,Int}) 143 function LazyTensors.apply(d::NormalDerivative, v::AbstractArray, I::NTuple{2,Int})
144 i = I[dim(d.bId)] 144 i = I[dim(d.bId)]
145 j = I[3-dim(d.bId)] 145 j = I[3-dim(d.bId)]
146 N_i = size(d.grid)[dim(d.bId)] 146 N_i = size(d.grid)[dim(d.bId)]
147 h_inv = d.grid.inverse_spacing[dim(d.bId)] 147 h_inv = inverse_spacing(d.grid)[dim(d.bId)]
148 return apply_d(d.op, h_inv, v[j], N_i, i, region(d.bId)) 148 return apply_d(d.op, h_inv, v[j], N_i, i, region(d.bId))
149 end 149 end
150 150
151 function LazyTensors.apply_transpose(d::NormalDerivative, v::AbstractArray, I::NTuple{1,Int}) 151 function LazyTensors.apply_transpose(d::NormalDerivative, v::AbstractArray, I::NTuple{1,Int})
152 u = selectdim(v,3-dim(d.bId),I[1]) 152 u = selectdim(v,3-dim(d.bId),I[1])
153 return apply_d_T(d.op, d.grid.inverse_spacing[dim(d.bId)], u, region(d.bId)) 153 return apply_d_T(d.op, inverse_spacing(d.grid)[dim(d.bId)], u, region(d.bId))
154 end 154 end
155 155
156 """ 156 """
157 BoundaryQuadrature{T,N,M,K} <: TensorOperator{T,1} 157 BoundaryQuadrature{T,N,M,K} <: TensorOperator{T,1}
158 158
166 export BoundaryQuadrature 166 export BoundaryQuadrature
167 167
168 # TODO: Make this independent of dimension 168 # TODO: Make this independent of dimension
169 # TODO: Dispatch directly on Index{R}? 169 # TODO: Dispatch directly on Index{R}?
170 function LazyTensors.apply(q::BoundaryQuadrature{T}, v::AbstractArray{T,1}, I::NTuple{1,Int}) where T 170 function LazyTensors.apply(q::BoundaryQuadrature{T}, v::AbstractArray{T,1}, I::NTuple{1,Int}) where T
171 h = q.grid.spacing[3-dim(q.bId)] 171 h = spacing(q.grid)[3-dim(q.bId)]
172 N = size(v) 172 N = size(v)
173 return apply_quadrature(q.op, h, v[I[1]], I[1], N[1]) 173 return apply_quadrature(q.op, h, v[I[1]], I[1], N[1])
174 end 174 end
175 175
176 LazyTensors.apply_transpose(q::BoundaryQuadrature{T}, v::AbstractArray{T,1}, I::NTuple{1,Int}) where T = LazyTensors.apply(q,v,I) 176 LazyTensors.apply_transpose(q::BoundaryQuadrature{T}, v::AbstractArray{T,1}, I::NTuple{1,Int}) where T = LazyTensors.apply(q,v,I)