Coverage for pySDC/projects/DAE/sweepers/semiImplicitDAEMPI.py: 100%

1from mpi4py import MPI 

2 

3from pySDC.projects.DAE.sweepers.fullyImplicitDAEMPI import SweeperDAEMPI 

4from pySDC.projects.DAE.sweepers.semiImplicitDAE import SemiImplicitDAE 

5 

6 

7class SemiImplicitDAEMPI(SweeperDAEMPI, SemiImplicitDAE): 

8 r""" 

9 Custom sweeper class to implement the fully-implicit SDC parallelized across the nodes for solving fully-implicit DAE problems of the form 

10 

11 .. math:: 

12 u' = f(u, z, t), 

13 

14 .. math:: 

15 0 = g(u, z, t) 

16 

17 More detailed description can be found in ``semiImplicitDAE.py``. To parallelize SDC across the method the idea is to use a diagonal :math:`\mathbf{Q}_\Delta` to have solves of the implicit system on each node that can be done in parallel since they are fully decoupled from values of previous nodes. First such diagonal :math:`\mathbf{Q}_\Delta` were developed in [1]_. Years later intensive theory about the topic was developed in [2]_. For the DAE case these ideas were basically transferred. 

18 

19 Note 

20 ---- 

21 For more details of implementing a sweeper enabling parallelization across the method we refer to documentation in [generic_implicit_MPI.py](https://github.com/Parallel-in-Time/pySDC/blob/master/pySDC/implementations/sweeper_classes/generic_implicit_MPI.py) and [fullyImplicitDAEMPI.py](https://github.com/Parallel-in-Time/pySDC/blob/master/pySDC/projects/DAE/sweepers/fullyImplicitDAEMPI.py). 

22 

23 For parallelization across the method for semi-explicit DAEs, differential and algebraic parts have to be treated separately. In ``integrate()`` these parts needs to be collected separately, otherwise information would get lost. 

24 

25 Reference 

26 --------- 

27 .. [1] R. Speck. Parallelizing spectral deferred corrections across the method. Comput. Vis. Sci. 19, No. 3-4, 75-83 (2018). 

28 .. [2] G. Čaklović, T. Lunet, S. Götschel, D. Ruprecht. Improving Efficiency of Parallel Across the Method Spectral Deferred Corrections. Preprint, arXiv:2403.18641 [math.NA] (2024). 

29 """ 

30 

31 def integrate(self, last_only=False): 

32 """ 

33 Integrates the gradient. Note that here only the differential part is integrated, i.e., the 

34 integral over the algebraic part is zero. ``me`` serves as buffer, and root process ``m`` 

35 stores the result of integral at node :math:`\tau_m` there. 

36 

37 Parameters 

38 ---------- 

39 last_only : bool, optional 

40 Integrate only the last node for the residual or all of them. 

41 

42 Returns 

43 ------- 

44 me : list of dtype_u 

45 Containing the integral as values. 

46 """ 

47 

48 L = self.level 

49 P = L.prob 

50 

51 me = P.dtype_u(P.init, val=0.0) 

52 for m in [self.coll.num_nodes - 1] if last_only else range(self.coll.num_nodes): 

53 integral = P.dtype_u(P.init, val=0.0) 

54 integral.diff[:] = L.dt * self.coll.Qmat[m + 1, self.rank + 1] * L.f[self.rank + 1].diff[:] 

55 recvBuf = me[:] if m == self.rank else None 

56 self.comm.Reduce(integral, recvBuf, root=m, op=MPI.SUM) 

57 

58 return me 

59 

60 def update_nodes(self): 

61 r""" 

62 Updates values of ``u`` and ``f`` at collocation nodes. This correspond to a single iteration of the 

63 preconditioned Richardson iteration in **"ordinary"** SDC. 

64 """ 

65 

66 L = self.level 

67 P = L.prob 

68 

69 # only if the level has been touched before 

70 assert L.status.unlocked 

71 

72 integral = self.integrate() 

73 integral.diff[:] -= L.dt * self.QI[self.rank + 1, self.rank + 1] * L.f[self.rank + 1].diff[:] 

74 integral.diff[:] += L.u[0].diff[:] 

75 

76 u_approx = P.dtype_u(integral) 

77 

78 u0 = P.dtype_u(P.init) 

79 u0.diff[:], u0.alg[:] = L.f[self.rank + 1].diff[:], L.u[self.rank + 1].alg[:] 

80 

81 u_new = P.solve_system( 

82 SemiImplicitDAE.F, 

83 u_approx, 

84 L.dt * self.QI[self.rank + 1, self.rank + 1], 

85 u0, 

86 L.time + L.dt * self.coll.nodes[self.rank], 

87 ) 

88 

89 L.f[self.rank + 1].diff[:] = u_new.diff[:] 

90 L.u[self.rank + 1].alg[:] = u_new.alg[:] 

91 

92 integral = self.integrate() 

93 L.u[self.rank + 1].diff[:] = L.u[0].diff[:] + integral.diff[:] 

94 

95 # indicate presence of new values at this level 

96 L.status.updated = True 

97 

98 return None