Coverage for pySDC/projects/Monodomain/problem_classes/TestODE.py: 91%

1import logging 

2import numpy as np 

3from pySDC.core.problem import Problem 

4from pySDC.core.common import RegisterParams 

5from pySDC.implementations.datatype_classes.mesh import mesh 

6from pySDC.projects.Monodomain.datatype_classes.my_mesh import imexexp_mesh 

7 

8 

9""" 

10Here we define the problems classes for the multirate Dahlquist test equation y'=lambda_I*y + lambda_E*y + lambda_e*y 

11Things are done so that it is compatible witht the sweepers. 

12""" 

13 

14 

15class Parabolic(RegisterParams): 

16 def __init__(self, **problem_params): 

17 self._makeAttributeAndRegister(*problem_params.keys(), localVars=problem_params, readOnly=True) 

18 self.shape = (1,) 

19 self.init = ((1,), None, np.dtype("float64")) 

20 

21 

22class TestODE(Problem): 

23 def __init__(self, **problem_params): 

24 self.logger = logging.getLogger("step") 

25 

26 self.parabolic = Parabolic(**problem_params) 

27 self.size = 1 # one state variable 

28 self.init = ((self.size, *self.parabolic.init[0]), self.parabolic.init[1], self.parabolic.init[2]) 

29 

30 # invoke super init 

31 super(TestODE, self).__init__(self.init) 

32 # store all problem params dictionary values as attributes 

33 self._makeAttributeAndRegister(*problem_params.keys(), localVars=problem_params, readOnly=True) 

34 

35 # initial and end time 

36 self.t0 = 0.0 

37 self.Tend = 1.0 if self.end_time < 0.0 else self.end_time 

38 

39 # set lambdas, if not provided by user 

40 if not hasattr(self, 'lmbda_laplacian'): 

41 self.lmbda_laplacian = -5.0 

42 if not hasattr(self, 'lmbda_gating'): 

43 self.lmbda_gating = -10.0 

44 if not hasattr(self, 'lmbda_others'): 

45 self.lmbda_others = -1.0 

46 

47 self.dtype_u = mesh 

48 self.dtype_f = mesh 

49 

50 def init_exp_extruded(self, new_dim_shape): 

51 return ((*new_dim_shape, 1, self.init[0][1]), self.init[1], self.init[2]) 

52 

53 def initial_value(self): 

54 u0 = self.dtype_u(self.init, val=1.0) 

55 

56 return u0 

57 

58 def eval_f(self, u, t, fh=None): 

59 if fh is None: 

60 fh = self.dtype_f(init=self.init, val=0.0) 

61 

62 fh[0] = (self.lmbda_laplacian + self.lmbda_gating + self.lmbda_others) * u[0] 

63 

64 return fh 

65 

66 

67class MultiscaleTestODE(TestODE): 

68 def __init__(self, **problem_params): 

69 super(MultiscaleTestODE, self).__init__(**problem_params) 

70 

71 self.dtype_f = imexexp_mesh 

72 

73 self.rhs_stiff_indeces = [0] 

74 self.rhs_stiff_args = [0] 

75 self.rhs_nonstiff_indeces = [0] 

76 self.rhs_nonstiff_args = [0] 

77 self.rhs_exp_args = [0] 

78 self.rhs_exp_indeces = [0] 

79 self.rhs_non_exp_indeces = [] 

80 

81 self.constant_lambda_and_phi = True 

82 

83 self.one = self.dtype_u(init=self.init, val=1.0) 

84 

85 def solve_system(self, rhs, factor, u0, t, u_sol=None): 

86 if u_sol is None: 

87 u_sol = self.dtype_u(init=self.init, val=0.0) 

88 

89 u_sol[0] = rhs[0] / (1 - factor * self.lmbda_laplacian) 

90 

91 return u_sol 

92 

93 def eval_f(self, u, t, eval_impl=True, eval_expl=True, eval_exp=True, fh=None, zero_untouched_indeces=True): 

94 

95 if fh is None: 

96 fh = self.dtype_f(init=self.init, val=0.0) 

97 

98 if eval_expl: 

99 fh.expl[0] = self.lmbda_others * u[0] 

100 

101 if eval_impl: 

102 fh.impl[0] = self.lmbda_laplacian * u[0] 

103 

104 if eval_exp: 

105 fh.exp[0] = self.lmbda_gating * u[0] 

106 

107 return fh 

108 

109 def eval_lmbda_yinf_exp(self, u, lmbda, yinf): 

110 lmbda[0] = self.lmbda_gating 

111 yinf[0] = 0.0 

112 

113 def lmbda_eval(self, u, t, lmbda=None): 

114 if lmbda is None: 

115 lmbda = self.dtype_u(init=self.init, val=0.0) 

116 

117 lmbda[0] = self.lmbda_gating 

118 

119 return lmbda