Coverage for pySDC/projects/DAE/misc/problemDAE.py: 100%

1import numpy as np 

2from scipy.optimize import root 

3 

4from pySDC.core.problem import Problem, WorkCounter 

5from pySDC.projects.DAE.misc.meshDAE import MeshDAE 

6 

7 

8class ProblemDAE(Problem): 

9 r""" 

10 This class implements a generic DAE class and illustrates the interface class for DAE problems. 

11 It ensures that all parameters are passed that are needed by DAE sweepers. 

12 

13 Parameters 

14 ---------- 

15 nvars : int 

16 Number of unknowns of the problem class. 

17 newton_tol : float 

18 Tolerance for the nonlinear solver. 

19 

20 Attributes 

21 ---------- 

22 work_counters : WorkCounter 

23 Counts the work, here the number of function calls during the nonlinear solve is logged and stored 

24 in work_counters['newton']. The number of each function class of the right-hand side is then stored 

25 in work_counters['rhs'] 

26 """ 

27 

28 dtype_u = MeshDAE 

29 dtype_f = MeshDAE 

30 

31 def __init__(self, nvars, newton_tol): 

32 """Initialization routine""" 

33 super().__init__((nvars, None, np.dtype('float64'))) 

34 self._makeAttributeAndRegister('nvars', 'newton_tol', localVars=locals(), readOnly=True) 

35 

36 self.work_counters['newton'] = WorkCounter() 

37 self.work_counters['rhs'] = WorkCounter() 

38 

39 def solve_system(self, impl_sys, u_approx, factor, u0, t): 

40 r""" 

41 Solver for nonlinear implicit system (defined in sweeper). 

42 

43 Parameters 

44 ---------- 

45 impl_sys : callable 

46 Implicit system to be solved. 

47 u_approx : dtype_u 

48 Approximation of solution :math:`u` which is needed to solve 

49 the implicit system. 

50 factor : float 

51 Abbrev. for the node-to-node stepsize. 

52 u0 : dtype_u 

53 Initial guess for solver. 

54 t : float 

55 Current time :math:`t`. 

56 

57 Returns 

58 ------- 

59 me : dtype_u 

60 Numerical solution. 

61 """ 

62 me = self.dtype_u(self.init) 

63 

64 def implSysFlatten(unknowns, **kwargs): 

65 sys = impl_sys(unknowns.reshape(me.shape).view(type(u0)), self, factor, u_approx, t, **kwargs) 

66 return sys.flatten() 

67 

68 opt = root( 

69 implSysFlatten, 

70 u0, 

71 method='hybr', 

72 tol=self.newton_tol, 

73 ) 

74 me[:] = opt.x.reshape(me.shape) 

75 self.work_counters['newton'].niter += opt.nfev 

76 return me 

77 

78 def du_exact(self, t): 

79 r""" 

80 Routine for the derivative of the exact solution at time :math:`t \leq 1`. 

81 For this problem, the exact solution is piecewise. 

82 

83 Parameters 

84 ---------- 

85 t : float 

86 Time of the exact solution. 

87 

88 Returns 

89 ------- 

90 me : dtype_u 

91 Derivative of exact solution. 

92 """ 

93 

94 raise NotImplementedError('ERROR: problem has to implement du_exact(self, t)!')