comparison SbpOperators/src/BoundaryValue.jl @ 306:f8a4850caed2

Move BoundaryValue from Laplace to separate file. Currently WIP
author Vidar Stiernström <vidar.stiernstrom@it.uu.se>
date Wed, 09 Sep 2020 21:06:27 +0200
parents
children
comparison
equal deleted inserted replaced
305:bd09d67ebb22 306:f8a4850caed2
1 """
2 BoundaryValue{T,N,M,K} <: TensorMapping{T,2,1}
3
4 Implements the boundary operator `e` as a TensorMapping
5 """
6 struct BoundaryValue{T,N,M,K} <: TensorMapping{T,2,1}
7 eClosure::Stencil{T,M}
8 bId::CartesianBoundary
9 end
10 export BoundaryValue
11
12 # TODO: This is obviouly strange. Is domain_size just discarded? Is there a way to avoid storing grid in BoundaryValue?
13 # Can we give special treatment to TensorMappings that go to a higher dim?
14 function LazyTensors.range_size(e::BoundaryValue{T}, domain_size::NTuple{1,Integer}) where T
15 if dim(e.bId) == 1
16 return (UnknownDim, domain_size[1])
17 elseif dim(e.bId) == 2
18 return (domain_size[1], UnknownDim)
19 end
20 end
21 LazyTensors.domain_size(e::BoundaryValue{T}, range_size::NTuple{2,Integer}) where T = (range_size[3-dim(e.bId)],)
22 # TODO: Make a nicer solution for 3-dim(e.bId)
23
24 # TODO: Make this independent of dimension
25 function LazyTensors.apply(e::BoundaryValue{T}, v::AbstractArray{T}, I::NTuple{2,Index}) where T
26 i = I[dim(e.bId)]
27 j = I[3-dim(e.bId)]
28 N_i = size(e.grid)[dim(e.bId)]
29 return apply_boundary_value(e.op, v[j], i, N_i, region(e.bId))
30 end
31
32 function LazyTensors.apply_transpose(e::BoundaryValue{T}, v::AbstractArray{T}, I::NTuple{1,Index}) where T
33 u = selectdim(v,3-dim(e.bId),Int(I[1]))
34 return apply_boundary_value_transpose(e.op, u, region(e.bId))
35 end
36
37 function apply_boundary_value_transpose(op::ConstantStencilOperator, v::AbstractVector, ::Type{Lower})
38 @boundscheck if length(v) < closuresize(op)
39 throw(BoundsError())
40 end
41 apply_stencil(op.eClosure,v,1)
42 end
43
44 function apply_boundary_value_transpose(op::ConstantStencilOperator, v::AbstractVector, ::Type{Upper})
45 @boundscheck if length(v) < closuresize(op)
46 throw(BoundsError())
47 end
48 apply_stencil_backwards(op.eClosure,v,length(v))
49 end
50 export apply_boundary_value_transpose
51
52 function apply_boundary_value(op::ConstantStencilOperator, v::Number, i::Index, N::Integer, ::Type{Lower})
53 @boundscheck if !(0<length(Int(i)) <= N)
54 throw(BoundsError())
55 end
56 op.eClosure[Int(i)-1]*v
57 end
58
59 function apply_boundary_value(op::ConstantStencilOperator, v::Number, i::Index, N::Integer, ::Type{Upper})
60 @boundscheck if !(0<length(Int(i)) <= N)
61 throw(BoundsError())
62 end
63 op.eClosure[N-Int(i)]*v
64 end
65 export apply_boundary_value
66
67
68 """
69 BoundaryValue{T,N,M,K} <: TensorMapping{T,2,1}
70
71 Implements the boundary operator `e` as a TensorMapping
72 """
73 struct BoundaryValue{D,T,M,R} <: TensorMapping{T,D,1}
74 e:BoundaryOperator{T,M,R}
75 bId::CartesianBoundary
76 end
77
78 function LazyTensors.apply_transpose(bv::BoundaryValue{T,M,Lower}, v::AbstractVector{T}, i::Index) where T
79 u = selectdim(v,3-dim(bv.bId),Int(I[1]))
80 return apply_transpose(bv.e, u, I)
81 end
82
83
84 """
85 BoundaryOperator{T,N,R} <: TensorMapping{T,1,1}
86
87 Implements the boundary operator `e` as a TensorMapping
88 """
89 export BoundaryOperator
90 struct BoundaryOperator{T,M,R<:Region} <: TensorMapping{T,1,1}
91 closure::Stencil{T,M}
92 end
93
94 function LazyTensors.range_size(e::BoundaryOperator, domain_size::NTuple{1,Integer})
95 return UnknownDim
96 end
97
98 LazyTensors.domain_size(e::BoundaryOperator{T}, range_size::NTuple{1,Integer}) where T = range_size
99
100 function LazyTensors.apply_transpose(e::BoundaryOperator{T,M,Lower}, v::AbstractVector{T}, i::Index{Lower}) where T
101 @boundscheck if length(v) < closuresize(e) #TODO: Use domain_size here?
102 throw(BoundsError())
103 end
104 apply_stencil(e.closure,v,Int(i))
105 end
106
107 function LazyTensors.apply_transpose(e::BoundaryOperator{T,M,Upper}}, v::AbstractVector{T}, i::Index{Upper}) where T
108 @boundscheck if length(v) < closuresize(e) #TODO: Use domain_size here?
109 throw(BoundsError())
110 end
111 apply_stencil_backwards(e.closure,v,Int(i))
112 end
113
114 function LazyTensors.apply_transpose(e::BoundaryOperator{T}, v::AbstractVector{T}, i::Index) where T
115 @boundscheck if length(v) < closuresize(e) #TODO: Use domain_size here?
116 throw(BoundsError())
117 end
118 return eltype(v)(0)
119 end
120
121 #TODO: Implement apply in a meaningful way. Should it return a vector or a single value (perferable?) Should fit into the
122 function LazyTensors.apply(e::BoundaryOperator, v::AbstractVector, i::Index)
123 @boundscheck if !(0<length(Int(i)) <= length(v))
124 throw(BoundsError())
125 end
126 return e.closure[Int(i)].*v
127 end