Mercurial > repos > public > sbplib_julia

comparison test/testSbpOperators.jl @ 565:15423a868d28 feature/boundary_ops

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Restructure and extend tests for BoundaryRestriction
author Vidar Stiernström <vidar.stiernstrom@it.uu.se>
date Mon, 30 Nov 2020 23:19:20 +0100
parents 8f7919a9b398
children fe026b4f99ec
comparison
equal deleted inserted replaced
564:ccb41095def6 565:15423a868d28
173 @test Qinv*v == Qinv'*v 173 @test Qinv*v == Qinv'*v
174 end 174 end
175 175
176 @testset "BoundaryRestrictrion" begin 176 @testset "BoundaryRestrictrion" begin
177 op = readOperator(sbp_operators_path()*"d2_4th.txt",sbp_operators_path()*"h_4th.txt") 177 op = readOperator(sbp_operators_path()*"d2_4th.txt",sbp_operators_path()*"h_4th.txt")
178 178 g_1D = EquidistantGrid(4, 0.0, 1.0)
179 g = EquidistantGrid(4, 0.0, 1.0) 179 g_2D = EquidistantGrid((4,5), (0.0, 0.0), (1.0,1.0))
180 180
181 e_l = BoundaryRestriction(g,op.eClosure,Lower()) 181 @testset "Constructors" begin
182 e_r = BoundaryRestriction(g,op.eClosure,Upper()) 182 # 1D
183 183 e_l = BoundaryRestriction{Float64,4,Lower}(op.eClosure,size(g_1D))
184 v = evalOn(g,x->1+x^2) 184 @test e_l == BoundaryRestriction(g_1D,op.eClosure,Lower())
185 u = fill(3.124) 185 @test e_l == boundary_restriction(g_1D,op.eClosure,CartesianBoundary{1,Lower}())
186 186 @test e_l isa TensorMapping{T,0,1} where T
187 @test (e_l*v)[] == v[1] 187
188 @test (e_r*v)[] == v[end] 188 e_r = BoundaryRestriction{Float64,4,Upper}(op.eClosure,size(g_1D))
189 @test e_l'*u == [u[], 0, 0, 0] 189 @test e_r == BoundaryRestriction(g_1D,op.eClosure,Upper())
190 @test e_r'*u == [0, 0, 0, u[]] 190 @test e_r == boundary_restriction(g_1D,op.eClosure,CartesianBoundary{1,Upper}())
191 @test_throws BoundsError (e_l*v)[Index{Lower}(3)] 191 @test e_r isa TensorMapping{T,0,1} where T
192 @test_throws BoundsError (e_r*v)[Index{Upper}(3)] 192
193 193 # 2D
194 g = EquidistantGrid((4,5), (0.0, 0.0), (1.0,1.0)) 194 e_w = boundary_restriction(g_2D,op.eClosure,CartesianBoundary{1,Upper}())
195 195 @test e_w isa InflatedTensorMapping
196 e_w = boundary_restriction(g, op.eClosure, CartesianBoundary{1,Lower}()) 196 @test e_w isa TensorMapping{T,1,2} where T
197 e_e = boundary_restriction(g, op.eClosure, CartesianBoundary{1,Upper}())
198 e_s = boundary_restriction(g, op.eClosure, CartesianBoundary{2,Lower}())
199 e_n = boundary_restriction(g, op.eClosure, CartesianBoundary{2,Upper}())
200
201 v = zeros(Float64, 4, 5)
202 v[:,5] = [1, 2, 3,4]
203 v[:,4] = [1, 2, 3,4]
204 v[:,3] = [4, 5, 6, 7]
205 v[:,2] = [7, 8, 9, 10]
206 v[:,1] = [10, 11, 12, 13]
207
208 @test e_w isa TensorMapping{T,1,2} where T
209 @test e_w' isa TensorMapping{T,2,1} where T
210
211 @test domain_size(e_w) == (4,5)
212 @test domain_size(e_e) == (4,5)
213 @test domain_size(e_s) == (4,5)
214 @test domain_size(e_n) == (4,5)
215
216 @test range_size(e_w) == (5,)
217 @test range_size(e_e) == (5,)
218 @test range_size(e_s) == (4,)
219 @test range_size(e_n) == (4,)
220
221 I_w = [(Index{Lower}(1),),
222 (Index{Interior}(2),),
223 (Index{Interior}(3),),
224 (Index{Interior}(4),),
225 (Index{Upper}(5),)]
226 v_w = [10,7,4,1.0,1];
227 for i = 1:length(I_w)
228 @test (e_w*v)[I_w[i]...] == v_w[i];
229 end 197 end
230 @test e_w*v == [10,7,4,1.0,1] 198
231 @test e_e*v == [13,10,7,4,4.0] 199 e_l = boundary_restriction(g_1D, op.eClosure, CartesianBoundary{1,Lower}())
232 @test e_s*v == [10,11,12,13.0] 200 e_r = boundary_restriction(g_1D, op.eClosure, CartesianBoundary{1,Upper}())
233 @test e_n*v == [1,2,3,4.0] 201
234 202 e_w = boundary_restriction(g_2D, op.eClosure, CartesianBoundary{1,Lower}())
235 g_x = [1,2,3,4.0] 203 e_e = boundary_restriction(g_2D, op.eClosure, CartesianBoundary{1,Upper}())
236 g_y = [5,4,3,2,1.0] 204 e_s = boundary_restriction(g_2D, op.eClosure, CartesianBoundary{2,Lower}())
237 205 e_n = boundary_restriction(g_2D, op.eClosure, CartesianBoundary{2,Upper}())
238 G_w = zeros(Float64, (4,5)) 206
239 G_w[1,:] = g_y 207 @testset "Sizes" begin
240 208 # 1D
241 G_e = zeros(Float64, (4,5)) 209 @test domain_size(e_l) == (4,)
242 G_e[4,:] = g_y 210 @test domain_size(e_r) == (4,)
243 211
244 G_s = zeros(Float64, (4,5)) 212 @test range_size(e_l) == ()
245 G_s[:,1] = g_x 213 @test range_size(e_r) == ()
246 214
247 G_n = zeros(Float64, (4,5)) 215 # 2D
248 G_n[:,5] = g_x 216 @test domain_size(e_w) == (4,5)
249 217 @test domain_size(e_e) == (4,5)
250 @test e_w'*g_y == G_w 218 @test domain_size(e_s) == (4,5)
251 @test e_e'*g_y == G_e 219 @test domain_size(e_n) == (4,5)
252 @test e_s'*g_x == G_s 220
253 @test e_n'*g_x == G_n 221 @test range_size(e_w) == (5,)
222 @test range_size(e_e) == (5,)
223 @test range_size(e_s) == (4,)
224 @test range_size(e_n) == (4,)
225 end
226
227
228 @testset "Application" begin
229 # 1D
230 v = evalOn(g_1D,x->1+x^2)
231 u = fill(3.124)
232 @test (e_l*v)[] == v[1]
233 @test (e_r*v)[] == v[end]
234 @test (e_r*v)[1] == v[end]
235 @test e_l'*u == [u[], 0, 0, 0]
236 @test e_r'*u == [0, 0, 0, u[]]
237 @test_throws BoundsError (e_l*v)[2]
238 @test_throws BoundsError (e_l'*u)[5]
239
240 # 2D
241 v = zeros(Float64, 4, 5)
242 v[:,5] = [1, 2, 3,4]
243 v[:,4] = [1, 2, 3,4]
244 v[:,3] = [4, 5, 6, 7]
245 v[:,2] = [7, 8, 9, 10]
246 v[:,1] = [10, 11, 12, 13]
247
248 @test e_w*v == [10,7,4,1.0,1]
249 @test e_e*v == [13,10,7,4,4.0]
250 @test e_s*v == [10,11,12,13.0]
251 @test e_n*v == [1,2,3,4.0]
252
253 I_w = [(Index{Lower}(1),),
254 (Index{Interior}(2),),
255 (Index{Interior}(3),),
256 (Index{Interior}(4),),
257 (Index{Upper}(5),)]
258 for i = 1:length(I_w)
259 @test (e_w*v)[I_w[i]...] == [10,7,4,1.0,1][i];
260 end
261
262 g_x = [1,2,3,4.0]
263 g_y = [5,4,3,2,1.0]
264
265 G_w = zeros(Float64, (4,5))
266 G_w[1,:] = g_y
267
268 G_e = zeros(Float64, (4,5))
269 G_e[4,:] = g_y
270
271 G_s = zeros(Float64, (4,5))
272 G_s[:,1] = g_x
273
274 G_n = zeros(Float64, (4,5))
275 G_n[:,5] = g_x
276
277 @test e_w'*g_y == G_w
278 @test e_e'*g_y == G_e
279 @test e_s'*g_x == G_s
280 @test e_n'*g_x == G_n
281 end
282
283 @testset "Inferred" begin
284 # 1D
285 v = ones(Float64, 4)
286 u = fill(1.)
287 @inferred (e_l*v)[] == 1
288 @inferred (e_r*v)[] == 1
289 @inferred e_l'*u == [1., 0., 0., 0.]
290 @inferred e_r'*u == [0., 0., 0., 1.]
291
292 # 2D
293 v = ones(Float64, 4, 5)
294 @inferred e_w*v == ones(Float64, 5)
295 @inferred e_e*v == ones(Float64, 5)
296 @inferred e_s*v == ones(Float64, 4)
297 @inferred e_n*v == ones(Float64, 4)
298
299 g_x = ones(Float64,4)
300 g_y = ones(Float64,5)
301
302 G_w = zeros(Float64, (4,5))
303 G_w[1,:] = g_y
304
305 G_e = zeros(Float64, (4,5))
306 G_e[4,:] = g_y
307
308 G_s = zeros(Float64, (4,5))
309 G_s[:,1] = g_x
310
311 G_n = zeros(Float64, (4,5))
312 G_n[:,5] = g_x
313
314 @inferred e_w'*g_y == G_w
315 @inferred e_e'*g_y == G_e
316 @inferred e_s'*g_x == G_s
317 @inferred e_n'*g_x == G_n
318 end
319
254 end 320 end
255 # 321 #
256 # @testset "NormalDerivative" begin 322 # @testset "NormalDerivative" begin
257 # op = readOperator(sbp_operators_path()*"d2_4th.txt",sbp_operators_path()*"h_4th.txt") 323 # op = readOperator(sbp_operators_path()*"d2_4th.txt",sbp_operators_path()*"h_4th.txt")
258 # g = EquidistantGrid((5,6), (0.0, 0.0), (4.0,5.0)) 324 # g = EquidistantGrid((5,6), (0.0, 0.0), (4.0,5.0))