Mercurial > repos > public > sbplib_julia

comparison SbpOperators/src/d2.jl @ 269:ccef055233a2 boundary_conditions

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Move general methods for D2 to ConstantStencilOperator and increase verbosity of method names for clarity
author Vidar Stiernström <vidar.stiernstrom@it.uu.se>
date Thu, 05 Dec 2019 10:48:31 +0100
parents f89718833620
children 51e7de109c25
comparison
equal deleted inserted replaced
268:f67ce2eb6019 269:ccef055233a2
1 export D2, closuresize, readOperator, apply_e, apply_d, apply_e_T, apply_d_T 1 export D2, closuresize, readOperator
2 2
3 @enum Parity begin 3 @enum Parity begin
4 odd = -1 4 odd = -1
5 even = 1 5 even = 1
6 end 6 end
16 end 16 end
17 17
18 function closuresize(D::D2)::Int 18 function closuresize(D::D2)::Int
19 return length(D.quadratureClosure) 19 return length(D.quadratureClosure)
20 end 20 end
21
22 # TODO: Dispatch on Index{R}?
23 apply_quadrature(op::D2{T}, h::Real, v::T, i::Integer, N::Integer, ::Type{Lower}) where T = v*h*op.quadratureClosure[i]
24 apply_quadrature(op::D2{T}, h::Real, v::T, i::Integer, N::Integer, ::Type{Upper}) where T = v*h*op.quadratureClosure[N-i+1]
25 apply_quadrature(op::D2{T}, h::Real, v::T, i::Integer, N::Integer, ::Type{Interior}) where T = v*h
26
27 # TODO: Avoid branching in inner loops
28 function apply_quadrature(op::D2{T}, h::Real, v::T, i::Integer, N::Integer) where T
29 r = getregion(i, closuresize(op), N)
30 return apply_quadrature(op, h, v, i, N, r)
31 end
32 export apply_quadrature
33
34 # TODO: Dispatch on Index{R}?
35 apply_inverse_quadrature(op::D2{T}, h_inv::Real, v::T, i::Integer, N::Integer, ::Type{Lower}) where T = v*h_inv*op.inverseQuadratureClosure[i]
36 apply_inverse_quadrature(op::D2{T}, h_inv::Real, v::T, i::Integer, N::Integer, ::Type{Upper}) where T = v*h_inv*op.inverseQuadratureClosure[N-i+1]
37 apply_inverse_quadrature(op::D2{T}, h_inv::Real, v::T, i::Integer, N::Integer, ::Type{Interior}) where T = v*h_inv
38
39 # TODO: Avoid branching in inner loops
40 function apply_inverse_quadrature(op::D2{T}, h_inv::Real, v::T, i::Integer, N::Integer) where T
41 r = getregion(i, closuresize(op), N)
42 return apply_inverse_quadrature(op, h_inv, v, i, N, r)
43 end
44 export apply_inverse_quadrature
45
46 function apply_e_T(op::D2, v::AbstractVector, ::Type{Lower})
47 @boundscheck if length(v) < closuresize(op)
48 throw(BoundsError())
49 end
50 apply(op.eClosure,v,1)
51 end
52
53 function apply_e_T(op::D2, v::AbstractVector, ::Type{Upper})
54 @boundscheck if length(v) < closuresize(op)
55 throw(BoundsError())
56 end
57 apply(flip(op.eClosure),v,length(v))
58 end
59
60
61 function apply_e(op::D2, v::Number, N::Integer, i::Integer, ::Type{Lower})
62 @boundscheck if !(0<length(i) <= N)
63 throw(BoundsError())
64 end
65 op.eClosure[i-1]*v
66 end
67
68 function apply_e(op::D2, v::Number, N::Integer, i::Integer, ::Type{Upper})
69 @boundscheck if !(0<length(i) <= N)
70 throw(BoundsError())
71 end
72 op.eClosure[N-i]*v
73 end
74
75 function apply_d_T(op::D2, h_inv::Real, v::AbstractVector, ::Type{Lower})
76 @boundscheck if length(v) < closuresize(op)
77 throw(BoundsError())
78 end
79 h_inv*apply(op.dClosure,v,1)
80 end
81
82 function apply_d_T(op::D2, h_inv::Real, v::AbstractVector, ::Type{Upper})
83 @boundscheck if length(v) < closuresize(op)
84 throw(BoundsError())
85 end
86 -h_inv*apply(flip(op.dClosure),v,length(v))
87 end
88
89 function apply_d(op::D2, h_inv::Real, v::Number, N::Integer, i::Integer, ::Type{Lower})
90 @boundscheck if !(0<length(i) <= N)
91 throw(BoundsError())
92 end
93 h_inv*op.dClosure[i-1]*v
94 end
95
96 function apply_d(op::D2, h_inv::Real, v::Number, N::Integer, i::Integer, ::Type{Upper})
97 @boundscheck if !(0<length(i) <= N)
98 throw(BoundsError())
99 end
100 -h_inv*op.dClosure[N-i]*v
101 end