annotate src/SbpOperators/stencil.jl @ 1063:f98893154e22 feature/nested_stencils

Fix type instability for stencil application when weights and vector didn't match (grafted from a7f898b1ce1ea3adfe21201d3ccafdffb09c5be6)
author Jonatan Werpers <jonatan@werpers.com>
date Mon, 14 Feb 2022 14:28:48 +0100
parents a76830879c63
children a3bc90c59e8e
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
671
e14627e79a54 Add stencil constructor for centered stencils and change from tuple to vararg in stencil constructor taking cneter
Jonatan Werpers <jonatan@werpers.com>
parents: 595
diff changeset
1 export CenteredStencil
1052
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
2 export CenteredNestedStencil
671
e14627e79a54 Add stencil constructor for centered stencils and change from tuple to vararg in stencil constructor taking cneter
Jonatan Werpers <jonatan@werpers.com>
parents: 595
diff changeset
3
806
9fc6d38da03f Remove type requirement for stencil weights
Jonatan Werpers <jonatan@werpers.com>
parents: 672
diff changeset
4 struct Stencil{T,N}
1062
a76830879c63 Fix bugs from stencil refactor
Jonatan Werpers <jonatan@werpers.com>
parents: 1061
diff changeset
5 range::UnitRange{Int64}
84
48079bd39969 Change to using tuples in stencils and ops
Jonatan Werpers <jonatan@werpers.com>
parents: 69
diff changeset
6 weights::NTuple{N,T}
126
66c239678a21 Add Assertion in stencil constructor
Ylva Rydin <ylva.rydin@telia.com>
parents: 85
diff changeset
7
1060
ff0e819f2075 Refactor code for regular stencils to use fewer type parameters and allow promotion
Jonatan Werpers <jonatan@werpers.com>
parents: 1059
diff changeset
8 function Stencil(range::UnitRange,weights::NTuple{N,T}) where {T, N}
ff0e819f2075 Refactor code for regular stencils to use fewer type parameters and allow promotion
Jonatan Werpers <jonatan@werpers.com>
parents: 1059
diff changeset
9 @assert length(range) == N
126
66c239678a21 Add Assertion in stencil constructor
Ylva Rydin <ylva.rydin@telia.com>
parents: 85
diff changeset
10 new{T,N}(range,weights)
66c239678a21 Add Assertion in stencil constructor
Ylva Rydin <ylva.rydin@telia.com>
parents: 85
diff changeset
11 end
8
433008d3b7d3 Move stencil to its own file
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
12 end
433008d3b7d3 Move stencil to its own file
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
13
584
4aa7fe13a984 Add scale() and eltype() methods for stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 333
diff changeset
14 """
595
03ef4d4740ab Add a constructor for Stencil where you can specify the center of the stencil
Jonatan Werpers <jonatan@werpers.com>
parents: 584
diff changeset
15 Stencil(weights::NTuple; center::Int)
03ef4d4740ab Add a constructor for Stencil where you can specify the center of the stencil
Jonatan Werpers <jonatan@werpers.com>
parents: 584
diff changeset
16
03ef4d4740ab Add a constructor for Stencil where you can specify the center of the stencil
Jonatan Werpers <jonatan@werpers.com>
parents: 584
diff changeset
17 Create a stencil with the given weights with element `center` as the center of the stencil.
03ef4d4740ab Add a constructor for Stencil where you can specify the center of the stencil
Jonatan Werpers <jonatan@werpers.com>
parents: 584
diff changeset
18 """
1060
ff0e819f2075 Refactor code for regular stencils to use fewer type parameters and allow promotion
Jonatan Werpers <jonatan@werpers.com>
parents: 1059
diff changeset
19 function Stencil(weights...; center::Int)
ff0e819f2075 Refactor code for regular stencils to use fewer type parameters and allow promotion
Jonatan Werpers <jonatan@werpers.com>
parents: 1059
diff changeset
20 weights = promote(weights...)
595
03ef4d4740ab Add a constructor for Stencil where you can specify the center of the stencil
Jonatan Werpers <jonatan@werpers.com>
parents: 584
diff changeset
21 N = length(weights)
1060
ff0e819f2075 Refactor code for regular stencils to use fewer type parameters and allow promotion
Jonatan Werpers <jonatan@werpers.com>
parents: 1059
diff changeset
22 range = (1:N) .- center
595
03ef4d4740ab Add a constructor for Stencil where you can specify the center of the stencil
Jonatan Werpers <jonatan@werpers.com>
parents: 584
diff changeset
23
03ef4d4740ab Add a constructor for Stencil where you can specify the center of the stencil
Jonatan Werpers <jonatan@werpers.com>
parents: 584
diff changeset
24 return Stencil(range, weights)
03ef4d4740ab Add a constructor for Stencil where you can specify the center of the stencil
Jonatan Werpers <jonatan@werpers.com>
parents: 584
diff changeset
25 end
03ef4d4740ab Add a constructor for Stencil where you can specify the center of the stencil
Jonatan Werpers <jonatan@werpers.com>
parents: 584
diff changeset
26
1060
ff0e819f2075 Refactor code for regular stencils to use fewer type parameters and allow promotion
Jonatan Werpers <jonatan@werpers.com>
parents: 1059
diff changeset
27 Stencil{T,N}(s::Stencil{S,N}) where {T,S,N} = Stencil(s.range, T.(s.weights))
ff0e819f2075 Refactor code for regular stencils to use fewer type parameters and allow promotion
Jonatan Werpers <jonatan@werpers.com>
parents: 1059
diff changeset
28 Stencil{T}(s::Stencil) where T = Stencil{T,length(s)}(s)
826
4433be383840 Add stencil constructor to change the type of the weights along with a convert method
Jonatan Werpers <jonatan@werpers.com>
parents: 807
diff changeset
29
1060
ff0e819f2075 Refactor code for regular stencils to use fewer type parameters and allow promotion
Jonatan Werpers <jonatan@werpers.com>
parents: 1059
diff changeset
30 Base.convert(::Type{Stencil{T1,N}}, s::Stencil{T2,N}) where {T1,T2,N} = Stencil{T1,N}(s)
ff0e819f2075 Refactor code for regular stencils to use fewer type parameters and allow promotion
Jonatan Werpers <jonatan@werpers.com>
parents: 1059
diff changeset
31 Base.convert(::Type{Stencil{T1}}, s::Stencil{T2,N}) where {T1,T2,N} = Stencil{T1,N}(s)
826
4433be383840 Add stencil constructor to change the type of the weights along with a convert method
Jonatan Werpers <jonatan@werpers.com>
parents: 807
diff changeset
32
1060
ff0e819f2075 Refactor code for regular stencils to use fewer type parameters and allow promotion
Jonatan Werpers <jonatan@werpers.com>
parents: 1059
diff changeset
33 Base.promote_rule(::Type{Stencil{T1,N}}, ::Type{Stencil{T2,N}}) where {T1,T2,N} = Stencil{promote_type(T1,T2),N}
ff0e819f2075 Refactor code for regular stencils to use fewer type parameters and allow promotion
Jonatan Werpers <jonatan@werpers.com>
parents: 1059
diff changeset
34
ff0e819f2075 Refactor code for regular stencils to use fewer type parameters and allow promotion
Jonatan Werpers <jonatan@werpers.com>
parents: 1059
diff changeset
35 function CenteredStencil(weights...)
671
e14627e79a54 Add stencil constructor for centered stencils and change from tuple to vararg in stencil constructor taking cneter
Jonatan Werpers <jonatan@werpers.com>
parents: 595
diff changeset
36 if iseven(length(weights))
e14627e79a54 Add stencil constructor for centered stencils and change from tuple to vararg in stencil constructor taking cneter
Jonatan Werpers <jonatan@werpers.com>
parents: 595
diff changeset
37 throw(ArgumentError("a centered stencil must have an odd number of weights."))
e14627e79a54 Add stencil constructor for centered stencils and change from tuple to vararg in stencil constructor taking cneter
Jonatan Werpers <jonatan@werpers.com>
parents: 595
diff changeset
38 end
e14627e79a54 Add stencil constructor for centered stencils and change from tuple to vararg in stencil constructor taking cneter
Jonatan Werpers <jonatan@werpers.com>
parents: 595
diff changeset
39
e14627e79a54 Add stencil constructor for centered stencils and change from tuple to vararg in stencil constructor taking cneter
Jonatan Werpers <jonatan@werpers.com>
parents: 595
diff changeset
40 r = length(weights) ÷ 2
e14627e79a54 Add stencil constructor for centered stencils and change from tuple to vararg in stencil constructor taking cneter
Jonatan Werpers <jonatan@werpers.com>
parents: 595
diff changeset
41
1060
ff0e819f2075 Refactor code for regular stencils to use fewer type parameters and allow promotion
Jonatan Werpers <jonatan@werpers.com>
parents: 1059
diff changeset
42 return Stencil(-r:r, weights)
671
e14627e79a54 Add stencil constructor for centered stencils and change from tuple to vararg in stencil constructor taking cneter
Jonatan Werpers <jonatan@werpers.com>
parents: 595
diff changeset
43 end
e14627e79a54 Add stencil constructor for centered stencils and change from tuple to vararg in stencil constructor taking cneter
Jonatan Werpers <jonatan@werpers.com>
parents: 595
diff changeset
44
672
59a81254fefc Formatting
Jonatan Werpers <jonatan@werpers.com>
parents: 671
diff changeset
45
595
03ef4d4740ab Add a constructor for Stencil where you can specify the center of the stencil
Jonatan Werpers <jonatan@werpers.com>
parents: 584
diff changeset
46 """
584
4aa7fe13a984 Add scale() and eltype() methods for stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 333
diff changeset
47 scale(s::Stencil, a)
4aa7fe13a984 Add scale() and eltype() methods for stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 333
diff changeset
48
4aa7fe13a984 Add scale() and eltype() methods for stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 333
diff changeset
49 Scale the weights of the stencil `s` with `a` and return a new stencil.
4aa7fe13a984 Add scale() and eltype() methods for stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 333
diff changeset
50 """
4aa7fe13a984 Add scale() and eltype() methods for stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 333
diff changeset
51 function scale(s::Stencil, a)
4aa7fe13a984 Add scale() and eltype() methods for stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 333
diff changeset
52 return Stencil(s.range, a.*s.weights)
4aa7fe13a984 Add scale() and eltype() methods for stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 333
diff changeset
53 end
4aa7fe13a984 Add scale() and eltype() methods for stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 333
diff changeset
54
1057
0728e84af0dc Add length method for stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 1056
diff changeset
55 Base.eltype(::Stencil{T,N}) where {T,N} = T
0728e84af0dc Add length method for stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 1056
diff changeset
56 Base.length(::Stencil{T,N}) where {T,N} = N
584
4aa7fe13a984 Add scale() and eltype() methods for stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 333
diff changeset
57
8
433008d3b7d3 Move stencil to its own file
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
58 function flip(s::Stencil)
433008d3b7d3 Move stencil to its own file
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
59 range = (-s.range[2], -s.range[1])
84
48079bd39969 Change to using tuples in stencils and ops
Jonatan Werpers <jonatan@werpers.com>
parents: 69
diff changeset
60 return Stencil(range, reverse(s.weights))
8
433008d3b7d3 Move stencil to its own file
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
61 end
433008d3b7d3 Move stencil to its own file
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
62
433008d3b7d3 Move stencil to its own file
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
63 # Provides index into the Stencil based on offset for the root element
128
7c0b9bb7ab4d Improve stencil application code to make it more friendly to compiler optimizations
Jonatan Werpers <jonatan@werpers.com>
parents: 122
diff changeset
64 @inline function Base.getindex(s::Stencil, i::Int)
1060
ff0e819f2075 Refactor code for regular stencils to use fewer type parameters and allow promotion
Jonatan Werpers <jonatan@werpers.com>
parents: 1059
diff changeset
65 @boundscheck if i ∉ s.range
584
4aa7fe13a984 Add scale() and eltype() methods for stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 333
diff changeset
66 return zero(eltype(s))
8
433008d3b7d3 Move stencil to its own file
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
67 end
71
18d0d794d3bb Make stencils respond to @ inbounds
Jonatan Werpers <jonatan@werpers.com>
parents: 67
diff changeset
68 return s.weights[1 + i - s.range[1]]
8
433008d3b7d3 Move stencil to its own file
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
69 end
433008d3b7d3 Move stencil to its own file
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
70
1060
ff0e819f2075 Refactor code for regular stencils to use fewer type parameters and allow promotion
Jonatan Werpers <jonatan@werpers.com>
parents: 1059
diff changeset
71 Base.@propagate_inbounds @inline function apply_stencil(s::Stencil, v::AbstractVector, i::Int)
1063
f98893154e22 Fix type instability for stencil application when weights and vector didn't match
Jonatan Werpers <jonatan@werpers.com>
parents: 1062
diff changeset
72 w = zero(promote_type(eltype(s),eltype(v)))
1060
ff0e819f2075 Refactor code for regular stencils to use fewer type parameters and allow promotion
Jonatan Werpers <jonatan@werpers.com>
parents: 1059
diff changeset
73 @simd for k ∈ 1:length(s)
ff0e819f2075 Refactor code for regular stencils to use fewer type parameters and allow promotion
Jonatan Werpers <jonatan@werpers.com>
parents: 1059
diff changeset
74 w += s.weights[k]*v[i + s.range[k]]
8
433008d3b7d3 Move stencil to its own file
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
75 end
1060
ff0e819f2075 Refactor code for regular stencils to use fewer type parameters and allow promotion
Jonatan Werpers <jonatan@werpers.com>
parents: 1059
diff changeset
76
8
433008d3b7d3 Move stencil to its own file
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
77 return w
433008d3b7d3 Move stencil to its own file
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
78 end
122
6c6979ff17f4 Introduce and use apply_backwards for stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 85
diff changeset
79
1060
ff0e819f2075 Refactor code for regular stencils to use fewer type parameters and allow promotion
Jonatan Werpers <jonatan@werpers.com>
parents: 1059
diff changeset
80 Base.@propagate_inbounds @inline function apply_stencil_backwards(s::Stencil, v::AbstractVector, i::Int)
ff0e819f2075 Refactor code for regular stencils to use fewer type parameters and allow promotion
Jonatan Werpers <jonatan@werpers.com>
parents: 1059
diff changeset
81 w = zero(eltype(v))
ff0e819f2075 Refactor code for regular stencils to use fewer type parameters and allow promotion
Jonatan Werpers <jonatan@werpers.com>
parents: 1059
diff changeset
82 @simd for k ∈ length(s):-1:1
ff0e819f2075 Refactor code for regular stencils to use fewer type parameters and allow promotion
Jonatan Werpers <jonatan@werpers.com>
parents: 1059
diff changeset
83 w += s.weights[k]*v[i - s.range[k]]
122
6c6979ff17f4 Introduce and use apply_backwards for stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 85
diff changeset
84 end
6c6979ff17f4 Introduce and use apply_backwards for stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 85
diff changeset
85 return w
6c6979ff17f4 Introduce and use apply_backwards for stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 85
diff changeset
86 end
1052
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
87
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
88
1056
fe83ea1db953 Add argument restrictions to stencil constructors to get clearer error messages
Jonatan Werpers <jonatan@werpers.com>
parents: 1055
diff changeset
89 struct NestedStencil{T,N,M}
fe83ea1db953 Add argument restrictions to stencil constructors to get clearer error messages
Jonatan Werpers <jonatan@werpers.com>
parents: 1055
diff changeset
90 s::Stencil{Stencil{T,N},M}
1052
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
91 end
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
92
1055
df498ce0cf52 Add methods for CenteredNestedStencil with tuple input
Jonatan Werpers <jonatan@werpers.com>
parents: 1053
diff changeset
93 # Stencil input
1061
d1ddcdc098bf Add promotion functionality for NestedTuples
Jonatan Werpers <jonatan@werpers.com>
parents: 1060
diff changeset
94 NestedStencil(s::Vararg{Stencil}; center) = NestedStencil(Stencil(s... ; center))
d1ddcdc098bf Add promotion functionality for NestedTuples
Jonatan Werpers <jonatan@werpers.com>
parents: 1060
diff changeset
95 CenteredNestedStencil(s::Vararg{Stencil}) = NestedStencil(CenteredStencil(s...))
1052
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
96
1055
df498ce0cf52 Add methods for CenteredNestedStencil with tuple input
Jonatan Werpers <jonatan@werpers.com>
parents: 1053
diff changeset
97 # Tuple input
1061
d1ddcdc098bf Add promotion functionality for NestedTuples
Jonatan Werpers <jonatan@werpers.com>
parents: 1060
diff changeset
98 function NestedStencil(weights::Vararg{NTuple{N,Any}}; center) where N
1053
cb0a5f41be02 Add some tests and fix some bugs
Jonatan Werpers <jonatan@werpers.com>
parents: 1052
diff changeset
99 inner_stencils = map(w -> Stencil(w...; center), weights)
1052
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
100 return NestedStencil(Stencil(inner_stencils... ; center))
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
101 end
1061
d1ddcdc098bf Add promotion functionality for NestedTuples
Jonatan Werpers <jonatan@werpers.com>
parents: 1060
diff changeset
102 function CenteredNestedStencil(weights::Vararg{NTuple{N,Any}}) where N
1055
df498ce0cf52 Add methods for CenteredNestedStencil with tuple input
Jonatan Werpers <jonatan@werpers.com>
parents: 1053
diff changeset
103 inner_stencils = map(w->CenteredStencil(w...), weights)
df498ce0cf52 Add methods for CenteredNestedStencil with tuple input
Jonatan Werpers <jonatan@werpers.com>
parents: 1053
diff changeset
104 return CenteredNestedStencil(inner_stencils...)
df498ce0cf52 Add methods for CenteredNestedStencil with tuple input
Jonatan Werpers <jonatan@werpers.com>
parents: 1053
diff changeset
105 end
1052
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
106
1058
c4f71d6f2d63 Add functions for converting the element type of nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 1057
diff changeset
107
c4f71d6f2d63 Add functions for converting the element type of nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 1057
diff changeset
108 # Conversion
1061
d1ddcdc098bf Add promotion functionality for NestedTuples
Jonatan Werpers <jonatan@werpers.com>
parents: 1060
diff changeset
109 function NestedStencil{T,N,M}(ns::NestedStencil{S,N,M}) where {T,S,N,M}
1058
c4f71d6f2d63 Add functions for converting the element type of nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 1057
diff changeset
110 return NestedStencil(Stencil{Stencil{T}}(ns.s))
c4f71d6f2d63 Add functions for converting the element type of nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 1057
diff changeset
111 end
c4f71d6f2d63 Add functions for converting the element type of nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 1057
diff changeset
112
1061
d1ddcdc098bf Add promotion functionality for NestedTuples
Jonatan Werpers <jonatan@werpers.com>
parents: 1060
diff changeset
113 function NestedStencil{T}(ns::NestedStencil{S,N,M}) where {T,S,N,M}
d1ddcdc098bf Add promotion functionality for NestedTuples
Jonatan Werpers <jonatan@werpers.com>
parents: 1060
diff changeset
114 NestedStencil{T,N,M}(ns)
d1ddcdc098bf Add promotion functionality for NestedTuples
Jonatan Werpers <jonatan@werpers.com>
parents: 1060
diff changeset
115 end
d1ddcdc098bf Add promotion functionality for NestedTuples
Jonatan Werpers <jonatan@werpers.com>
parents: 1060
diff changeset
116
d1ddcdc098bf Add promotion functionality for NestedTuples
Jonatan Werpers <jonatan@werpers.com>
parents: 1060
diff changeset
117 function Base.convert(::Type{NestedStencil{T,N,M}}, s::NestedStencil{S,N,M}) where {T,S,N,M}
d1ddcdc098bf Add promotion functionality for NestedTuples
Jonatan Werpers <jonatan@werpers.com>
parents: 1060
diff changeset
118 return NestedStencil{T,N,M}(s)
d1ddcdc098bf Add promotion functionality for NestedTuples
Jonatan Werpers <jonatan@werpers.com>
parents: 1060
diff changeset
119 end
1058
c4f71d6f2d63 Add functions for converting the element type of nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 1057
diff changeset
120 Base.convert(::Type{NestedStencil{T}}, stencil) where T = NestedStencil{T}(stencil)
c4f71d6f2d63 Add functions for converting the element type of nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 1057
diff changeset
121
1061
d1ddcdc098bf Add promotion functionality for NestedTuples
Jonatan Werpers <jonatan@werpers.com>
parents: 1060
diff changeset
122 function Base.promote_rule(::Type{NestedStencil{T,N,M}}, ::Type{NestedStencil{S,N,M}}) where {T,S,N,M}
d1ddcdc098bf Add promotion functionality for NestedTuples
Jonatan Werpers <jonatan@werpers.com>
parents: 1060
diff changeset
123 return NestedStencil{promote_type(T,S),N,M}
d1ddcdc098bf Add promotion functionality for NestedTuples
Jonatan Werpers <jonatan@werpers.com>
parents: 1060
diff changeset
124 end
d1ddcdc098bf Add promotion functionality for NestedTuples
Jonatan Werpers <jonatan@werpers.com>
parents: 1060
diff changeset
125
1052
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
126 Base.eltype(::NestedStencil{T}) where T = T
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
127
1059
4d06642174ec Add scale method for nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 1058
diff changeset
128 function scale(ns::NestedStencil, a)
4d06642174ec Add scale method for nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 1058
diff changeset
129 range = ns.s.range
4d06642174ec Add scale method for nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 1058
diff changeset
130 weights = ns.s.weights
4d06642174ec Add scale method for nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 1058
diff changeset
131
4d06642174ec Add scale method for nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 1058
diff changeset
132 return NestedStencil(Stencil(range, scale.(weights,a)))
4d06642174ec Add scale method for nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 1058
diff changeset
133 end
4d06642174ec Add scale method for nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 1058
diff changeset
134
1052
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
135 function flip(ns::NestedStencil)
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
136 s_flip = flip(ns.s)
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
137 return NestedStencil(Stencil(s_flip.range, flip.(s_flip.weights)))
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
138 end
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
139
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
140 Base.getindex(ns::NestedStencil, i::Int) = ns.s[i]
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
141
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
142 "Apply inner stencils to `c` and get a concrete stencil"
1053
cb0a5f41be02 Add some tests and fix some bugs
Jonatan Werpers <jonatan@werpers.com>
parents: 1052
diff changeset
143 Base.@propagate_inbounds function apply_inner_stencils(ns::NestedStencil, c::AbstractVector, i::Int)
1052
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
144 weights = apply_stencil.(ns.s.weights, Ref(c), i)
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
145 return Stencil(ns.s.range, weights)
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
146 end
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
147
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
148 "Apply the whole nested stencil"
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
149 Base.@propagate_inbounds function apply_stencil(ns::NestedStencil, c::AbstractVector, v::AbstractVector, i::Int)
1053
cb0a5f41be02 Add some tests and fix some bugs
Jonatan Werpers <jonatan@werpers.com>
parents: 1052
diff changeset
150 s = apply_inner_stencils(ns,c,i)
1052
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
151 return apply_stencil(s, v, i)
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
152 end
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
153
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
154 "Apply inner stencils backwards to `c` and get a concrete stencil"
1053
cb0a5f41be02 Add some tests and fix some bugs
Jonatan Werpers <jonatan@werpers.com>
parents: 1052
diff changeset
155 Base.@propagate_inbounds @inline function apply_inner_stencils_backwards(ns::NestedStencil, c::AbstractVector, i::Int)
1052
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
156 weights = apply_stencil_backwards.(ns.s.weights, Ref(c), i)
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
157 return Stencil(ns.s.range, weights)
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
158 end
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
159
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
160 "Apply the whole nested stencil backwards"
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
161 Base.@propagate_inbounds @inline function apply_stencil_backwards(ns::NestedStencil, c::AbstractVector, v::AbstractVector, i::Int)
1053
cb0a5f41be02 Add some tests and fix some bugs
Jonatan Werpers <jonatan@werpers.com>
parents: 1052
diff changeset
162 s = apply_inner_stencils_backwards(ns,c,i)
1052
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
163 return apply_stencil_backwards(s, v, i)
ca718fd4e816 Start implementing nested stencils
Jonatan Werpers <jonatan@werpers.com>
parents: 826
diff changeset
164 end