sbplib_julia: test/SbpOperators/volumeops/laplace/laplace

annotate test/SbpOperators/volumeops/laplace/laplace_test.jl @ 1598:19cdec9c21cb feature/boundary_conditions

Implement and test sat_tensors for Dirichlet and Neumann conditions

author	Vidar Stiernström <vidar.stiernstrom@gmail.com>
date	Sun, 26 May 2024 18:19:02 -0700
parents	d68d02dd882f
children	fca4a01d60c9

rev	line source
728 45966c77cb20 Split tests for SbpOperators over several files Jonatan Werpers <jonatan@werpers.com> parents: diff changeset	1 using Test
45966c77cb20 Split tests for SbpOperators over several files Jonatan Werpers <jonatan@werpers.com> parents: diff changeset	2
45966c77cb20 Split tests for SbpOperators over several files Jonatan Werpers <jonatan@werpers.com> parents: diff changeset	3 using Sbplib.SbpOperators
45966c77cb20 Split tests for SbpOperators over several files Jonatan Werpers <jonatan@werpers.com> parents: diff changeset	4 using Sbplib.Grids
732 6114274447f5 Add missing using and imports in test files Jonatan Werpers <jonatan@werpers.com> parents: 728 diff changeset	5 using Sbplib.LazyTensors
1485 e96ee7d7ac9c Test sat_tensors for Laplace w. Neumann conditions Vidar Stiernström <vidar.stiernstrom@it.uu.se> parents: 1291 diff changeset	6 using Sbplib.BoundaryConditions
728 45966c77cb20 Split tests for SbpOperators over several files Jonatan Werpers <jonatan@werpers.com> parents: diff changeset	7
1291 356ec6a72974 Implement changes in SbpOperators Jonatan Werpers <jonatan@werpers.com> parents: 1285 diff changeset	8 @testset "Laplace" begin
1285 7d52c4835d15 Skip broken testsets Jonatan Werpers <jonatan@werpers.com> parents: 1040 diff changeset	9 # Default stencils (4th order)
7d52c4835d15 Skip broken testsets Jonatan Werpers <jonatan@werpers.com> parents: 1040 diff changeset	10 operator_path = sbp_operators_path()*"standard_diagonal.toml"
7d52c4835d15 Skip broken testsets Jonatan Werpers <jonatan@werpers.com> parents: 1040 diff changeset	11 stencil_set = read_stencil_set(operator_path; order=4)
1529 43aaf710463e Change to signature of equidistant_grid to same style as many array methods. Jonatan Werpers <jonatan@werpers.com> parents: 1291 diff changeset	12 g_1D = equidistant_grid(0.0, 1., 101)
43aaf710463e Change to signature of equidistant_grid to same style as many array methods. Jonatan Werpers <jonatan@werpers.com> parents: 1291 diff changeset	13 g_3D = equidistant_grid((0.0, -1.0, 0.0), (1., 1., 1.), 51, 101, 52)
728 45966c77cb20 Split tests for SbpOperators over several files Jonatan Werpers <jonatan@werpers.com> parents: diff changeset	14
45966c77cb20 Split tests for SbpOperators over several files Jonatan Werpers <jonatan@werpers.com> parents: diff changeset	15 @testset "Constructors" begin
750 f88b2117dc69 Merge in default Vidar Stiernström <vidar.stiernstrom@it.uu.se> parents: 732 diff changeset	16 @testset "1D" begin
1291 356ec6a72974 Implement changes in SbpOperators Jonatan Werpers <jonatan@werpers.com> parents: 1285 diff changeset	17 @test Laplace(g_1D, stencil_set) == Laplace(laplace(g_1D, stencil_set), stencil_set)
356ec6a72974 Implement changes in SbpOperators Jonatan Werpers <jonatan@werpers.com> parents: 1285 diff changeset	18 @test Laplace(g_1D, stencil_set) isa LazyTensor{Float64,1,1}
750 f88b2117dc69 Merge in default Vidar Stiernström <vidar.stiernstrom@it.uu.se> parents: 732 diff changeset	19 end
f88b2117dc69 Merge in default Vidar Stiernström <vidar.stiernstrom@it.uu.se> parents: 732 diff changeset	20 @testset "3D" begin
1291 356ec6a72974 Implement changes in SbpOperators Jonatan Werpers <jonatan@werpers.com> parents: 1285 diff changeset	21 @test Laplace(g_3D, stencil_set) == Laplace(laplace(g_3D, stencil_set),stencil_set)
356ec6a72974 Implement changes in SbpOperators Jonatan Werpers <jonatan@werpers.com> parents: 1285 diff changeset	22 @test Laplace(g_3D, stencil_set) isa LazyTensor{Float64,3,3}
728 45966c77cb20 Split tests for SbpOperators over several files Jonatan Werpers <jonatan@werpers.com> parents: diff changeset	23 end
45966c77cb20 Split tests for SbpOperators over several files Jonatan Werpers <jonatan@werpers.com> parents: diff changeset	24 end
45966c77cb20 Split tests for SbpOperators over several files Jonatan Werpers <jonatan@werpers.com> parents: diff changeset	25
45966c77cb20 Split tests for SbpOperators over several files Jonatan Werpers <jonatan@werpers.com> parents: diff changeset	26 # Exact differentiation is measured point-wise. In other cases
45966c77cb20 Split tests for SbpOperators over several files Jonatan Werpers <jonatan@werpers.com> parents: diff changeset	27 # the error is measured in the l2-norm.
45966c77cb20 Split tests for SbpOperators over several files Jonatan Werpers <jonatan@werpers.com> parents: diff changeset	28 @testset "Accuracy" begin
1291 356ec6a72974 Implement changes in SbpOperators Jonatan Werpers <jonatan@werpers.com> parents: 1285 diff changeset	29 l2(v) = sqrt(prod(spacing.(g_3D.grids))*sum(v.^2));
728 45966c77cb20 Split tests for SbpOperators over several files Jonatan Werpers <jonatan@werpers.com> parents: diff changeset	30 polynomials = ()
45966c77cb20 Split tests for SbpOperators over several files Jonatan Werpers <jonatan@werpers.com> parents: diff changeset	31 maxOrder = 4;
45966c77cb20 Split tests for SbpOperators over several files Jonatan Werpers <jonatan@werpers.com> parents: diff changeset	32 for i = 0:maxOrder-1
45966c77cb20 Split tests for SbpOperators over several files Jonatan Werpers <jonatan@werpers.com> parents: diff changeset	33 f_i(x,y,z) = 1/factorial(i)*(y^i + x^i + z^i)
1291 356ec6a72974 Implement changes in SbpOperators Jonatan Werpers <jonatan@werpers.com> parents: 1285 diff changeset	34 polynomials = (polynomials...,eval_on(g_3D,f_i))
728 45966c77cb20 Split tests for SbpOperators over several files Jonatan Werpers <jonatan@werpers.com> parents: diff changeset	35 end
1291 356ec6a72974 Implement changes in SbpOperators Jonatan Werpers <jonatan@werpers.com> parents: 1285 diff changeset	36 # v = eval_on(g_3D, (x,y,z) -> sin(x) + cos(y) + exp(z))
356ec6a72974 Implement changes in SbpOperators Jonatan Werpers <jonatan@werpers.com> parents: 1285 diff changeset	37 # Δv = eval_on(g_3D,(x,y,z) -> -sin(x) - cos(y) + exp(z))
356ec6a72974 Implement changes in SbpOperators Jonatan Werpers <jonatan@werpers.com> parents: 1285 diff changeset	38
356ec6a72974 Implement changes in SbpOperators Jonatan Werpers <jonatan@werpers.com> parents: 1285 diff changeset	39 v = eval_on(g_3D, x̄ -> sin(x̄[1]) + cos(x̄[2]) + exp(x̄[3]))
356ec6a72974 Implement changes in SbpOperators Jonatan Werpers <jonatan@werpers.com> parents: 1285 diff changeset	40 Δv = eval_on(g_3D, x̄ -> -sin(x̄[1]) - cos(x̄[2]) + exp(x̄[3]))
356ec6a72974 Implement changes in SbpOperators Jonatan Werpers <jonatan@werpers.com> parents: 1285 diff changeset	41 @inferred v[1,2,3]
728 45966c77cb20 Split tests for SbpOperators over several files Jonatan Werpers <jonatan@werpers.com> parents: diff changeset	42
45966c77cb20 Split tests for SbpOperators over several files Jonatan Werpers <jonatan@werpers.com> parents: diff changeset	43 # 2nd order interior stencil, 1st order boundary stencil,
45966c77cb20 Split tests for SbpOperators over several files Jonatan Werpers <jonatan@werpers.com> parents: diff changeset	44 # implies that L*v should be exact for binomials up to order 2.
45966c77cb20 Split tests for SbpOperators over several files Jonatan Werpers <jonatan@werpers.com> parents: diff changeset	45 @testset "2nd order" begin
1018 5ec49dd2c7c4 Reintroduce read_stencil_set Vidar Stiernström <vidar.stiernstrom@it.uu.se> parents: 990 diff changeset	46 stencil_set = read_stencil_set(operator_path; order=2)
924 12e8e431b43c Start restructuring Laplace making it more minimal. Vidar Stiernström <vidar.stiernstrom@it.uu.se> parents: 872 diff changeset	47 Δ = Laplace(g_3D, stencil_set)
12e8e431b43c Start restructuring Laplace making it more minimal. Vidar Stiernström <vidar.stiernstrom@it.uu.se> parents: 872 diff changeset	48 @test Δ*polynomials[1] ≈ zeros(Float64, size(g_3D)...) atol = 5e-9
12e8e431b43c Start restructuring Laplace making it more minimal. Vidar Stiernström <vidar.stiernstrom@it.uu.se> parents: 872 diff changeset	49 @test Δ*polynomials[2] ≈ zeros(Float64, size(g_3D)...) atol = 5e-9
12e8e431b43c Start restructuring Laplace making it more minimal. Vidar Stiernström <vidar.stiernstrom@it.uu.se> parents: 872 diff changeset	50 @test Δ*polynomials[3] ≈ polynomials[1] atol = 5e-9
12e8e431b43c Start restructuring Laplace making it more minimal. Vidar Stiernström <vidar.stiernstrom@it.uu.se> parents: 872 diff changeset	51 @test Δ*v ≈ Δv rtol = 5e-2 norm = l2
728 45966c77cb20 Split tests for SbpOperators over several files Jonatan Werpers <jonatan@werpers.com> parents: diff changeset	52 end
45966c77cb20 Split tests for SbpOperators over several files Jonatan Werpers <jonatan@werpers.com> parents: diff changeset	53
45966c77cb20 Split tests for SbpOperators over several files Jonatan Werpers <jonatan@werpers.com> parents: diff changeset	54 # 4th order interior stencil, 2nd order boundary stencil,
45966c77cb20 Split tests for SbpOperators over several files Jonatan Werpers <jonatan@werpers.com> parents: diff changeset	55 # implies that L*v should be exact for binomials up to order 3.
45966c77cb20 Split tests for SbpOperators over several files Jonatan Werpers <jonatan@werpers.com> parents: diff changeset	56 @testset "4th order" begin
1018 5ec49dd2c7c4 Reintroduce read_stencil_set Vidar Stiernström <vidar.stiernstrom@it.uu.se> parents: 990 diff changeset	57 stencil_set = read_stencil_set(operator_path; order=4)
924 12e8e431b43c Start restructuring Laplace making it more minimal. Vidar Stiernström <vidar.stiernstrom@it.uu.se> parents: 872 diff changeset	58 Δ = Laplace(g_3D, stencil_set)
728 45966c77cb20 Split tests for SbpOperators over several files Jonatan Werpers <jonatan@werpers.com> parents: diff changeset	59 # NOTE: high tolerances for checking the "exact" differentiation
45966c77cb20 Split tests for SbpOperators over several files Jonatan Werpers <jonatan@werpers.com> parents: diff changeset	60 # due to accumulation of round-off errors/cancellation errors?
924 12e8e431b43c Start restructuring Laplace making it more minimal. Vidar Stiernström <vidar.stiernstrom@it.uu.se> parents: 872 diff changeset	61 @test Δ*polynomials[1] ≈ zeros(Float64, size(g_3D)...) atol = 5e-9
12e8e431b43c Start restructuring Laplace making it more minimal. Vidar Stiernström <vidar.stiernstrom@it.uu.se> parents: 872 diff changeset	62 @test Δ*polynomials[2] ≈ zeros(Float64, size(g_3D)...) atol = 5e-9
12e8e431b43c Start restructuring Laplace making it more minimal. Vidar Stiernström <vidar.stiernstrom@it.uu.se> parents: 872 diff changeset	63 @test Δ*polynomials[3] ≈ polynomials[1] atol = 5e-9
12e8e431b43c Start restructuring Laplace making it more minimal. Vidar Stiernström <vidar.stiernstrom@it.uu.se> parents: 872 diff changeset	64 @test Δ*polynomials[4] ≈ polynomials[2] atol = 5e-9
12e8e431b43c Start restructuring Laplace making it more minimal. Vidar Stiernström <vidar.stiernstrom@it.uu.se> parents: 872 diff changeset	65 @test Δ*v ≈ Δv rtol = 5e-4 norm = l2
728 45966c77cb20 Split tests for SbpOperators over several files Jonatan Werpers <jonatan@werpers.com> parents: diff changeset	66 end
45966c77cb20 Split tests for SbpOperators over several files Jonatan Werpers <jonatan@werpers.com> parents: diff changeset	67 end
45966c77cb20 Split tests for SbpOperators over several files Jonatan Werpers <jonatan@werpers.com> parents: diff changeset	68 end
924 12e8e431b43c Start restructuring Laplace making it more minimal. Vidar Stiernström <vidar.stiernstrom@it.uu.se> parents: 872 diff changeset	69
1291 356ec6a72974 Implement changes in SbpOperators Jonatan Werpers <jonatan@werpers.com> parents: 1285 diff changeset	70 @testset "laplace" begin
1285 7d52c4835d15 Skip broken testsets Jonatan Werpers <jonatan@werpers.com> parents: 1040 diff changeset	71 operator_path = sbp_operators_path()*"standard_diagonal.toml"
7d52c4835d15 Skip broken testsets Jonatan Werpers <jonatan@werpers.com> parents: 1040 diff changeset	72 stencil_set = read_stencil_set(operator_path; order=4)
1529 43aaf710463e Change to signature of equidistant_grid to same style as many array methods. Jonatan Werpers <jonatan@werpers.com> parents: 1291 diff changeset	73 g_1D = equidistant_grid(0.0, 1., 101)
43aaf710463e Change to signature of equidistant_grid to same style as many array methods. Jonatan Werpers <jonatan@werpers.com> parents: 1291 diff changeset	74 g_3D = equidistant_grid((0.0, -1.0, 0.0), (1., 1., 1.), 51, 101, 52)
1285 7d52c4835d15 Skip broken testsets Jonatan Werpers <jonatan@werpers.com> parents: 1040 diff changeset	75
924 12e8e431b43c Start restructuring Laplace making it more minimal. Vidar Stiernström <vidar.stiernstrom@it.uu.se> parents: 872 diff changeset	76 @testset "1D" begin
1291 356ec6a72974 Implement changes in SbpOperators Jonatan Werpers <jonatan@werpers.com> parents: 1285 diff changeset	77 Δ = laplace(g_1D, stencil_set)
356ec6a72974 Implement changes in SbpOperators Jonatan Werpers <jonatan@werpers.com> parents: 1285 diff changeset	78 @test Δ == second_derivative(g_1D, stencil_set)
356ec6a72974 Implement changes in SbpOperators Jonatan Werpers <jonatan@werpers.com> parents: 1285 diff changeset	79 @test Δ isa LazyTensor{Float64,1,1}
924 12e8e431b43c Start restructuring Laplace making it more minimal. Vidar Stiernström <vidar.stiernstrom@it.uu.se> parents: 872 diff changeset	80 end
12e8e431b43c Start restructuring Laplace making it more minimal. Vidar Stiernström <vidar.stiernstrom@it.uu.se> parents: 872 diff changeset	81 @testset "3D" begin
1291 356ec6a72974 Implement changes in SbpOperators Jonatan Werpers <jonatan@werpers.com> parents: 1285 diff changeset	82 Δ = laplace(g_3D, stencil_set)
356ec6a72974 Implement changes in SbpOperators Jonatan Werpers <jonatan@werpers.com> parents: 1285 diff changeset	83 @test Δ isa LazyTensor{Float64,3,3}
356ec6a72974 Implement changes in SbpOperators Jonatan Werpers <jonatan@werpers.com> parents: 1285 diff changeset	84 Dxx = second_derivative(g_3D, stencil_set, 1)
356ec6a72974 Implement changes in SbpOperators Jonatan Werpers <jonatan@werpers.com> parents: 1285 diff changeset	85 Dyy = second_derivative(g_3D, stencil_set, 2)
356ec6a72974 Implement changes in SbpOperators Jonatan Werpers <jonatan@werpers.com> parents: 1285 diff changeset	86 Dzz = second_derivative(g_3D, stencil_set, 3)
924 12e8e431b43c Start restructuring Laplace making it more minimal. Vidar Stiernström <vidar.stiernstrom@it.uu.se> parents: 872 diff changeset	87 @test Δ == Dxx + Dyy + Dzz
1291 356ec6a72974 Implement changes in SbpOperators Jonatan Werpers <jonatan@werpers.com> parents: 1285 diff changeset	88 @test Δ isa LazyTensor{Float64,3,3}
924 12e8e431b43c Start restructuring Laplace making it more minimal. Vidar Stiernström <vidar.stiernstrom@it.uu.se> parents: 872 diff changeset	89 end
12e8e431b43c Start restructuring Laplace making it more minimal. Vidar Stiernström <vidar.stiernstrom@it.uu.se> parents: 872 diff changeset	90 end
12e8e431b43c Start restructuring Laplace making it more minimal. Vidar Stiernström <vidar.stiernstrom@it.uu.se> parents: 872 diff changeset	91
1485 e96ee7d7ac9c Test sat_tensors for Laplace w. Neumann conditions Vidar Stiernström <vidar.stiernstrom@it.uu.se> parents: 1291 diff changeset	92 @testset "sat_tensors" begin
e96ee7d7ac9c Test sat_tensors for Laplace w. Neumann conditions Vidar Stiernström <vidar.stiernstrom@it.uu.se> parents: 1291 diff changeset	93 operator_path = sbp_operators_path()*"standard_diagonal.toml"
1598 19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	94 orders = (2,4)
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	95 tols = (5e-2,5e-4)
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	96 sz = (201,401)
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	97 g = equidistant_grid((0.,0.), (1.,1.), sz...)
1485 e96ee7d7ac9c Test sat_tensors for Laplace w. Neumann conditions Vidar Stiernström <vidar.stiernstrom@it.uu.se> parents: 1291 diff changeset	98
1598 19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	99 # Verify implementation of sat_tesnors by testing accuracy and symmetry (TODO)
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	100 # of the operator D = Δ + SAT, where SAT is the tensor composition of the
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	101 # operators from sat_tensor. Note that SAT*u should approximate 0 for the
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	102 # conditions chosen.
1485 e96ee7d7ac9c Test sat_tensors for Laplace w. Neumann conditions Vidar Stiernström <vidar.stiernstrom@it.uu.se> parents: 1291 diff changeset	103
1598 19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	104 @testset "Dirichlet" begin
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	105 for (o, tol) ∈ zip(orders,tols)
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	106 stencil_set = read_stencil_set(operator_path; order=o)
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	107 Δ = Laplace(g, stencil_set)
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	108 H = inner_product(g, stencil_set)
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	109 u = collect(eval_on(g, (x,y) -> sin(πx)sin(2π*y)))
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	110 Δu = collect(eval_on(g, (x,y) -> -5π^2sin(πx)sin(2π*y)))
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	111 D = Δ
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	112 for id ∈ boundary_identifiers(g)
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	113 D = D + foldl(∘, sat_tensors(Δ, g, DirichletCondition(0., id)))
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	114 end
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	115 e = D*u .- Δu
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	116 # Accuracy
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	117 @test sqrt(sum(H*e.^2)) ≈ 0 atol = tol
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	118 # Symmetry
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	119 # TODO: # Consider generating the matrices to H and D and test D'H == H'D
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	120 r = randn(size(u))
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	121 @test_broken (D'∘H - H∘D)*r .≈ 0 atol = 1e-13 # TODO: Need to implement apply_transpose for D.
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	122 end
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	123 end
1485 e96ee7d7ac9c Test sat_tensors for Laplace w. Neumann conditions Vidar Stiernström <vidar.stiernstrom@it.uu.se> parents: 1291 diff changeset	124
e96ee7d7ac9c Test sat_tensors for Laplace w. Neumann conditions Vidar Stiernström <vidar.stiernstrom@it.uu.se> parents: 1291 diff changeset	125 @testset "Neumann" begin
1598 19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	126 @testset "Dirichlet" begin
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	127 for (o, tol) ∈ zip(orders,tols)
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	128 stencil_set = read_stencil_set(operator_path; order=o)
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	129 Δ = Laplace(g, stencil_set)
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	130 H = inner_product(g, stencil_set)
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	131 u = collect(eval_on(g, (x,y) -> cos(πx)cos(2π*y)))
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	132 Δu = collect(eval_on(g, (x,y) -> -5π^2cos(πx)cos(2π*y)))
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	133 op = Δ
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	134 for id ∈ boundary_identifiers(g)
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	135 op = op + foldl(∘, sat_tensors(Δ, g, NeumannCondition(0., id)))
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	136 end
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	137 e = op*u .- Δu
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	138 # Accuracy
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	139 @test sqrt(sum(H*e.^2)) ≈ 0 atol = tol
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	140 # Symmetry
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	141 # TODO: # Consider generating the matrices to H and D and test D'H == H'D
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	142 r = randn(size(u))
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	143 @test_broken (D'∘H - H∘D)*r .≈ 0 atol = 1e-13 # TODO: Need to implement apply_transpose for D.
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	144 end
19cdec9c21cb Implement and test sat_tensors for Dirichlet and Neumann conditions Vidar Stiernström <vidar.stiernstrom@gmail.com> parents: 1594 diff changeset	145 end
1485 e96ee7d7ac9c Test sat_tensors for Laplace w. Neumann conditions Vidar Stiernström <vidar.stiernstrom@it.uu.se> parents: 1291 diff changeset	146 end
e96ee7d7ac9c Test sat_tensors for Laplace w. Neumann conditions Vidar Stiernström <vidar.stiernstrom@it.uu.se> parents: 1291 diff changeset	147 end
e96ee7d7ac9c Test sat_tensors for Laplace w. Neumann conditions Vidar Stiernström <vidar.stiernstrom@it.uu.se> parents: 1291 diff changeset	148

Mercurial > repos > public > sbplib_julia

annotate test/SbpOperators/volumeops/laplace/laplace_test.jl @ 1598:19cdec9c21cb feature/boundary_conditions