Coverage for pySDC/implementations/convergence_controller_classes/estimate_embedded_error.py: 79%

104 statements  

« prev     ^ index     » next       coverage.py v7.6.1, created at 2024-09-09 14:59 +0000

1import numpy as np 

2 

3from pySDC.core.convergence_controller import ConvergenceController, Pars, Status 

4from pySDC.implementations.convergence_controller_classes.store_uold import StoreUOld 

5 

6from pySDC.implementations.sweeper_classes.Runge_Kutta import RungeKutta 

7 

8 

9class EstimateEmbeddedError(ConvergenceController): 

10 """ 

11 The embedded error is obtained by computing two solutions of different accuracy and pretending the more accurate 

12 one is an exact solution from the point of view of the less accurate solution. In practice, we like to compute the 

13 solutions with different order methods, meaning that in SDC we can just subtract two consecutive sweeps, as long as 

14 you make sure your preconditioner is compatible, which you have to just try out... 

15 """ 

16 

17 @classmethod 

18 def get_implementation(cls, flavor='standard', useMPI=False): 

19 """ 

20 Retrieve the implementation for a specific flavor of this class. 

21 

22 Args: 

23 flavor (str): The implementation that you want 

24 

25 Returns: 

26 cls: The child class that implements the desired flavor 

27 """ 

28 if flavor == 'standard': 

29 return cls 

30 elif flavor == 'linearized': 

31 if useMPI: 

32 return EstimateEmbeddedErrorLinearizedMPI 

33 else: 

34 return EstimateEmbeddedErrorLinearizedNonMPI 

35 elif flavor == 'collocation': 

36 return EstimateEmbeddedErrorCollocation 

37 else: 

38 raise NotImplementedError(f'Flavor {flavor} of EstimateEmbeddedError is not implemented!') 

39 

40 def setup(self, controller, params, description, **kwargs): 

41 """ 

42 Add a default value for control order to the parameters and check if we are using a Runge-Kutta sweeper 

43 

44 Args: 

45 controller (pySDC.Controller): The controller 

46 params (dict): Parameters for the convergence controller 

47 description (dict): The description object used to instantiate the controller 

48 

49 Returns: 

50 dict: Updated parameters 

51 """ 

52 sweeper_type = 'SDC' 

53 if RungeKutta in description['sweeper_class'].__mro__: 

54 sweeper_type = 'RK' 

55 elif 'SweeperMPI' in [me.__name__ for me in description['sweeper_class'].__mro__]: 

56 sweeper_type = 'MPI' 

57 return { 

58 "control_order": -80, 

59 "sweeper_type": sweeper_type, 

60 **super().setup(controller, params, description, **kwargs), 

61 } 

62 

63 def dependencies(self, controller, description, **kwargs): 

64 """ 

65 Load the convergence controller that stores the solution of the last sweep unless we are doing Runge-Kutta. 

66 Add the hook for recording the error. 

67 

68 Args: 

69 controller (pySDC.Controller): The controller 

70 description (dict): The description object used to instantiate the controller 

71 

72 Returns: 

73 None 

74 """ 

75 if RungeKutta not in description["sweeper_class"].__bases__: 

76 controller.add_convergence_controller(StoreUOld, description=description) 

77 

78 from pySDC.implementations.hooks.log_embedded_error_estimate import LogEmbeddedErrorEstimate 

79 

80 controller.add_hook(LogEmbeddedErrorEstimate) 

81 return None 

82 

83 def estimate_embedded_error_serial(self, L): 

84 """ 

85 Estimate the serial embedded error, which may need to be modified for a parallel estimate. 

86 

87 Depending on the type of sweeper, the lower order solution is stored in a different place. 

88 

89 Args: 

90 L (pySDC.level): The level 

91 

92 Returns: 

93 dtype_u: The embedded error estimate 

94 """ 

95 if self.params.sweeper_type == "RK": 

96 # lower order solution is stored in the second to last entry of L.u 

97 return abs(L.u[-2] - L.u[-1]) 

98 elif self.params.sweeper_type == "SDC": 

99 # order rises by one between sweeps, making this so ridiculously easy 

100 return abs(L.uold[-1] - L.u[-1]) 

101 elif self.params.sweeper_type == 'MPI': 

102 comm = L.sweep.comm 

103 return comm.bcast(abs(L.uold[comm.rank + 1] - L.u[comm.rank + 1]), root=comm.size - 1) 

104 else: 

105 raise NotImplementedError( 

106 f" \ 

107 Don't know how to estimate embedded error for sweeper type \ 

108\"{self.params.sweeper_type}\"" 

109 ) 

