Coverage for pySDC/implementations/convergence_controller_classes/estimate_embedded_error.py: 82%
122 statements
« prev ^ index » next coverage.py v7.5.0, created at 2024-04-29 09:02 +0000
1import numpy as np
3from pySDC.core.ConvergenceController import ConvergenceController, Pars, Status
4from pySDC.implementations.convergence_controller_classes.store_uold import StoreUOld
6from pySDC.implementations.sweeper_classes.Runge_Kutta import RungeKutta
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 """
17 @classmethod
18 def get_implementation(cls, flavor='standard', useMPI=False):
19 """
20 Retrieve the implementation for a specific flavor of this class.
22 Args:
23 flavor (str): The implementation that you want
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!')
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
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
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 }
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.
68 Args:
69 controller (pySDC.Controller): The controller
70 description (dict): The description object used to instantiate the controller
72 Returns:
73 None
74 """
75 if RungeKutta not in description["sweeper_class"].__bases__:
76 controller.add_convergence_controller(StoreUOld, description=description)
78 from pySDC.implementations.hooks.log_embedded_error_estimate import LogEmbeddedErrorEstimate
80 controller.add_hook(LogEmbeddedErrorEstimate)
81 return None
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.
87 Depending on the type of sweeper, the lower order solution is stored in a different place.
89 Args:
90 L (pySDC.level): The level
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 )
111 def setup_status_variables(self, controller, **kwargs):
112 """
113 Add the embedded error variable to the error function.
115 Args:
116 controller (pySDC.Controller): The controller
117 """
118 if 'comm' in kwargs.keys():
119 steps = [controller.S]
120 else:
121 if 'active_slots' in kwargs.keys():
122 steps = [controller.MS[i] for i in kwargs['active_slots']]
123 else:
124 steps = controller.MS
125 where = ["levels", "status"]
126 for S in steps:
127 self.add_variable(S, name='error_embedded_estimate', where=where, init=None)
129 def reset_status_variables(self, controller, **kwargs):
130 self.setup_status_variables(controller, **kwargs)
132 def post_iteration_processing(self, controller, S, **kwargs):
133 """
134 Estimate the local error here.
136 If you are doing MSSDC, this is the global error within the block in Gauss-Seidel mode.
137 In Jacobi mode, I haven't thought about what this is.
139 Args:
140 controller (pySDC.Controller): The controller
141 S (pySDC.Step): The current step
143 Returns:
144 None
145 """
147 if S.status.iter > 0 or self.params.sweeper_type == "RK":
148 for L in S.levels:
149 L.status.error_embedded_estimate = max([self.estimate_embedded_error_serial(L), np.finfo(float).eps])
150 self.debug(f'L.status.error_embedded_estimate={L.status.error_embedded_estimate:.5e}', S)
152 return None
155class EstimateEmbeddedErrorLinearizedNonMPI(EstimateEmbeddedError):
156 def __init__(self, controller, params, description, **kwargs):
157 """
158 Initialisation routine. Add the buffers for communication.
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 super().__init__(controller, params, description, **kwargs)
166 self.buffers = Pars({'e_em_last': 0.0})
168 def setup(self, controller, params, description, **kwargs):
169 """
170 Add option for averaging the local errors.
172 Args:
173 controller (pySDC.Controller): The controller
174 params (dict): Parameters for the convergence controller
175 description (dict): The description object used to instantiate the controller
177 Returns:
178 dict: Updated parameters
179 """
180 return {
181 'averaged': False,
182 **super().setup(controller, params, description, **kwargs),
183 }
185 def reset_buffers_nonMPI(self, controller, **kwargs):
186 """
187 Reset buffers for imitated communication.
189 Args:
190 controller (pySDC.controller): The controller
192 Returns:
193 None
194 """
195 self.buffers.e_em_last = 0.0
196 return None
198 def post_iteration_processing(self, controller, S, **kwargs):
199 """
200 Compute embedded error estimate on the last node of each level
201 In serial this is the local error, but in block Gauss-Seidel MSSDC this is a semi-global error in each block
203 Args:
204 controller (pySDC.Controller): The controller
205 S (pySDC.Step): The current step
207 Returns:
208 None
209 """
210 if len(S.levels) > 1 and len(controller.MS) > 1:
211 raise NotImplementedError(
212 "Embedded error estimate only works for serial multi-level or parallel single \
213level"
214 )
216 if S.status.iter > 0 or self.params.sweeper_type == "RK":
217 if self.params.averaged:
218 averaging = float(S.status.slot + 1)
219 else:
220 averaging = 1.0
222 for L in S.levels:
223 temp = self.estimate_embedded_error_serial(L)
224 L.status.error_embedded_estimate = max(
225 [abs(temp - self.buffers.e_em_last) / averaging, np.finfo(float).eps]
226 )
228 if not self.params.averaged:
229 self.buffers.e_em_last = temp * 1.0
231 return None
234class EstimateEmbeddedErrorLinearizedMPI(EstimateEmbeddedError):
235 def __init__(self, controller, params, description, **kwargs):
236 """
237 Initialisation routine. Add the buffers for communication.
239 Args:
240 controller (pySDC.Controller): The controller
241 params (dict): Parameters for the convergence controller
242 description (dict): The description object used to instantiate the controller
243 """
244 super().__init__(controller, params, description, **kwargs)
245 self.buffers = Pars({'e_em_last': 0.0})
247 def post_iteration_processing(self, controller, S, **kwargs):
248 """
249 Compute embedded error estimate on the last node of each level
250 In serial this is the local error, but in block Gauss-Seidel MSSDC this is a semi-global error in each block
252 Args:
253 controller (pySDC.Controller): The controller
254 S (pySDC.Step): The current step
256 Returns:
257 None
258 """
259 comm = kwargs['comm']
261 if S.status.iter > 0 or self.params.sweeper_type == "RK":
262 for L in S.levels:
263 # get accumulated local errors from previous steps
264 if not S.status.first:
265 if not S.status.prev_done:
266 self.buffers.e_em_last = self.recv(comm, S.status.slot - 1)
267 else:
268 self.buffers.e_em_last = 0.0
270 # estimate accumulated local error
271 temp = self.estimate_embedded_error_serial(L)
273 # estimate local error as difference of accumulated errors
274 L.status.error_embedded_estimate = max([abs(temp - self.buffers.e_em_last), np.finfo(float).eps])
276 # send the accumulated local errors forward
277 if not S.status.last:
278 self.send(comm, dest=S.status.slot + 1, data=temp, blocking=True)
280 return None
283class EstimateEmbeddedErrorCollocation(ConvergenceController):
284 """
285 Estimates an embedded error based on changing the underlying quadrature rule. The error estimate is stored as
286 `error_embedded_estimate_collocation` in the status of the level. Note that we only compute the estimate on the
287 finest level. The error is stored as a tuple with the first index denoting to which iteration it belongs. This
288 is useful since the error estimate is not available immediately after, but only when the next collocation problem
289 is converged to make sure the two solutions are of different accuracy.
291 Changing the collocation method between iterations happens using the `AdaptiveCollocation` convergence controller.
292 Please refer to that for documentation on how to use this. Just pass the parameters for that convergence controller
293 as `adaptive_coll_params` to the parameters for this one and they will be passed on when the `AdaptiveCollocation`
294 convergence controller is automatically added while loading dependencies.
295 """
297 def setup(self, controller, params, description, **kwargs):
298 """
299 Add a default value for control order to the parameters
301 Args:
302 controller (pySDC.Controller): The controller
303 params (dict): Parameters for the convergence controller
304 description (dict): The description object used to instantiate the controller
306 Returns:
307 dict: Updated parameters
308 """
309 defaults = {
310 "control_order": 210,
311 "adaptive_coll_params": {},
312 **super().setup(controller, params, description, **kwargs),
313 }
314 return defaults
316 def dependencies(self, controller, description, **kwargs):
317 """
318 Load the `AdaptiveCollocation` convergence controller to switch between collocation problems between iterations.
320 Args:
321 controller (pySDC.Controller): The controller
322 description (dict): The description object used to instantiate the controller
323 """
324 from pySDC.implementations.convergence_controller_classes.adaptive_collocation import AdaptiveCollocation
326 controller.add_convergence_controller(
327 AdaptiveCollocation, params=self.params.adaptive_coll_params, description=description
328 )
329 from pySDC.implementations.hooks.log_embedded_error_estimate import LogEmbeddedErrorEstimate
331 controller.add_hook(LogEmbeddedErrorEstimate)
333 def post_iteration_processing(self, controller, step, **kwargs):
334 """
335 Compute the embedded error as the difference between the interpolated and the current solution on the finest
336 level.
338 Args:
339 controller (pySDC.Controller.controller): The controller
340 step (pySDC.Step.step): The current step
341 """
342 if step.status.done:
343 lvl = step.levels[0]
344 lvl.sweep.compute_end_point()
345 self.status.u += [lvl.uend]
346 self.status.iter += [step.status.iter]
348 if len(self.status.u) > 1:
349 lvl.status.error_embedded_estimate_collocation = (
350 self.status.iter[-2],
351 max([np.finfo(float).eps, abs(self.status.u[-1] - self.status.u[-2])]),
352 )
354 def setup_status_variables(self, controller, **kwargs):
355 """
356 Add the embedded error variable to the levels and add a status variable for previous steps.
358 Args:
359 controller (pySDC.Controller): The controller
360 """
361 self.status = Status(['u', 'iter'])
362 self.status.u = [] # the solutions of converged collocation problems
363 self.status.iter = [] # the iteration in which the solution converged
365 if 'comm' in kwargs.keys():
366 steps = [controller.S]
367 else:
368 if 'active_slots' in kwargs.keys():
369 steps = [controller.MS[i] for i in kwargs['active_slots']]
370 else:
371 steps = controller.MS
372 where = ["levels", "status"]
373 for S in steps:
374 self.add_variable(S, name='error_embedded_estimate_collocation', where=where, init=None)
376 def reset_status_variables(self, controller, **kwargs):
377 self.setup_status_variables(controller, **kwargs)