Mercurial > repos > public > sbplib_julia
diff test/SbpOperators/volumeops/derivatives/second_derivative_test.jl @ 1858:4a9be96f2569 feature/documenter_logo
Merge default
author | Jonatan Werpers <jonatan@werpers.com> |
---|---|
date | Sun, 12 Jan 2025 21:18:44 +0100 |
parents | 471a948cd2b2 |
children |
line wrap: on
line diff
--- a/test/SbpOperators/volumeops/derivatives/second_derivative_test.jl Fri Jan 21 15:23:08 2022 +0100 +++ b/test/SbpOperators/volumeops/derivatives/second_derivative_test.jl Sun Jan 12 21:18:44 2025 +0100 @@ -1,32 +1,32 @@ using Test -using Sbplib.SbpOperators -using Sbplib.Grids -using Sbplib.LazyTensors +using Diffinitive.SbpOperators +using Diffinitive.Grids +using Diffinitive.LazyTensors -import Sbplib.SbpOperators.VolumeOperator +import Diffinitive.SbpOperators.VolumeOperator + +# TODO: Refactor these test to look more like the tests in first_derivative_test.jl. @testset "SecondDerivative" begin - stencil_set = read_stencil_set(sbp_operators_path()*"standard_diagonal.toml"; order=4) + operator_path = sbp_operators_path()*"standard_diagonal.toml" + stencil_set = read_stencil_set(operator_path; order=4) inner_stencil = parse_stencil(stencil_set["D2"]["inner_stencil"]) closure_stencils = parse_stencil.(stencil_set["D2"]["closure_stencils"]) Lx = 3.5 Ly = 3. - g_1D = EquidistantGrid(121, 0.0, Lx) - g_2D = EquidistantGrid((121,123), (0.0, 0.0), (Lx, Ly)) + g_1D = equidistant_grid(0.0, Lx, 121) + g_2D = equidistant_grid((0.0, 0.0), (Lx, Ly), 121, 123) @testset "Constructors" begin @testset "1D" begin - Dₓₓ = second_derivative(g_1D,inner_stencil,closure_stencils) - @test Dₓₓ == second_derivative(g_1D,inner_stencil,closure_stencils,1) - @test Dₓₓ isa VolumeOperator + Dₓₓ = second_derivative(g_1D, stencil_set) + @test Dₓₓ == second_derivative(g_1D, inner_stencil, closure_stencils) + @test Dₓₓ isa LazyTensor{Float64,1,1} end @testset "2D" begin - Dₓₓ = second_derivative(g_2D,inner_stencil,closure_stencils,1) - D2 = second_derivative(g_1D,inner_stencil,closure_stencils) - I = IdentityMapping{Float64}(size(g_2D)[2]) - @test Dₓₓ == D2⊗I - @test Dₓₓ isa TensorMapping{T,2,2} where T + Dₓₓ = second_derivative(g_2D,stencil_set,1) + @test Dₓₓ isa LazyTensor{Float64,2,2} end end @@ -39,18 +39,16 @@ maxOrder = 4; for i = 0:maxOrder-1 f_i(x) = 1/factorial(i)*x^i - monomials = (monomials...,evalOn(g_1D,f_i)) + monomials = (monomials...,eval_on(g_1D,f_i)) end - v = evalOn(g_1D,x -> sin(x)) - vₓₓ = evalOn(g_1D,x -> -sin(x)) + v = eval_on(g_1D,x -> sin(x)) + vₓₓ = eval_on(g_1D,x -> -sin(x)) # 2nd order interior stencil, 1nd order boundary stencil, # implies that L*v should be exact for monomials up to order 2. @testset "2nd order" begin - stencil_set = read_stencil_set(sbp_operators_path()*"standard_diagonal.toml"; order=2) - inner_stencil = parse_stencil(stencil_set["D2"]["inner_stencil"]) - closure_stencils = parse_stencil.(stencil_set["D2"]["closure_stencils"]) - Dₓₓ = second_derivative(g_1D,inner_stencil,closure_stencils) + stencil_set = read_stencil_set(operator_path; order=2) + Dₓₓ = second_derivative(g_1D,stencil_set) @test Dₓₓ*monomials[1] ≈ zeros(Float64,size(g_1D)...) atol = 5e-10 @test Dₓₓ*monomials[2] ≈ zeros(Float64,size(g_1D)...) atol = 5e-10 @test Dₓₓ*monomials[3] ≈ monomials[1] atol = 5e-10 @@ -60,10 +58,8 @@ # 4th order interior stencil, 2nd order boundary stencil, # implies that L*v should be exact for monomials up to order 3. @testset "4th order" begin - stencil_set = read_stencil_set(sbp_operators_path()*"standard_diagonal.toml"; order=4) - inner_stencil = parse_stencil(stencil_set["D2"]["inner_stencil"]) - closure_stencils = parse_stencil.(stencil_set["D2"]["closure_stencils"]) - Dₓₓ = second_derivative(g_1D,inner_stencil,closure_stencils) + stencil_set = read_stencil_set(operator_path; order=4) + Dₓₓ = second_derivative(g_1D,stencil_set) # NOTE: high tolerances for checking the "exact" differentiation # due to accumulation of round-off errors/cancellation errors? @test Dₓₓ*monomials[1] ≈ zeros(Float64,size(g_1D)...) atol = 5e-10 @@ -75,42 +71,38 @@ end @testset "2D" begin - l2(v) = sqrt(prod(spacing(g_2D))*sum(v.^2)); + l2(v) = sqrt(prod(spacing.(g_2D.grids))*sum(v.^2)); binomials = () maxOrder = 4; for i = 0:maxOrder-1 f_i(x,y) = 1/factorial(i)*y^i + x^i - binomials = (binomials...,evalOn(g_2D,f_i)) + binomials = (binomials...,eval_on(g_2D,f_i)) end - v = evalOn(g_2D, (x,y) -> sin(x)+cos(y)) - v_yy = evalOn(g_2D,(x,y) -> -cos(y)) + v = eval_on(g_2D, (x,y) -> sin(x)+cos(y)) + v_yy = eval_on(g_2D,(x,y) -> -cos(y)) # 2nd order interior stencil, 1st order boundary stencil, # implies that L*v should be exact for binomials up to order 2. @testset "2nd order" begin - stencil_set = read_stencil_set(sbp_operators_path()*"standard_diagonal.toml"; order=2) - inner_stencil = parse_stencil(stencil_set["D2"]["inner_stencil"]) - closure_stencils = parse_stencil.(stencil_set["D2"]["closure_stencils"]) - Dyy = second_derivative(g_2D,inner_stencil,closure_stencils,2) + stencil_set = read_stencil_set(operator_path; order=2) + Dyy = second_derivative(g_2D,stencil_set,2) @test Dyy*binomials[1] ≈ zeros(Float64,size(g_2D)...) atol = 5e-9 @test Dyy*binomials[2] ≈ zeros(Float64,size(g_2D)...) atol = 5e-9 - @test Dyy*binomials[3] ≈ evalOn(g_2D,(x,y)->1.) atol = 5e-9 + @test Dyy*binomials[3] ≈ eval_on(g_2D,(x,y)->1.) atol = 5e-9 @test Dyy*v ≈ v_yy rtol = 5e-2 norm = l2 end # 4th order interior stencil, 2nd order boundary stencil, # implies that L*v should be exact for binomials up to order 3. @testset "4th order" begin - stencil_set = read_stencil_set(sbp_operators_path()*"standard_diagonal.toml"; order=4) - inner_stencil = parse_stencil(stencil_set["D2"]["inner_stencil"]) - closure_stencils = parse_stencil.(stencil_set["D2"]["closure_stencils"]) - Dyy = second_derivative(g_2D,inner_stencil,closure_stencils,2) + stencil_set = read_stencil_set(operator_path; order=4) + Dyy = second_derivative(g_2D,stencil_set,2) # NOTE: high tolerances for checking the "exact" differentiation # due to accumulation of round-off errors/cancellation errors? @test Dyy*binomials[1] ≈ zeros(Float64,size(g_2D)...) atol = 5e-9 @test Dyy*binomials[2] ≈ zeros(Float64,size(g_2D)...) atol = 5e-9 - @test Dyy*binomials[3] ≈ evalOn(g_2D,(x,y)->1.) atol = 5e-9 - @test Dyy*binomials[4] ≈ evalOn(g_2D,(x,y)->y) atol = 5e-9 + @test Dyy*binomials[3] ≈ eval_on(g_2D,(x,y)->1.) atol = 5e-9 + @test Dyy*binomials[4] ≈ eval_on(g_2D,(x,y)->y) atol = 5e-9 @test Dyy*v ≈ v_yy rtol = 5e-4 norm = l2 end end