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
3from pySDC.core.convergence_controller 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 self.add_status_variable_to_level('error_embedded_estimate')
120 def post_iteration_processing(self, controller, S, **kwargs):
121 """
122 Estimate the local error here.
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.
127 Args:
128 controller (pySDC.Controller): The controller
129 S (pySDC.Step): The current step
131 Returns:
132 None
133 """
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)
140 return None
143class EstimateEmbeddedErrorLinearizedNonMPI(EstimateEmbeddedError):
144 def __init__(self, controller, params, description, **kwargs):
145 """
146 Initialisation routine. Add the buffers for communication.
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})
156 def setup(self, controller, params, description, **kwargs):
157 """
158 Add option for averaging the local errors.
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
165 Returns:
166 dict: Updated parameters
167 """
168 return {
169 'averaged': False,
170 **super().setup(controller, params, description, **kwargs),
171 }
173 def reset_buffers_nonMPI(self, controller, **kwargs):
174 """
175 Reset buffers for imitated communication.
177 Args:
178 controller (pySDC.controller): The controller
180 Returns:
181 None
182 """
183 self.buffers.e_em_last = 0.0
184 return None
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
191 Args:
192 controller (pySDC.Controller): The controller
193 S (pySDC.Step): The current step
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 )
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
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 )
216 if not self.params.averaged:
217 self.buffers.e_em_last = temp * 1.0
219 return None
222class EstimateEmbeddedErrorLinearizedMPI(EstimateEmbeddedError):
223 def __init__(self, controller, params, description, **kwargs):
224 """
225 Initialisation routine. Add the buffers for communication.
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})
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
240 Args:
241 controller (pySDC.Controller): The controller
242 S (pySDC.Step): The current step
244 Returns:
245 None
246 """
247 comm = kwargs['comm']
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
258 # estimate accumulated local error
259 temp = self.estimate_embedded_error_serial(L)
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])
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)
268 return None
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.
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 """
285 def setup(self, controller, params, description, **kwargs):
286 """
287 Add a default value for control order to the parameters
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
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
304 def dependencies(self, controller, description, **kwargs):
305 """
306 Load the `AdaptiveCollocation` convergence controller to switch between collocation problems between iterations.
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
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
319 controller.add_hook(LogEmbeddedErrorEstimate)
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.
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]
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 )
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.
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
353 self.add_status_variable_to_level('error_embedded_estimate_collocation')