Mercurial > repos > public > sbplib_julia
changeset 1021:ee5a641a8277 feature/dissipation_operators
Add some tests and start implementing dissipation()
author | Jonatan Werpers <jonatan@werpers.com> |
---|---|
date | Fri, 18 Mar 2022 15:44:03 +0100 |
parents | 3f5137ce3aa1 |
children | e74c41c4b60e |
files | src/SbpOperators/SbpOperators.jl src/SbpOperators/volumeops/derivatives/dissipation.jl test/SbpOperators/volumeops/derivatives/dissipation_test.jl |
diffstat | 3 files changed, 50 insertions(+), 2 deletions(-) [+] |
line wrap: on
line diff
--- a/src/SbpOperators/SbpOperators.jl Fri Mar 18 12:24:10 2022 +0100 +++ b/src/SbpOperators/SbpOperators.jl Fri Mar 18 15:44:03 2022 +0100 @@ -9,6 +9,7 @@ export normal_derivative export first_derivative export second_derivative +export dissipation using Sbplib.RegionIndices using Sbplib.LazyTensors
--- a/src/SbpOperators/volumeops/derivatives/dissipation.jl Fri Mar 18 12:24:10 2022 +0100 +++ b/src/SbpOperators/volumeops/derivatives/dissipation.jl Fri Mar 18 15:44:03 2022 +0100 @@ -1,3 +1,12 @@ +function dissipation(g::EquidistantGrid, p, direction) + h_inv = inverse_spacing(g)[direction] + + # D = volume_operator(g,CenteredStencil(1),(CenteredStencil(1)), ) + return nothing, nothing +end + +dissipation(g::EquidistantGrid{1}, p) = dissipation(g, p, 1) + function dissipation_interior_weights(p) if p == 0 return (1,)
--- a/test/SbpOperators/volumeops/derivatives/dissipation_test.jl Fri Mar 18 12:24:10 2022 +0100 +++ b/test/SbpOperators/volumeops/derivatives/dissipation_test.jl Fri Mar 18 15:44:03 2022 +0100 @@ -1,8 +1,8 @@ using Test using Sbplib.SbpOperators -# using Sbplib.Grids -# using Sbplib.LazyTensors +using Sbplib.Grids +using Sbplib.LazyTensors using Sbplib.SbpOperators: Stencil @@ -13,6 +13,44 @@ using Sbplib.SbpOperators: dissipation_lower_closure_stencils,dissipation_upper_closure_stencils using Sbplib.SbpOperators: dissipation_transpose_lower_closure_stencils, dissipation_transpose_upper_closure_stencils +""" + monomial(x,k) + +Evaluates ``x^k/k!` with the convetion that it is ``0`` for all ``k<0``. +Has the property that ``d/dx monomial(x,k) = monomial(x,k-1)`` +""" +function monomial(x,k) + if k < 0 + return zero(x) + end + x^k/factorial(k) +end + +@testset "dissipation" begin + g = EquidistantGrid(20, 0., 11.) + D,Dᵀ = dissipation(g, 1) + + @test D isa TensorMapping{Float64,1,1} where T + @test Dᵀ isa TensorMapping{Float64,1,1} where T + + @testset "Accuracy conditions" begin + N = 20 + g = EquidistantGrid(N, 0//1,2//1) + @testset "D_$p" for p ∈ [1,2,3,4] + D,Dᵀ = dissipation(g, p) + + @testset "x^$k" for k ∈ 0:1 + v = evalOn(g, x->monomial(x,k)) + + x, = points(g)[10] + @test (D*v)[10] == monomial(x,k-1) + end + + # Test Dᵀ works backwards and interior forwards + end + end +end + @testset "dissipation_interior_weights" begin @test dissipation_interior_weights(1) == (-1, 1) @test dissipation_interior_weights(2) == (1,-2, 1)