Mercurial > repos > public > sbplib_julia

comparison test/SbpOperators/volumeops/derivatives/second_derivative_test.jl @ 980:f885e1de6dc4 feature/variable_derivatives

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Merge default
author Jonatan Werpers <jonatan@werpers.com>
date Tue, 15 Mar 2022 21:38:55 +0100
parents 2ae62dbaf839
children 7bf3121c6864 1ba8a398af9c
comparison
equal deleted inserted replaced
971:bc12be1b1ae5 980:f885e1de6dc4
5 using Sbplib.LazyTensors 5 using Sbplib.LazyTensors
6 6
7 import Sbplib.SbpOperators.VolumeOperator 7 import Sbplib.SbpOperators.VolumeOperator
8 8
9 @testset "SecondDerivative" begin 9 @testset "SecondDerivative" begin
10 stencil_set = read_stencil_set(sbp_operators_path()*"standard_diagonal.toml"; order=4) 10 operator_path = sbp_operators_path()*"standard_diagonal.toml"
11 stencil_set = read_stencil_set(operator_path; order=4)
11 inner_stencil = parse_stencil(stencil_set["D2"]["inner_stencil"]) 12 inner_stencil = parse_stencil(stencil_set["D2"]["inner_stencil"])
12 closure_stencils = parse_stencil.(stencil_set["D2"]["closure_stencils"]) 13 closure_stencils = parse_stencil.(stencil_set["D2"]["closure_stencils"])
13 Lx = 3.5 14 Lx = 3.5
14 Ly = 3. 15 Ly = 3.
15 g_1D = EquidistantGrid(121, 0.0, Lx) 16 g_1D = EquidistantGrid(121, 0.0, Lx)
16 g_2D = EquidistantGrid((121,123), (0.0, 0.0), (Lx, Ly)) 17 g_2D = EquidistantGrid((121,123), (0.0, 0.0), (Lx, Ly))
17 18
18 @testset "Constructors" begin 19 @testset "Constructors" begin
19 @testset "1D" begin 20 @testset "1D" begin
20 Dₓₓ = second_derivative(g_1D,inner_stencil,closure_stencils) 21 Dₓₓ = second_derivative(g_1D,inner_stencil,closure_stencils,1)
21 @test Dₓₓ == second_derivative(g_1D,inner_stencil,closure_stencils,1) 22 @test Dₓₓ == second_derivative(g_1D,inner_stencil,closure_stencils)
23 @test Dₓₓ == second_derivative(g_1D,stencil_set,1)
22 @test Dₓₓ isa VolumeOperator 24 @test Dₓₓ isa VolumeOperator
23 end 25 end
24 @testset "2D" begin 26 @testset "2D" begin
25 Dₓₓ = second_derivative(g_2D,inner_stencil,closure_stencils,1) 27 Dₓₓ = second_derivative(g_2D,inner_stencil,closure_stencils,1)
26 D2 = second_derivative(g_1D,inner_stencil,closure_stencils) 28 D2 = second_derivative(g_1D,inner_stencil,closure_stencils)
27 I = IdentityMapping{Float64}(size(g_2D)[2]) 29 I = IdentityMapping{Float64}(size(g_2D)[2])
28 @test Dₓₓ == D2⊗I 30 @test Dₓₓ == D2⊗I
31 @test Dₓₓ == second_derivative(g_2D,stencil_set,1)
29 @test Dₓₓ isa TensorMapping{T,2,2} where T 32 @test Dₓₓ isa TensorMapping{T,2,2} where T
30 end 33 end
31 end 34 end
32 35
33 # Exact differentiation is measured point-wise. In other cases 36 # Exact differentiation is measured point-wise. In other cases
45 vₓₓ = evalOn(g_1D,x -> -sin(x)) 48 vₓₓ = evalOn(g_1D,x -> -sin(x))
46 49
47 # 2nd order interior stencil, 1nd order boundary stencil, 50 # 2nd order interior stencil, 1nd order boundary stencil,
48 # implies that L*v should be exact for monomials up to order 2. 51 # implies that L*v should be exact for monomials up to order 2.
49 @testset "2nd order" begin 52 @testset "2nd order" begin
50 stencil_set = read_stencil_set(sbp_operators_path()*"standard_diagonal.toml"; order=2) 53 stencil_set = read_stencil_set(operator_path; order=2)
51 inner_stencil = parse_stencil(stencil_set["D2"]["inner_stencil"]) 54 inner_stencil = parse_stencil(stencil_set["D2"]["inner_stencil"])
52 closure_stencils = parse_stencil.(stencil_set["D2"]["closure_stencils"]) 55 closure_stencils = parse_stencil.(stencil_set["D2"]["closure_stencils"])
53 Dₓₓ = second_derivative(g_1D,inner_stencil,closure_stencils) 56 Dₓₓ = second_derivative(g_1D,inner_stencil,closure_stencils)
54 @test Dₓₓ*monomials[1] ≈ zeros(Float64,size(g_1D)...) atol = 5e-10 57 @test Dₓₓ*monomials[1] ≈ zeros(Float64,size(g_1D)...) atol = 5e-10
55 @test Dₓₓ*monomials[2] ≈ zeros(Float64,size(g_1D)...) atol = 5e-10 58 @test Dₓₓ*monomials[2] ≈ zeros(Float64,size(g_1D)...) atol = 5e-10
58 end 61 end
59 62
60 # 4th order interior stencil, 2nd order boundary stencil, 63 # 4th order interior stencil, 2nd order boundary stencil,
61 # implies that L*v should be exact for monomials up to order 3. 64 # implies that L*v should be exact for monomials up to order 3.
62 @testset "4th order" begin 65 @testset "4th order" begin
63 stencil_set = read_stencil_set(sbp_operators_path()*"standard_diagonal.toml"; order=4) 66 stencil_set = read_stencil_set(operator_path; order=4)
64 inner_stencil = parse_stencil(stencil_set["D2"]["inner_stencil"]) 67 inner_stencil = parse_stencil(stencil_set["D2"]["inner_stencil"])
65 closure_stencils = parse_stencil.(stencil_set["D2"]["closure_stencils"]) 68 closure_stencils = parse_stencil.(stencil_set["D2"]["closure_stencils"])
66 Dₓₓ = second_derivative(g_1D,inner_stencil,closure_stencils) 69 Dₓₓ = second_derivative(g_1D,inner_stencil,closure_stencils)
67 # NOTE: high tolerances for checking the "exact" differentiation 70 # NOTE: high tolerances for checking the "exact" differentiation
68 # due to accumulation of round-off errors/cancellation errors? 71 # due to accumulation of round-off errors/cancellation errors?
86 v_yy = evalOn(g_2D,(x,y) -> -cos(y)) 89 v_yy = evalOn(g_2D,(x,y) -> -cos(y))
87 90
88 # 2nd order interior stencil, 1st order boundary stencil, 91 # 2nd order interior stencil, 1st order boundary stencil,
89 # implies that L*v should be exact for binomials up to order 2. 92 # implies that L*v should be exact for binomials up to order 2.
90 @testset "2nd order" begin 93 @testset "2nd order" begin
91 stencil_set = read_stencil_set(sbp_operators_path()*"standard_diagonal.toml"; order=2) 94 stencil_set = read_stencil_set(operator_path; order=2)
92 inner_stencil = parse_stencil(stencil_set["D2"]["inner_stencil"]) 95 inner_stencil = parse_stencil(stencil_set["D2"]["inner_stencil"])
93 closure_stencils = parse_stencil.(stencil_set["D2"]["closure_stencils"]) 96 closure_stencils = parse_stencil.(stencil_set["D2"]["closure_stencils"])
94 Dyy = second_derivative(g_2D,inner_stencil,closure_stencils,2) 97 Dyy = second_derivative(g_2D,inner_stencil,closure_stencils,2)
95 @test Dyy*binomials[1] ≈ zeros(Float64,size(g_2D)...) atol = 5e-9 98 @test Dyy*binomials[1] ≈ zeros(Float64,size(g_2D)...) atol = 5e-9
96 @test Dyy*binomials[2] ≈ zeros(Float64,size(g_2D)...) atol = 5e-9 99 @test Dyy*binomials[2] ≈ zeros(Float64,size(g_2D)...) atol = 5e-9
99 end 102 end
100 103
101 # 4th order interior stencil, 2nd order boundary stencil, 104 # 4th order interior stencil, 2nd order boundary stencil,
102 # implies that L*v should be exact for binomials up to order 3. 105 # implies that L*v should be exact for binomials up to order 3.
103 @testset "4th order" begin 106 @testset "4th order" begin
104 stencil_set = read_stencil_set(sbp_operators_path()*"standard_diagonal.toml"; order=4) 107 stencil_set = read_stencil_set(operator_path; order=4)
105 inner_stencil = parse_stencil(stencil_set["D2"]["inner_stencil"]) 108 inner_stencil = parse_stencil(stencil_set["D2"]["inner_stencil"])
106 closure_stencils = parse_stencil.(stencil_set["D2"]["closure_stencils"]) 109 closure_stencils = parse_stencil.(stencil_set["D2"]["closure_stencils"])
107 Dyy = second_derivative(g_2D,inner_stencil,closure_stencils,2) 110 Dyy = second_derivative(g_2D,inner_stencil,closure_stencils,2)
108 # NOTE: high tolerances for checking the "exact" differentiation 111 # NOTE: high tolerances for checking the "exact" differentiation
109 # due to accumulation of round-off errors/cancellation errors? 112 # due to accumulation of round-off errors/cancellation errors?