110 

111 def setup_status_variables(self, controller, **kwargs): 

112 """ 

113 Add the embedded error variable to the error function. 

114 

115 Args: 

116 controller (pySDC.Controller): The controller 

117 """ 

118 self.add_status_variable_to_level('error_embedded_estimate') 

119 

120 def post_iteration_processing(self, controller, S, **kwargs): 

121 """ 

122 Estimate the local error here. 

123 

124 If you are doing MSSDC, this is the global error within the block in Gauss-Seidel mode. 

125 In Jacobi mode, I haven't thought about what this is. 

126 

127 Args: 

128 controller (pySDC.Controller): The controller 

129 S (pySDC.Step): The current step 

130 

131 Returns: 

132 None 

133 """ 

134 

135 if S.status.iter > 0 or self.params.sweeper_type == "RK": 

136 for L in S.levels: 

137 L.status.error_embedded_estimate = max([self.estimate_embedded_error_serial(L), np.finfo(float).eps]) 

138 self.debug(f'L.status.error_embedded_estimate={L.status.error_embedded_estimate:.5e}', S) 

139 

140 return None 

141 

142 

143class EstimateEmbeddedErrorLinearizedNonMPI(EstimateEmbeddedError): 

144 def __init__(self, controller, params, description, **kwargs): 

145 """ 

146 Initialisation routine. Add the buffers for communication. 

147 

148 Args: 

149 controller (pySDC.Controller): The controller 

150 params (dict): Parameters for the convergence controller 

151 description (dict): The description object used to instantiate the controller 

152 """ 

153 super().__init__(controller, params, description, **kwargs) 

154 self.buffers = Pars({'e_em_last': 0.0}) 

155 

156 def setup(self, controller, params, description, **kwargs): 

157 """ 

158 Add option for averaging the local errors. 

159 

160 Args: 

161 controller (pySDC.Controller): The controller 

162 params (dict): Parameters for the convergence controller 

163 description (dict): The description object used to instantiate the controller 

164 

165 Returns: 

166 dict: Updated parameters 

167 """ 

168 return { 

169 'averaged': False, 

170 **super().setup(controller, params, description, **kwargs), 

171 } 

172 

173 def reset_buffers_nonMPI(self, controller, **kwargs): 

174 """ 

175 Reset buffers for imitated communication. 

176 

177 Args: 

178 controller (pySDC.controller): The controller 

179 

180 Returns: 

181 None 

182 """ 

183 self.buffers.e_em_last = 0.0 

184 return None 

185 

186 def post_iteration_processing(self, controller, S, **kwargs): 

187 """ 

188 Compute embedded error estimate on the last node of each level 

189 In serial this is the local error, but in block Gauss-Seidel MSSDC this is a semi-global error in each block 

190 

191 Args: 

192 controller (pySDC.Controller): The controller 

193 S (pySDC.Step): The current step 

194 

195 Returns: 

196 None 

197 """ 

198 if len(S.levels) > 1 and len(controller.MS) > 1: 

199 raise NotImplementedError( 

200 "Embedded error estimate only works for serial multi-level or parallel single \ 

201level" 

202 ) 

203 

204 if S.status.iter > 0 or self.params.sweeper_type == "RK": 

205 if self.params.averaged: 

206 averaging = float(S.status.slot + 1) 

207 else: 

208 averaging = 1.0 

209 

210 for L in S.levels: 

211 temp = self.estimate_embedded_error_serial(L) 

212 L.status.error_embedded_estimate = max( 

213 [abs(temp - self.buffers.e_em_last) / averaging, np.finfo(float).eps] 

214 ) 

215 

216 if not self.params.averaged: 

217 self.buffers.e_em_last = temp * 1.0 

218 

219 return None 

220 

221 

222class EstimateEmbeddedErrorLinearizedMPI(EstimateEmbeddedError): 

223 def __init__(self, controller, params, description, **kwargs): 

224 """ 

225 Initialisation routine. Add the buffers for communication. 

226 

227 Args: 

228 controller (pySDC.Controller): The controller 

229 params (dict): Parameters for the convergence controller 

230 description (dict): The description object used to instantiate the controller 

231 """ 

232 super().__init__(controller, params, description, **kwargs) 

233 self.buffers = Pars({'e_em_last': 0.0}) 

234 

235 def post_iteration_processing(self, controller, S, **kwargs): 

236 """ 

237 Compute embedded error estimate on the last node of each level 

238 In serial this is the local error, but in block Gauss-Seidel MSSDC this is a semi-global error in each block 

239 

240 Args: 

241 controller (pySDC.Controller): The controller 

242 S (pySDC.Step): The current step 

243 

244 Returns: 

245 None 

246 """ 

