comparison test/SbpOperators/volumeops/derivatives/dissipation_test.jl @ 1858:4a9be96f2569 feature/documenter_logo

Merge default
author Jonatan Werpers <jonatan@werpers.com>
date Sun, 12 Jan 2025 21:18:44 +0100
parents 471a948cd2b2
children
comparison
equal deleted inserted replaced
1857:ffde7dad9da5 1858:4a9be96f2569
1 using Test
2
3 using Diffinitive.SbpOperators
4 using Diffinitive.Grids
5 using Diffinitive.LazyTensors
6
7 using Diffinitive.SbpOperators: Stencil
8
9 using Diffinitive.SbpOperators: dissipation_interior_weights
10 using Diffinitive.SbpOperators: dissipation_interior_stencil, dissipation_transpose_interior_stencil
11 using Diffinitive.SbpOperators: midpoint, midpoint_transpose
12 using Diffinitive.SbpOperators: dissipation_lower_closure_size, dissipation_upper_closure_size
13 using Diffinitive.SbpOperators: dissipation_lower_closure_stencils,dissipation_upper_closure_stencils
14 using Diffinitive.SbpOperators: dissipation_transpose_lower_closure_stencils, dissipation_transpose_upper_closure_stencils
15
16
17 @testset "undivided_skewed04" begin
18 monomial(x,k) = k < 0 ? zero(x) : x^k/factorial(k)
19 g = equidistant_grid(0., 11., 20)
20 D,Dᵀ = undivided_skewed04(g, 1)
21
22 @test D isa LazyTensor{Float64,1,1}
23 @test Dᵀ isa LazyTensor{Float64,1,1}
24
25 @testset "Accuracy conditions" begin
26 N = 20
27 g = equidistant_grid(0//1, 2//1, N)
28 h = only(spacing(g))
29 @testset "D_$p" for p ∈ [1,2,3,4]
30 D,Dᵀ = undivided_skewed04(g, p)
31
32 @testset "x^$k" for k ∈ 0:p
33 v = eval_on(g, x->monomial(x,k))
34 vₚₓ = eval_on(g, x->monomial(x,k-p))
35
36 @test D*v == h^p * vₚₓ
37 end
38 end
39 end
40
41 @testset "transpose equality" begin
42 function get_matrix(D)
43 N = only(range_size(D))
44 M = only(domain_size(D))
45
46 Dmat = zeros(N,M)
47 e = zeros(M)
48 for i ∈ 1:M
49 if i > 1
50 e[i-1] = 0.
51 end
52 e[i] = 1.
53 Dmat[:,i] = D*e
54 end
55
56 return Dmat
57 end
58
59 g = equidistant_grid(0., 1., 11)
60 @testset "D_$p" for p ∈ [1,2,3,4]
61 D,Dᵀ = undivided_skewed04(g, p)
62
63 D̄ = get_matrix(D)
64 D̄ᵀ = get_matrix(Dᵀ)
65
66 @test D̄ == D̄ᵀ'
67 end
68 end
69
70 @testset "2D" begin
71 N = 20
72 g = equidistant_grid((0,0), (2,1), N, 2N)
73 h = spacing.(g.grids)
74
75 D,Dᵀ = undivided_skewed04(g, 3, 2)
76
77 v = eval_on(g, x->monomial(x[1],4)*monomial(x[2],3))
78 d³vdy³ = eval_on(g, x->monomial(x[1],4)*monomial(x[2],0))
79
80 @test D*v ≈ h[2]^3*d³vdy³
81 end
82 end
83
84 @testset "dissipation_interior_weights" begin
85 @test dissipation_interior_weights(1) == (-1, 1)
86 @test dissipation_interior_weights(2) == (1,-2, 1)
87 @test dissipation_interior_weights(3) == (-1, 3,-3, 1)
88 @test dissipation_interior_weights(4) == (1, -4, 6, -4, 1)
89 end
90
91 @testset "dissipation_interior_stencil" begin
92 @test dissipation_interior_stencil(dissipation_interior_weights(1)) == Stencil(-1, 1, center=2)
93 @test dissipation_interior_stencil(dissipation_interior_weights(2)) == Stencil( 1,-2, 1, center=2)
94 @test dissipation_interior_stencil(dissipation_interior_weights(3)) == Stencil(-1, 3,-3, 1, center=3)
95 @test dissipation_interior_stencil(dissipation_interior_weights(4)) == Stencil( 1,-4, 6,-4, 1, center=3)
96 end
97
98 @testset "dissipation_transpose_interior_stencil" begin
99 @test dissipation_transpose_interior_stencil(dissipation_interior_weights(1)) == Stencil(1,-1, center=1)
100 @test dissipation_transpose_interior_stencil(dissipation_interior_weights(2)) == Stencil(1,-2, 1, center=2)
101 @test dissipation_transpose_interior_stencil(dissipation_interior_weights(3)) == Stencil(1,-3, 3,-1, center=2)
102 @test dissipation_transpose_interior_stencil(dissipation_interior_weights(4)) == Stencil(1,-4, 6,-4, 1, center=3)
103 end
104
105 @testset "midpoint" begin
106 @test midpoint((1,1)) == 2
107 @test midpoint((1,1,1)) == 2
108 @test midpoint((1,1,1,1)) == 3
109 @test midpoint((1,1,1,1,1)) == 3
110 end
111
112 @testset "midpoint_transpose" begin
113 @test midpoint_transpose((1,1)) == 1
114 @test midpoint_transpose((1,1,1)) == 2
115 @test midpoint_transpose((1,1,1,1)) == 2
116 @test midpoint_transpose((1,1,1,1,1)) == 3
117 end
118
119 @testset "dissipation_lower_closure_size" begin
120 @test dissipation_lower_closure_size((1,1)) == 1
121 @test dissipation_lower_closure_size((1,1,1)) == 1
122 @test dissipation_lower_closure_size((1,1,1,1)) == 2
123 @test dissipation_lower_closure_size((1,1,1,1,1)) == 2
124 end
125
126 @testset "dissipation_upper_closure_size" begin
127 @test dissipation_upper_closure_size((1,1)) == 0
128 @test dissipation_upper_closure_size((1,1,1)) == 1
129 @test dissipation_upper_closure_size((1,1,1,1)) == 1
130 @test dissipation_upper_closure_size((1,1,1,1,1)) == 2
131 end
132
133 @testset "dissipation_lower_closure_stencils" begin
134 cases = (
135 (-1,1) => (
136 Stencil(-1, 1, center=1),
137 ),
138 (1,-2,1) => (
139 Stencil( 1,-2, 1, center=1),
140 ),
141 (-1,3,-3,1) => (
142 Stencil(-1,3,-3,1, center=1),
143 Stencil(-1,3,-3,1, center=2),
144 ),
145 (1, -4, 6, -4, 1) => (
146 Stencil(1, -4, 6, -4, 1, center=1),
147 Stencil(1, -4, 6, -4, 1, center=2),
148 )
149 )
150 @testset "interior_weights = $w" for (w, closure_stencils) ∈ cases
151 @test dissipation_lower_closure_stencils(w) == closure_stencils
152 end
153 end
154
155 @testset "dissipation_upper_closure_stencils" begin
156 cases = (
157 (-1,1) => (),
158 (1,-2,1) => (
159 Stencil( 1,-2, 1, center=3),
160 ),
161 (-1,3,-3,1) => (
162 Stencil(-1,3,-3,1, center=4),
163 ),
164 (1, -4, 6, -4, 1) => (
165 Stencil(1, -4, 6, -4, 1, center=4),
166 Stencil(1, -4, 6, -4, 1, center=5),
167 )
168 )
169 @testset "interior_weights = $w" for (w, closure_stencils) ∈ cases
170 @test dissipation_upper_closure_stencils(w) == closure_stencils
171 end
172 end
173
174
175 @testset "dissipation_transpose_lower_closure_stencils" begin
176 cases = (
177 (-1,1) => (
178 Stencil(-1,-1, 0, center=1),
179 Stencil( 1, 1,-1, center=2),
180 ),
181 (1,-2,1) => (
182 Stencil( 1, 1, 0, 0, center=1),
183 Stencil(-2,-2, 1, 0, center=2),
184 Stencil( 1, 1,-2, 1, center=3),
185 ),
186 (-1,3,-3,1) => (
187 Stencil(-1,-1,-1, 0, 0, 0, center=1),
188 Stencil( 3, 3, 3,-1, 0, 0, center=2),
189 Stencil(-3,-3,-3, 3,-1, 0, center=3),
190 Stencil( 1, 1, 1,-3, 3,-1, center=4),
191 ),
192 )
193 @testset "interior_weights = $w" for (w, closure_stencils) ∈ cases
194 @test dissipation_transpose_lower_closure_stencils(w) == closure_stencils
195 end
196 end
197
198 @testset "dissipation_transpose_upper_closure_stencils" begin
199 cases = (
200 (-1,1) => (
201 Stencil( 1,-1, center = 1),
202 Stencil( 0, 1, center = 2),
203 ),
204 (1,-2,1) => (
205 Stencil( 1,-2, 1, 1, center=2),
206 Stencil( 0, 1,-2,-2, center=3),
207 Stencil( 0, 0, 1, 1, center=4),
208 ),
209 (-1,3,-3,1) => (
210 Stencil( 1,-3, 3,-1,-1, center=2),
211 Stencil( 0, 1,-3, 3, 3, center=3),
212 Stencil( 0, 0, 1,-3,-3, center=4),
213 Stencil( 0, 0, 0, 1, 1, center=5),
214 ),
215 )
216 @testset "interior_weights = $w" for (w, closure_stencils) ∈ cases
217 @test dissipation_transpose_upper_closure_stencils(w) == closure_stencils
218 end
219 end