Coverage for pySDC/implementations/convergence_controller_classes/estimate_embedded_error.py: 79%
107 statements
« prev ^ index » next coverage.py v7.6.9, created at 2024-12-20 14:51 +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 L.sweep.compute_end_point()
97 return abs(L.uend - L.sweep.u_secondary)
98 elif self.params.sweeper_type == "SDC":
99 # order rises by one between sweeps
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"Don't know how to estimate embedded error for sweeper type \
107\"{self.params.sweeper_type}\""
108 )
110 def setup_status_variables(self, controller, **kwargs):
111 """
112 Add the embedded error to the level status
114 Args:
115 controller (pySDC.Controller): The controller
116 """
117 self.add_status_variable_to_level('error_embedded_estimate')
118 self.add_status_variable_to_level('increment')
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 L.status.increment = L.status.error_embedded_estimate * 1
139 self.debug(f'L.status.error_embedded_estimate={L.status.error_embedded_estimate:.5e}', S)
141 return None
144class EstimateEmbeddedErrorLinearizedNonMPI(EstimateEmbeddedError):
145 def __init__(self, controller, params, description, **kwargs):
146 """
147 Initialisation routine. Add the buffers for communication.
149 Args:
150 controller (pySDC.Controller): The controller
151 params (dict): Parameters for the convergence controller
152 description (dict): The description object used to instantiate the controller
153 """
154 super().__init__(controller, params, description, **kwargs)
155 self.buffers = Pars({'e_em_last': 0.0})
157 def setup(self, controller, params, description, **kwargs):
158 """
159 Add option for averaging the local errors.
161 Args:
162 controller (pySDC.Controller): The controller
163 params (dict): Parameters for the convergence controller
164 description (dict): The description object used to instantiate the controller
166 Returns:
167 dict: Updated parameters
168 """
169 return {
170 'averaged': False,
171 **super().setup(controller, params, description, **kwargs),
172 }
174 def reset_buffers_nonMPI(self, controller, **kwargs):
175 """
176 Reset buffers for imitated communication.
178 Args:
179 controller (pySDC.controller): The controller
181 Returns:
182 None
183 """
184 self.buffers.e_em_last = 0.0
185 return None
187 def post_iteration_processing(self, controller, S, **kwargs):
188 """
189 Compute embedded error estimate on the last node of each level
190 In serial this is the local error, but in block Gauss-Seidel MSSDC this is a semi-global error in each block
192 Args:
193 controller (pySDC.Controller): The controller
194 S (pySDC.Step): The current step
196 Returns:
197 None
198 """
199 if len(S.levels) > 1 and len(controller.MS) > 1:
200 raise NotImplementedError(
201 "Embedded error estimate only works for serial multi-level or parallel single \
202level"
203 )
205 if S.status.iter > 0 or self.params.sweeper_type == "RK":
206 if self.params.averaged:
207 averaging = float(S.status.slot + 1)
208 else:
209 averaging = 1.0
211 for L in S.levels:
212 temp = self.estimate_embedded_error_serial(L)
213 L.status.error_embedded_estimate = max(
214 [abs(temp - self.buffers.e_em_last) / averaging, np.finfo(float).eps]
215 )
217 if not self.params.averaged:
218 self.buffers.e_em_last = temp * 1.0
220 return None
223class EstimateEmbeddedErrorLinearizedMPI(EstimateEmbeddedError):
224 def __init__(self, controller, params, description, **kwargs):
225 """
226 Initialisation routine. Add the buffers for communication.
228 Args:
229 controller (pySDC.Controller): The controller
230 params (dict): Parameters for the convergence controller
231 description (dict): The description object used to instantiate the controller
232 """
233 super().__init__(controller, params, description, **kwargs)
234 self.buffers = Pars({'e_em_last': 0.0})
236 def post_iteration_processing(self, controller, S, **kwargs):
237 """
238 Compute embedded error estimate on the last node of each level
239 In serial this is the local error, but in block Gauss-Seidel MSSDC this is a semi-global error in each block
241 Args:
242 controller (pySDC.Controller): The controller
243 S (pySDC.Step): The current step
245 Returns:
246 None
247 """
248 comm = kwargs['comm']
250 if S.status.iter > 0 or self.params.sweeper_type == "RK":
251 for L in S.levels:
252 # get accumulated local errors from previous steps
253 if not S.status.first:
254 if not S.status.prev_done:
255 self.buffers.e_em_last = self.recv(comm, S.status.slot - 1)
256 else:
257 self.buffers.e_em_last = 0.0
259 # estimate accumulated local error
260 temp = self.estimate_embedded_error_serial(L)
262 # estimate local error as difference of accumulated errors
263 L.status.error_embedded_estimate = max([abs(temp - self.buffers.e_em_last), np.finfo(float).eps])
265 # send the accumulated local errors forward
266 if not S.status.last:
267 self.send(comm, dest=S.status.slot + 1, data=temp, blocking=True)
269 return None
272class EstimateEmbeddedErrorCollocation(ConvergenceController):
273 """
274 Estimates an embedded error based on changing the underlying quadrature rule. The error estimate is stored as
275 `error_embedded_estimate_collocation` in the status of the level. Note that we only compute the estimate on the
276 finest level. The error is stored as a tuple with the first index denoting to which iteration it belongs. This
277 is useful since the error estimate is not available immediately after, but only when the next collocation problem
278 is converged to make sure the two solutions are of different accuracy.
280 Changing the collocation method between iterations happens using the `AdaptiveCollocation` convergence controller.
281 Please refer to that for documentation on how to use this. Just pass the parameters for that convergence controller
282 as `adaptive_coll_params` to the parameters for this one and they will be passed on when the `AdaptiveCollocation`
283 convergence controller is automatically added while loading dependencies.
284 """
286 def setup(self, controller, params, description, **kwargs):
287 """
288 Add a default value for control order to the parameters
290 Args:
291 controller (pySDC.Controller): The controller
292 params (dict): Parameters for the convergence controller
293 description (dict): The description object used to instantiate the controller
295 Returns:
296 dict: Updated parameters
297 """
298 defaults = {
299 "control_order": 210,
300 "adaptive_coll_params": {},
301 **super().setup(controller, params, description, **kwargs),
302 }
303 return defaults
305 def dependencies(self, controller, description, **kwargs):
306 """
307 Load the `AdaptiveCollocation` convergence controller to switch between collocation problems between iterations.
309 Args:
310 controller (pySDC.Controller): The controller
311 description (dict): The description object used to instantiate the controller
312 """
313 from pySDC.implementations.convergence_controller_classes.adaptive_collocation import AdaptiveCollocation
315 controller.add_convergence_controller(
316 AdaptiveCollocation, params=self.params.adaptive_coll_params, description=description
317 )
318 from pySDC.implementations.hooks.log_embedded_error_estimate import LogEmbeddedErrorEstimate
320 controller.add_hook(LogEmbeddedErrorEstimate)
322 def post_iteration_processing(self, controller, step, **kwargs):
323 """
324 Compute the embedded error as the difference between the interpolated and the current solution on the finest
325 level.
327 Args:
328 controller (pySDC.Controller.controller): The controller
329 step (pySDC.Step.step): The current step
330 """
331 if step.status.done:
332 lvl = step.levels[0]
333 lvl.sweep.compute_end_point()
334 self.status.u += [lvl.uend]
335 self.status.iter += [step.status.iter]
337 if len(self.status.u) > 1:
338 lvl.status.error_embedded_estimate_collocation = (
339 self.status.iter[-2],
340 max([np.finfo(float).eps, abs(self.status.u[-1] - self.status.u[-2])]),
341 )
343 def setup_status_variables(self, *args, **kwargs):
344 """
345 Add the embedded error variable to the levels and add a status variable for previous steps.
347 Args:
348 controller (pySDC.Controller): The controller
349 """
350 self.status = Status(['u', 'iter'])
351 self.status.u = [] # the solutions of converged collocation problems
352 self.status.iter = [] # the iteration in which the solution converged
354 self.add_status_variable_to_level('error_embedded_estimate_collocation')