comparison test/SbpOperators/boundaryops/boundary_operator_test.jl @ 1099:05a25a5063bb refactor/sbpoperators/inflation

Try to remove volume_operator and boundary_operator methods
author Jonatan Werpers <jonatan@werpers.com>
date Mon, 21 Mar 2022 12:51:39 +0100
parents 1ba8a398af9c
children 157a78959e5d
comparison
equal deleted inserted replaced
1024:5be17f647018 1099:05a25a5063bb
4 using Sbplib.SbpOperators 4 using Sbplib.SbpOperators
5 using Sbplib.Grids 5 using Sbplib.Grids
6 using Sbplib.RegionIndices 6 using Sbplib.RegionIndices
7 import Sbplib.SbpOperators.Stencil 7 import Sbplib.SbpOperators.Stencil
8 import Sbplib.SbpOperators.BoundaryOperator 8 import Sbplib.SbpOperators.BoundaryOperator
9 import Sbplib.SbpOperators.boundary_operator 9
10 # TODO: What should happen to all the commented tests? Deleted? Replicated for user code?
10 11
11 @testset "BoundaryOperator" begin 12 @testset "BoundaryOperator" begin
12 closure_stencil = Stencil((0,2), (2.,1.,3.)) 13 closure_stencil = Stencil((0,2), (2.,1.,3.))
13 g_1D = EquidistantGrid(11, 0.0, 1.0) 14 g_1D = EquidistantGrid(11, 0.0, 1.0)
14 g_2D = EquidistantGrid((11,15), (0.0, 0.0), (1.0,1.0)) 15 g_2D = EquidistantGrid((11,15), (0.0, 0.0), (1.0,1.0))
15 16
16 @testset "Constructors" begin 17 @testset "Constructors" begin
17 @testset "1D" begin 18 @testset "1D" begin # TODO: Remove these testsets
18 op_l = BoundaryOperator{Lower}(closure_stencil,size(g_1D)[1]) 19 op_l = BoundaryOperator{Lower}(closure_stencil,size(g_1D)[1])
19 @test op_l == BoundaryOperator(g_1D,closure_stencil,Lower()) 20 @test op_l == BoundaryOperator(g_1D,closure_stencil,Lower())
20 @test op_l == boundary_operator(g_1D,closure_stencil,CartesianBoundary{1,Lower}())
21 @test op_l isa LazyTensor{T,0,1} where T 21 @test op_l isa LazyTensor{T,0,1} where T
22 22
23 op_r = BoundaryOperator{Upper}(closure_stencil,size(g_1D)[1]) 23 op_r = BoundaryOperator{Upper}(closure_stencil,size(g_1D)[1]) # TBD: Is this constructor really needed? looks weird!
24 @test op_r == BoundaryOperator(g_1D,closure_stencil,Upper()) 24 @test op_r == BoundaryOperator(g_1D,closure_stencil,Upper())
25 @test op_r == boundary_operator(g_1D,closure_stencil,CartesianBoundary{1,Upper}())
26 @test op_r isa LazyTensor{T,0,1} where T 25 @test op_r isa LazyTensor{T,0,1} where T
27 end 26 end
28 27
29 @testset "2D" begin 28 # @testset "2D" begin
30 e_w = boundary_operator(g_2D,closure_stencil,CartesianBoundary{1,Upper}()) 29 # e_w = boundary_operator(g_2D,closure_stencil,CartesianBoundary{1,Upper}())
31 @test e_w isa InflatedLazyTensor 30 # @test e_w isa InflatedLazyTensor
32 @test e_w isa LazyTensor{T,1,2} where T 31 # @test e_w isa LazyTensor{T,1,2} where T
33 end 32 # end
34 end 33 end
35 op_l, op_r = boundary_operator.(Ref(g_1D), Ref(closure_stencil), boundary_identifiers(g_1D)) 34
36 op_w, op_e, op_s, op_n = boundary_operator.(Ref(g_2D), Ref(closure_stencil), boundary_identifiers(g_2D)) 35 op_l = BoundaryOperator(g_1D, closure_stencil, Lower())
36 op_r = BoundaryOperator(g_1D, closure_stencil, Upper())
37 # op_w, op_e, op_s, op_n = boundary_operator.(Ref(g_2D), Ref(closure_stencil), boundary_identifiers(g_2D))
37 38
38 @testset "Sizes" begin 39 @testset "Sizes" begin
39 @testset "1D" begin 40 @testset "1D" begin
40 @test domain_size(op_l) == (11,) 41 @test domain_size(op_l) == (11,)
41 @test domain_size(op_r) == (11,) 42 @test domain_size(op_r) == (11,)
42 43
43 @test range_size(op_l) == () 44 @test range_size(op_l) == ()
44 @test range_size(op_r) == () 45 @test range_size(op_r) == ()
45 end 46 end
46 47
47 @testset "2D" begin 48 # @testset "2D" begin
48 @test domain_size(op_w) == (11,15) 49 # @test domain_size(op_w) == (11,15)
49 @test domain_size(op_e) == (11,15) 50 # @test domain_size(op_e) == (11,15)
50 @test domain_size(op_s) == (11,15) 51 # @test domain_size(op_s) == (11,15)
51 @test domain_size(op_n) == (11,15) 52 # @test domain_size(op_n) == (11,15)
52 53
53 @test range_size(op_w) == (15,) 54 # @test range_size(op_w) == (15,)
54 @test range_size(op_e) == (15,) 55 # @test range_size(op_e) == (15,)
55 @test range_size(op_s) == (11,) 56 # @test range_size(op_s) == (11,)
56 @test range_size(op_n) == (11,) 57 # @test range_size(op_n) == (11,)
57 end 58 # end
58 end 59 end
59 60
60 @testset "Application" begin 61 @testset "Application" begin
61 @testset "1D" begin 62 @testset "1D" begin
62 v = evalOn(g_1D,x->1+x^2) 63 v = evalOn(g_1D,x->1+x^2)
74 @test (op_l'*u)[1] isa ComplexF64 75 @test (op_l'*u)[1] isa ComplexF64
75 @test (op_l'*u)[5] isa ComplexF64 76 @test (op_l'*u)[5] isa ComplexF64
76 @test (op_l'*u)[11] isa ComplexF64 77 @test (op_l'*u)[11] isa ComplexF64
77 end 78 end
78 79
79 @testset "2D" begin 80 # @testset "2D" begin
80 v = rand(size(g_2D)...) 81 # v = rand(size(g_2D)...)
81 u = fill(3.124) 82 # u = fill(3.124)
82 @test op_w*v ≈ 2*v[1,:] + v[2,:] + 3*v[3,:] rtol = 1e-14 83 # @test op_w*v ≈ 2*v[1,:] + v[2,:] + 3*v[3,:] rtol = 1e-14
83 @test op_e*v ≈ 2*v[end,:] + v[end-1,:] + 3*v[end-2,:] rtol = 1e-14 84 # @test op_e*v ≈ 2*v[end,:] + v[end-1,:] + 3*v[end-2,:] rtol = 1e-14
84 @test op_s*v ≈ 2*v[:,1] + v[:,2] + 3*v[:,3] rtol = 1e-14 85 # @test op_s*v ≈ 2*v[:,1] + v[:,2] + 3*v[:,3] rtol = 1e-14
85 @test op_n*v ≈ 2*v[:,end] + v[:,end-1] + 3*v[:,end-2] rtol = 1e-14 86 # @test op_n*v ≈ 2*v[:,end] + v[:,end-1] + 3*v[:,end-2] rtol = 1e-14
86 87
87 88
88 g_x = rand(size(g_2D)[1]) 89 # g_x = rand(size(g_2D)[1])
89 g_y = rand(size(g_2D)[2]) 90 # g_y = rand(size(g_2D)[2])
90 91
91 G_w = zeros(Float64, size(g_2D)...) 92 # G_w = zeros(Float64, size(g_2D)...)
92 G_w[1,:] = 2*g_y 93 # G_w[1,:] = 2*g_y
93 G_w[2,:] = g_y 94 # G_w[2,:] = g_y
94 G_w[3,:] = 3*g_y 95 # G_w[3,:] = 3*g_y
95 96
96 G_e = zeros(Float64, size(g_2D)...) 97 # G_e = zeros(Float64, size(g_2D)...)
97 G_e[end,:] = 2*g_y 98 # G_e[end,:] = 2*g_y
98 G_e[end-1,:] = g_y 99 # G_e[end-1,:] = g_y
99 G_e[end-2,:] = 3*g_y 100 # G_e[end-2,:] = 3*g_y
100 101
101 G_s = zeros(Float64, size(g_2D)...) 102 # G_s = zeros(Float64, size(g_2D)...)
102 G_s[:,1] = 2*g_x 103 # G_s[:,1] = 2*g_x
103 G_s[:,2] = g_x 104 # G_s[:,2] = g_x
104 G_s[:,3] = 3*g_x 105 # G_s[:,3] = 3*g_x
105 106
106 G_n = zeros(Float64, size(g_2D)...) 107 # G_n = zeros(Float64, size(g_2D)...)
107 G_n[:,end] = 2*g_x 108 # G_n[:,end] = 2*g_x
108 G_n[:,end-1] = g_x 109 # G_n[:,end-1] = g_x
109 G_n[:,end-2] = 3*g_x 110 # G_n[:,end-2] = 3*g_x
110 111
111 @test op_w'*g_y == G_w 112 # @test op_w'*g_y == G_w
112 @test op_e'*g_y == G_e 113 # @test op_e'*g_y == G_e
113 @test op_s'*g_x == G_s 114 # @test op_s'*g_x == G_s
114 @test op_n'*g_x == G_n 115 # @test op_n'*g_x == G_n
115 end 116 # end
116 117
117 @testset "Regions" begin 118 @testset "Regions" begin
118 u = fill(3.124) 119 u = fill(3.124)
119 @test (op_l'*u)[Index(1,Lower)] == 2*u[] 120 @test (op_l'*u)[Index(1,Lower)] == 2*u[]
120 @test (op_l'*u)[Index(2,Lower)] == u[] 121 @test (op_l'*u)[Index(2,Lower)] == u[]