comparison test/testSbpOperators.jl @ 629:22a0971f7f84 feature/volume_and_boundary_operators

Add tests for SecondDerivative
author Vidar Stiernström <vidar.stiernstrom@it.uu.se>
date Thu, 31 Dec 2020 08:26:00 +0100
parents 316dbfd31d35
children 9f0f1ace5101
comparison
equal deleted inserted replaced
628:316dbfd31d35 629:22a0971f7f84
256 Lx = 3.5 256 Lx = 3.5
257 Ly = 3. 257 Ly = 3.
258 g_1D = EquidistantGrid(121, 0.0, Lx) 258 g_1D = EquidistantGrid(121, 0.0, Lx)
259 g_2D = EquidistantGrid((121,123), (0.0, 0.0), (Lx, Ly)) 259 g_2D = EquidistantGrid((121,123), (0.0, 0.0), (Lx, Ly))
260 260
261 # TODO: These areant really constructors. Better name?
261 @testset "Constructors" begin 262 @testset "Constructors" begin
262 @testset "1D" begin 263 @testset "1D" begin
263 Dₓₓ = SecondDerivative(g_1D,op.innerStencil,op.closureStencils) 264 Dₓₓ = SecondDerivative(g_1D,op.innerStencil,op.closureStencils)
264 Dₓₓ == SecondDerivative(g_1D,op.innerStencil,op.closureStencils,1) 265 @test Dₓₓ == SecondDerivative(g_1D,op.innerStencil,op.closureStencils,1)
265 Dₓₓ isa VolumeOperator 266 @test Dₓₓ isa VolumeOperator
266 end 267 end
267 @testset "2D" begin 268 @testset "2D" begin
268 Dₓₓ = SecondDerivative(g_2D,op.innerStencil,op.closureStencils,1) 269 Dₓₓ = SecondDerivative(g_2D,op.innerStencil,op.closureStencils,1)
269 Dₓₓ isa TensorMappingComposition 270 D2 = SecondDerivative(g_1D,op.innerStencil,op.closureStencils)
270 Dₓₓ isa TensorMapping{2,2} 271 I = IdentityMapping{Float64}(size(g_2D)[2])
272 @test Dₓₓ == D2⊗I
273 @test Dₓₓ isa TensorMapping{T,2,2} where T
271 end 274 end
272 end 275 end
273 276
274 @testset "Accuracy" begin 277 @testset "Accuracy" begin
275 @testset "1D" begin 278 @testset "1D" begin
276 v0 = evalOn(g_1D,x->0.)
277 monomials = () 279 monomials = ()
278 maxOrder = 4; 280 maxOrder = 4;
279 for i = 0:maxOrder-1 281 for i = 0:maxOrder-1
280 f_i(x) = 1/factorial(i)*x^i 282 f_i(x) = 1/factorial(i)*x^i
281 monomials = (monomials...,evalOn(g_1D,f_i)) 283 monomials = (monomials...,evalOn(g_1D,f_i))
284 286
285 #TODO: Error when reading second order stencil! 287 #TODO: Error when reading second order stencil!
286 # @testset "2nd order" begin 288 # @testset "2nd order" begin
287 # op = read_D2_operator(sbp_operators_path()*"standard_diagonal.toml"; order=2) 289 # op = read_D2_operator(sbp_operators_path()*"standard_diagonal.toml"; order=2)
288 # Dₓₓ = SecondDerivative(g_1D,op.innerStencil,op.closureStencils) 290 # Dₓₓ = SecondDerivative(g_1D,op.innerStencil,op.closureStencils)
289 # @test Dₓₓ*monomials[1] ≈ v0 atol = 5e-13 291 # @test Dₓₓ*monomials[1] ≈ zeros(Float64,size(g_1D)...) atol = 5e-13
290 # @test Dₓₓ*monomials[2] ≈ v0 atol = 5e-13 292 # @test Dₓₓ*monomials[2] ≈ zeros(Float64,size(g_1D)...) atol = 5e-13
291 # end 293 # end
292 294
293 # 4th order interior stencil, 2nd order boundary stencil, 295 # 4th order interior stencil, 2nd order boundary stencil,
294 # implies that L*v should be exact for monomials up to order 3. 296 # implies that L*v should be exact for monomials up to order 3.
295 # Exact differentiation is measured point-wise. For other grid functions 297 # Exact differentiation is measured point-wise. For other grid functions
297 @testset "4th order" begin 299 @testset "4th order" begin
298 op = read_D2_operator(sbp_operators_path()*"standard_diagonal.toml"; order=4) 300 op = read_D2_operator(sbp_operators_path()*"standard_diagonal.toml"; order=4)
299 Dₓₓ = SecondDerivative(g_1D,op.innerStencil,op.closureStencils) 301 Dₓₓ = SecondDerivative(g_1D,op.innerStencil,op.closureStencils)
300 # TODO: high tolerances for checking the "exact" differentiation 302 # TODO: high tolerances for checking the "exact" differentiation
301 # due to accumulation of round-off errors/cancellation errors? 303 # due to accumulation of round-off errors/cancellation errors?
302 @test Dₓₓ*monomials[1] ≈ v0 atol = 5e-10 304 @test Dₓₓ*monomials[1] ≈ zeros(Float64,size(g_1D)...) atol = 5e-10
303 @test Dₓₓ*monomials[2] ≈ v0 atol = 5e-10 305 @test Dₓₓ*monomials[2] ≈ zeros(Float64,size(g_1D)...) atol = 5e-10
304 @test Dₓₓ*monomials[3] ≈ monomials[1] atol = 5e-10 306 @test Dₓₓ*monomials[3] ≈ monomials[1] atol = 5e-10
305 @test Dₓₓ*monomials[4] ≈ monomials[2] atol = 5e-10 307 @test Dₓₓ*monomials[4] ≈ monomials[2] atol = 5e-10
306 @test Dₓₓ*evalOn(g_1D,x -> sin(x)) ≈ evalOn(g_1D,x -> -sin(x)) rtol = 5e-4 norm = l2 308 @test Dₓₓ*evalOn(g_1D,x -> sin(x)) ≈ evalOn(g_1D,x -> -sin(x)) rtol = 5e-4 norm = l2
307 end 309 end
308 end 310 end
331 @testset "4th order" begin 333 @testset "4th order" begin
332 op = read_D2_operator(sbp_operators_path()*"standard_diagonal.toml"; order=4) 334 op = read_D2_operator(sbp_operators_path()*"standard_diagonal.toml"; order=4)
333 Dyy = SecondDerivative(g_2D,op.innerStencil,op.closureStencils,2) 335 Dyy = SecondDerivative(g_2D,op.innerStencil,op.closureStencils,2)
334 # TODO: high tolerances for checking the "exact" differentiation 336 # TODO: high tolerances for checking the "exact" differentiation
335 # due to accumulation of round-off errors/cancellation errors? 337 # due to accumulation of round-off errors/cancellation errors?
336 @test Dyy*binomials[1] ≈ evalOn(g_2D,(x,y)->0.) atol = 5e-9 338 @test Dyy*binomials[1] ≈ zeros(Float64,size(g_2D)...) atol = 5e-9
337 @test Dyy*binomials[2] ≈ evalOn(g_2D,(x,y)->0.) atol = 5e-9 339 @test Dyy*binomials[2] ≈ zeros(Float64,size(g_2D)...) atol = 5e-9
338 @test Dyy*binomials[3] ≈ evalOn(g_2D,(x,y)->1.) atol = 5e-9 340 @test Dyy*binomials[3] ≈ evalOn(g_2D,(x,y)->1.) atol = 5e-9
339 @test Dyy*binomials[4] ≈ evalOn(g_2D,(x,y)->y) atol = 5e-9 341 @test Dyy*binomials[4] ≈ evalOn(g_2D,(x,y)->y) atol = 5e-9
340 @test Dyy*evalOn(g_2D, (x,y) -> sin(x)+cos(y)) ≈ evalOn(g_2D,(x,y) -> -cos(y)) rtol = 5e-4 norm = l2 342 @test Dyy*evalOn(g_2D, (x,y) -> sin(x)+cos(y)) ≈ evalOn(g_2D,(x,y) -> -cos(y)) rtol = 5e-4 norm = l2
341 end 343 end
342 end 344 end