Mercurial > repos > public > sbplib_julia

comparison test/testSbpOperators.jl @ 504:21fba50cb5b0 feature/quadrature_as_outer_product

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Use LazyOuterProduct to construct multi-dimensional quadratures. This change allwed to: - Replace the types Quadrature and InverseQuadrature by functions returning outer products of the 1D operators. - Avoid convoluted naming of types. DiagonalInnerProduct is now renamed to DiagonalQuadrature, similarly for InverseDiagonalInnerProduct.
author Vidar Stiernström <vidar.stiernstrom@it.uu.se>
date Sat, 07 Nov 2020 13:07:31 +0100
parents 3cecbfb3d623
children 26485066394a
comparison
equal deleted inserted replaced
503:fbbb3733650c 504:21fba50cb5b0
108 l2(v) = sqrt(prod(h)*sum(v.^2)) 108 l2(v) = sqrt(prod(h)*sum(v.^2))
109 @test L*v4 ≈ v2 atol=5e-4 norm=l2 109 @test L*v4 ≈ v2 atol=5e-4 norm=l2
110 @test L*v5 ≈ v5ₓₓ atol=5e-4 norm=l2 110 @test L*v5 ≈ v5ₓₓ atol=5e-4 norm=l2
111 end 111 end
112 112
113 @testset "DiagonalInnerProduct" begin 113 @testset "DiagonalQuadrature" begin
114 op = readOperator(sbp_operators_path()*"d2_4th.txt",sbp_operators_path()*"h_4th.txt") 114 op = readOperator(sbp_operators_path()*"d2_4th.txt",sbp_operators_path()*"h_4th.txt")
115 L = 2.3 115 Lx = 2.3
116 g = EquidistantGrid(77, 0.0, L) 116 g = EquidistantGrid(77, 0.0, Lx)
117 H = DiagonalInnerProduct(g,op.quadratureClosure) 117 H = DiagonalQuadrature(g,op.quadratureClosure)
118 v = ones(Float64, size(g)) 118 v = ones(Float64, size(g))
119 119
120 @test H isa TensorMapping{T,1,1} where T 120 @test H isa TensorMapping{T,1,1} where T
121 @test H' isa TensorMapping{T,1,1} where T 121 @test H' isa TensorMapping{T,1,1} where T
122 @test sum(H*v) ≈ L 122 @test H == diagonal_quadrature(g,op.quadratureClosure)
123 @test sum(H*v) ≈ Lx
123 @test H*v == H'*v 124 @test H*v == H'*v
124 end 125 @inferred H*v
125 126
126 @testset "Quadrature" begin 127 # Test multiple dimension
128 Ly = 5.2
129 g = EquidistantGrid((77,66), (0.0, 0.0), (Lx,Ly))
130
131 H = diagonal_quadrature(g, op.quadratureClosure)
132
133 @test H isa TensorMapping{T,2,2} where T
134 @test H' isa TensorMapping{T,2,2} where T
135
136 v = ones(Float64, size(g))
137 @test sum(H*v) ≈ Lx*Ly
138
139 v = 2*ones(Float64, size(g))
140 @test sum(H*v) ≈ 2*Lx*Ly
141
142 @test_broken H*v == H'*v
143
144 @inferred H*v
145 end
146
147 @testset "InverseDiagonalQuadrature" begin
127 op = readOperator(sbp_operators_path()*"d2_4th.txt",sbp_operators_path()*"h_4th.txt") 148 op = readOperator(sbp_operators_path()*"d2_4th.txt",sbp_operators_path()*"h_4th.txt")
128 Lx = 2.3 149 Lx = 2.3
150 g = EquidistantGrid(77, 0.0, Lx)
151 H = DiagonalQuadrature(g, op.quadratureClosure)
152 Hi = InverseDiagonalQuadrature(g,op.quadratureClosure)
153 v = evalOn(g, x->sin(x))
154
155 @test Hi isa TensorMapping{T,1,1} where T
156 @test Hi' isa TensorMapping{T,1,1} where T
157 @test Hi == inverse_diagonal_quadrature(g,op.quadratureClosure)
158 @test Hi*H*v ≈ v
159 @test Hi*v == Hi'*v
160 @inferred Hi*v
161
162 # Test multiple dimension
129 Ly = 5.2 163 Ly = 5.2
130 g = EquidistantGrid((77,66), (0.0, 0.0), (Lx,Ly)) 164 g = EquidistantGrid((77,66), (0.0, 0.0), (Lx,Ly))
131 165
132 Q = Quadrature(g, op.quadratureClosure) 166 H = diagonal_quadrature(g, op.quadratureClosure)
133 167 Hi = inverse_diagonal_quadrature(g, op.quadratureClosure)
134 @test Q isa TensorMapping{T,2,2} where T 168 v = evalOn(g, (x,y)-> x^2 + (y-1)^2 + x*y)
135 @test Q' isa TensorMapping{T,2,2} where T 169
136 170 @test Hi isa TensorMapping{T,2,2} where T
137 v = ones(Float64, size(g)) 171 @test Hi' isa TensorMapping{T,2,2} where T
138 @test sum(Q*v) ≈ Lx*Ly
139
140 v = 2*ones(Float64, size(g))
141 @test_broken sum(Q*v) ≈ 2*Lx*Ly
142
143 @test Q*v == Q'*v
144 end
145
146 @testset "InverseDiagonalInnerProduct" begin
147 op = readOperator(sbp_operators_path()*"d2_4th.txt",sbp_operators_path()*"h_4th.txt")
148 L = 2.3
149 g = EquidistantGrid(77, 0.0, L)
150 H = DiagonalInnerProduct(g, op.quadratureClosure)
151 Hi = InverseDiagonalInnerProduct(g,op.quadratureClosure)
152 v = evalOn(g, x->sin(x))
153
154 @test Hi isa TensorMapping{T,1,1} where T
155 @test Hi' isa TensorMapping{T,1,1} where T
156 @test Hi*H*v ≈ v 172 @test Hi*H*v ≈ v
157 @test Hi*v == Hi'*v 173 @test_broken Hi*v == Hi'*v
158 end
159
160 @testset "InverseQuadrature" begin
161 op = readOperator(sbp_operators_path()*"d2_4th.txt",sbp_operators_path()*"h_4th.txt")
162 Lx = 7.3
163 Ly = 8.2
164 g = EquidistantGrid((77,66), (0.0, 0.0), (Lx,Ly))
165
166 Q = Quadrature(g, op.quadratureClosure)
167 Qinv = InverseQuadrature(g, op.quadratureClosure)
168 v = evalOn(g, (x,y)-> x^2 + (y-1)^2 + x*y)
169
170 @test Qinv isa TensorMapping{T,2,2} where T
171 @test Qinv' isa TensorMapping{T,2,2} where T
172 @test_broken Qinv*(Q*v) ≈ v
173 @test Qinv*v == Qinv'*v
174 end 174 end
175 # 175 #
176 # @testset "BoundaryValue" begin 176 # @testset "BoundaryValue" 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 # g = EquidistantGrid((4,5), (0.0, 0.0), (1.0,1.0)) 178 # g = EquidistantGrid((4,5), (0.0, 0.0), (1.0,1.0))