comparison test/BoundaryConditions/sat_test.jl @ 1482:329720b9ba0d feature/boundary_conditions

Add test for 1D grid
author Vidar Stiernström <vidar.stiernstrom@it.uu.se>
date Mon, 25 Dec 2023 19:23:49 +0100
parents b96858a50e35
children 330c39505a94
comparison
equal deleted inserted replaced
1481:ee242c3fe4af 1482:329720b9ba0d
6 using Sbplib.LazyTensors 6 using Sbplib.LazyTensors
7 using Sbplib.SbpOperators 7 using Sbplib.SbpOperators
8 8
9 stencil_set = read_stencil_set(sbp_operators_path()*"standard_diagonal.toml"; order = 4) 9 stencil_set = read_stencil_set(sbp_operators_path()*"standard_diagonal.toml"; order = 4)
10 10
11 struct MockOp 11 struct MockOp end
12 end
13 12
14 function BoundaryConditions.sat_tensors(op::MockOp, g::TensorGrid, bc::DirichletCondition) 13 function BoundaryConditions.sat_tensors(op::MockOp, g::Grid, bc::DirichletCondition)
15 e = boundary_restriction(g, stencil_set, id(bc)) 14 e = boundary_restriction(g, stencil_set, id(bc))
16 L = e 15 L = e
17 sat_op = e' 16 sat_op = e'
18 return sat_op, L 17 return sat_op, L
19 end 18 end
20 19
21 function BoundaryConditions.sat_tensors(op::MockOp, g::TensorGrid, bc::DirichletCondition, a) 20 function BoundaryConditions.sat_tensors(op::MockOp, g::Grid, bc::DirichletCondition, a)
22 e = boundary_restriction(g, stencil_set, id(bc)) 21 e = boundary_restriction(g, stencil_set, id(bc))
23 L = a*e 22 L = a*e
24 sat_op = e' 23 sat_op = e'
25 return sat_op, L 24 return sat_op, L
26 end 25 end
27 26
28 function BoundaryConditions.sat_tensors(op::MockOp, g::TensorGrid, bc::NeumannCondition) 27 function BoundaryConditions.sat_tensors(op::MockOp, g::Grid, bc::NeumannCondition)
29 e = boundary_restriction(g, stencil_set, id(bc)) 28 e = boundary_restriction(g, stencil_set, id(bc))
30 d = normal_derivative(g, stencil_set, id(bc)) 29 d = normal_derivative(g, stencil_set, id(bc))
31 L = d 30 L = d
32 sat_op = e' 31 sat_op = e'
33 return sat_op, L 32 return sat_op, L
34 end 33 end
35 34
36 @testset "sat" begin 35 @testset "sat" begin
37 op = MockOp() 36 op = MockOp()
38 grid = equidistant_grid((11,13), (0.,0.), (1.,1.)) 37 @testset "1D" begin
39 W, E, S, N = boundary_identifiers(grid) 38 grid = equidistant_grid(11, 0., 1.)
40 u = eval_on(grid, (x,y) -> x+y^2) 39 l, r = boundary_identifiers(grid)
40 u = eval_on(grid, x-> 1. + 2x^2)
41 dc = DirichletCondition(1.0, l)
42 g_l = discretize_data(grid, dc)
43 SAT_l = sat(op, grid, dc)
44 @test SAT_l(u, g_l) ≈ zeros((size(grid))) atol = 1e-13
45
46 nc = NeumannCondition(4.0, r)
47 g_r = discretize_data(grid, nc)
48 SAT_r = sat(op, grid, nc)
49 @test SAT_r(u, g_r) ≈ zeros((size(grid))) atol = 1e-13
50 end
51 @testset "2D" begin
52 grid = equidistant_grid((11,13), (0.,0.), (1.,1.))
53 W, E, S, N = boundary_identifiers(grid)
54 u = eval_on(grid, (x,y) -> x+y^2)
41 55
42 56 dc_W = DirichletCondition(1.0, W)
43 dc_W = DirichletCondition(1.0, W) 57 SAT_W = sat(op, grid, dc_W)
44 SAT_W = sat(op, grid, dc_W) 58 g_W = discretize_data(grid, dc_W)
45 g_W = discretize_data(grid, dc_W) 59 r_W = zeros(size(grid))
46 r_W = zeros(size(grid)) 60 r_W[1,:] .= map(y -> (y^2-1.), range(0., 1., length=13))
47 r_W[1,:] .= map(y -> (y^2-1.), range(0., 1., length=13)) 61 @test SAT_W(u, g_W) ≈ r_W atol = 1e-13
48 @test SAT_W(u, g_W) ≈ r_W atol = 1e-13
49 62
50 dc_E = DirichletCondition(2, E) 63 dc_E = DirichletCondition(2, E)
51 SAT_E = sat(op, grid, dc_E, 2.) 64 SAT_E = sat(op, grid, dc_E, 2.)
52 g_E = discretize_data(grid, dc_E) 65 g_E = discretize_data(grid, dc_E)
53 r_E = zeros(size(grid)) 66 r_E = zeros(size(grid))
54 r_E[end,:] .= map(y -> (2*(1. + y^2)-2.), range(0., 1., length=13)) 67 r_E[end,:] .= map(y -> (2*(1. + y^2)-2.), range(0., 1., length=13))
55 @test SAT_E(u, g_E) ≈ r_E atol = 1e-13 68 @test SAT_E(u, g_E) ≈ r_E atol = 1e-13
56 69
57 nc_S = NeumannCondition(.0, S) 70 nc_S = NeumannCondition(.0, S)
58 SAT_S = sat(op, grid, nc_S) 71 SAT_S = sat(op, grid, nc_S)
59 g_S = discretize_data(grid, nc_S) 72 g_S = discretize_data(grid, nc_S)
60 @test SAT_S(u, g_S) ≈ zeros(size(grid)) atol = 1e-13 73 @test SAT_S(u, g_S) ≈ zeros(size(grid)) atol = 1e-13
61 74
62 nc_N = NeumannCondition(2.0, N) 75 nc_N = NeumannCondition(2.0, N)
63 SAT_S = sat(op, grid, nc_N) 76 SAT_N = sat(op, grid, nc_N)
64 g_N = discretize_data(grid, nc_N) 77 g_N = discretize_data(grid, nc_N)
65 @test SAT_S(u, g_N) ≈ zeros(size(grid)) atol = 1e-13 78 @test SAT_N(u, g_N) ≈ zeros(size(grid)) atol = 1e-13
79 end
66 end 80 end