15
15
import unittest
16
16
17
17
18
- def errors (pts , function , extrap ):
18
+ def errors (pts , function , method_options ):
19
19
20
20
domain = "test"
21
21
x = pybamm .SpatialVariable ("x" , domain = domain )
@@ -26,10 +26,9 @@ def errors(pts, function, extrap):
26
26
var_pts = {x : pts }
27
27
mesh = pybamm .Mesh (geometry , submesh_types , var_pts )
28
28
29
- spatial_methods = {"test" : pybamm .FiniteVolume }
29
+ spatial_methods = {"test" : pybamm .FiniteVolume ( method_options ) }
30
30
disc = pybamm .Discretisation (mesh , spatial_methods )
31
31
32
- disc .spatial_methods ["test" ].extrapolation = extrap
33
32
var = pybamm .Variable ("var" , domain = "test" )
34
33
left_extrap = pybamm .BoundaryValue (var , "left" )
35
34
right_extrap = pybamm .BoundaryValue (var , "right" )
@@ -47,21 +46,23 @@ def errors(pts, function, extrap):
47
46
return l_error , r_error
48
47
49
48
50
- def get_errors (function , extrap , pts ):
49
+ def get_errors (function , method_options , pts ):
51
50
52
51
l_errors = np .zeros (pts .shape )
53
52
r_errors = np .zeros (pts .shape )
54
53
55
54
for i , pt in enumerate (pts ):
56
- l_errors [i ], r_errors [i ] = errors (pt , function , extrap )
55
+ l_errors [i ], r_errors [i ] = errors (pt , function , method_options )
57
56
58
57
return l_errors , r_errors
59
58
60
59
61
60
class TestExtrapolation (unittest .TestCase ):
62
- def test_quadratic_convergence (self ):
61
+ def test_convergence_without_bcs (self ):
63
62
64
63
# all tests are performed on x in [0, 1]
64
+ linear = {"extrapolation" : {"order" : "linear" }}
65
+ quad = {"extrapolation" : {"order" : "quadratic" }}
65
66
66
67
def x_squared (x ):
67
68
y = x ** 2
@@ -71,8 +72,9 @@ def x_squared(x):
71
72
72
73
pts = 10 ** np .arange (1 , 6 , 1 )
73
74
dx = 1 / pts
74
- l_errors_lin , r_errors_lin = get_errors (x_squared , "linear" , pts )
75
- l_errors_quad , r_errors_quad = get_errors (x_squared , "quadratic" , pts )
75
+
76
+ l_errors_lin , r_errors_lin = get_errors (x_squared , linear , pts )
77
+ l_errors_quad , r_errors_quad = get_errors (x_squared , quad , pts )
76
78
77
79
l_lin_rates = np .log (l_errors_lin [:- 1 ] / l_errors_lin [1 :]) / np .log (
78
80
dx [:- 1 ] / dx [1 :]
@@ -95,7 +97,7 @@ def x_cubed(x):
95
97
r_true = 1
96
98
return y , l_true , r_true
97
99
98
- l_errors_lin , r_errors_lin = get_errors (x_squared , " linear" , pts )
100
+ l_errors_lin , r_errors_lin = get_errors (x_squared , linear , pts )
99
101
100
102
l_lin_rates = np .log (l_errors_lin [:- 1 ] / l_errors_lin [1 :]) / np .log (
101
103
dx [:- 1 ] / dx [1 :]
@@ -111,7 +113,7 @@ def x_cubed(x):
111
113
# quadratic case
112
114
pts = 5 ** np .arange (1 , 7 , 1 )
113
115
dx = 1 / pts
114
- l_errors_quad , r_errors_quad = get_errors (x_cubed , "quadratic" , pts )
116
+ l_errors_quad , r_errors_quad = get_errors (x_cubed , quad , pts )
115
117
116
118
l_quad_rates = np .log (l_errors_quad [:- 1 ] / l_errors_quad [1 :]) / np .log (
117
119
dx [:- 1 ] / dx [1 :]
@@ -124,6 +126,9 @@ def x_cubed(x):
124
126
np .testing .assert_array_almost_equal (l_quad_rates , 3 )
125
127
np .testing .assert_array_almost_equal (r_quad_rates , 3 , decimal = 3 )
126
128
129
+ def test_extrapolation_with_bcs (self ):
130
+ # simple particle with a flux bc
131
+
127
132
128
133
if __name__ == "__main__" :
129
134
print ("Add -v for more debug output" )
0 commit comments