247 comm = kwargs['comm'] 

248 

249 if S.status.iter > 0 or self.params.sweeper_type == "RK": 

250 for L in S.levels: 

251 # get accumulated local errors from previous steps 

252 if not S.status.first: 

253 if not S.status.prev_done: 

254 self.buffers.e_em_last = self.recv(comm, S.status.slot - 1) 

255 else: 

256 self.buffers.e_em_last = 0.0 

257 

258 # estimate accumulated local error 

259 temp = self.estimate_embedded_error_serial(L) 

260 

261 # estimate local error as difference of accumulated errors 

262 L.status.error_embedded_estimate = max([abs(temp - self.buffers.e_em_last), np.finfo(float).eps]) 

263 

264 # send the accumulated local errors forward 

265 if not S.status.last: 

266 self.send(comm, dest=S.status.slot + 1, data=temp, blocking=True) 

267 

268 return None 

269 

270 

271class EstimateEmbeddedErrorCollocation(ConvergenceController): 

272 """ 

273 Estimates an embedded error based on changing the underlying quadrature rule. The error estimate is stored as 

274 `error_embedded_estimate_collocation` in the status of the level. Note that we only compute the estimate on the 

275 finest level. The error is stored as a tuple with the first index denoting to which iteration it belongs. This 

276 is useful since the error estimate is not available immediately after, but only when the next collocation problem 

277 is converged to make sure the two solutions are of different accuracy. 

278 

279 Changing the collocation method between iterations happens using the `AdaptiveCollocation` convergence controller. 

280 Please refer to that for documentation on how to use this. Just pass the parameters for that convergence controller 

281 as `adaptive_coll_params` to the parameters for this one and they will be passed on when the `AdaptiveCollocation` 

282 convergence controller is automatically added while loading dependencies. 

283 """ 

284 

285 def setup(self, controller, params, description, **kwargs): 

286 """ 

287 Add a default value for control order to the parameters 

288 

289 Args: 

290 controller (pySDC.Controller): The controller 

291 params (dict): Parameters for the convergence controller 

292 description (dict): The description object used to instantiate the controller 

293 

294 Returns: 

295 dict: Updated parameters 

296 """ 

297 defaults = { 

298 "control_order": 210, 

299 "adaptive_coll_params": {}, 

300 **super().setup(controller, params, description, **kwargs), 

301 } 

302 return defaults 

303 

304 def dependencies(self, controller, description, **kwargs): 

305 """ 

306 Load the `AdaptiveCollocation` convergence controller to switch between collocation problems between iterations. 

307 

308 Args: 

309 controller (pySDC.Controller): The controller 

310 description (dict): The description object used to instantiate the controller 

311 """ 

312 from pySDC.implementations.convergence_controller_classes.adaptive_collocation import AdaptiveCollocation 

313 

314 controller.add_convergence_controller( 

315 AdaptiveCollocation, params=self.params.adaptive_coll_params, description=description 

316 ) 

317 from pySDC.implementations.hooks.log_embedded_error_estimate import LogEmbeddedErrorEstimate 

318 

319 controller.add_hook(LogEmbeddedErrorEstimate) 

320 

321 def post_iteration_processing(self, controller, step, **kwargs): 

322 """ 

323 Compute the embedded error as the difference between the interpolated and the current solution on the finest 

324 level. 

325 

326 Args: 

327 controller (pySDC.Controller.controller): The controller 

328 step (pySDC.Step.step): The current step 

329 """ 

330 if step.status.done: 

331 lvl = step.levels[0] 

332 lvl.sweep.compute_end_point() 

333 self.status.u += [lvl.uend] 

334 self.status.iter += [step.status.iter] 

335 

336 if len(self.status.u) > 1: 

337 lvl.status.error_embedded_estimate_collocation = ( 

338 self.status.iter[-2], 

339 max([np.finfo(float).eps, abs(self.status.u[-1] - self.status.u[-2])]), 

340 ) 

341 

342 def setup_status_variables(self, *args, **kwargs): 

343 """ 

344 Add the embedded error variable to the levels and add a status variable for previous steps. 

345 

346 Args: 

347 controller (pySDC.Controller): The controller 

348 """ 

349 self.status = Status(['u', 'iter']) 

350 self.status.u = [] # the solutions of converged collocation problems 

351 self.status.iter = [] # the iteration in which the solution converged 

352 

353 self.add_status_variable_to_level('error_embedded_estimate_collocation')