Coverage for pySDC/implementations/convergence_controller_classes/adaptive_collocation.py: 99%
77 statements
« prev ^ index » next coverage.py v7.16.2, created at 2026-10-17 13:41 +0000
1import numpy as np
3from qmat.lagrange import LagrangeApproximation
4from pySDC.core.convergence_controller import ConvergenceController
7class AdaptiveCollocation(ConvergenceController):
8 """
9 This convergence controller allows to change the underlying quadrature between iterations.
11 Supplying multiple quadrature rules will result in a change to a new quadrature whenever the previous one is
12 converged until all methods given are converged, at which point the step is ended.
13 Whenever the quadrature is changed, the solution is interpolated from the old nodes to the new nodes to ensure
14 accelerated convergence compared to starting from the initial conditions.
16 Use this convergence controller by supplying parameters that the sweeper accepts as a list to the `params`.
17 For instance, supplying
19 .. code-block:: python
21 params = {
22 'num_nodes': [2, 3],
23 }
25 will use collocation methods like you passed to the `sweeper_params` in the `description` object, but will change
26 the number of nodes to 2 before the first iteration and to 3 as soon as the 2-node collocation problem is converged.
27 This will override whatever you set for the number of nodes in the `sweeper_params`, but you still need to set
28 something to allow instantiation of the levels before this convergence controller becomes active.
29 Make sure all lists you supply here have the same length.
31 Feel free to set `logger_level = 15` in the controller parameters to get comprehensive text output on what exactly
32 is happening.
34 This convergence controller has various applications.
35 - You could try to obtain speedup by doing some inexactness. It is currently possible to set various residual
36 tolerances, which will be passed to the levels, corresponding to the accuracy with which each collocation problem
37 is solved.
38 - You can compute multiple solutions to the same initial value problem with different order. This allows, for
39 instance, to do adaptive time stepping.
41 When trying to obtain speedup with this, be ware that the interpolation is not for free. In particular, it is
42 necessary to reevaluate the right hand side at all nodes afterwards.
43 """
45 def setup(self, controller, params, description, **kwargs):
46 """
47 Record what variables we want to vary.
49 Args:
50 controller (pySDC.Controller.controller): The controller
51 params (dict): The params passed for this specific convergence controller
52 description (dict): The description object used to instantiate the controller
54 Returns:
55 (dict): The updated params dictionary
56 """
58 defaults = {
59 'control_order': 300,
60 'num_colls': 0,
61 'sweeper_params': description['sweeper_params'],
62 'vary_keys_sweeper': [],
63 'vary_keys_level': [],
64 }
66 # only these keys can be changed by this convergence controller
67 self.allowed_sweeper_keys = ['quad_type', 'num_nodes', 'node_type', 'do_coll_update']
68 self.allowed_level_keys = ['restol']
70 # add the keys to lists so we know what we need to change later
71 for key in params.keys():
72 if type(params[key]) == list:
73 if key in self.allowed_sweeper_keys:
74 defaults['vary_keys_sweeper'] += [key]
75 elif key in self.allowed_level_keys:
76 defaults['vary_keys_level'] += [key]
77 else:
78 raise NotImplementedError(f'Don\'t know what to do with key {key} here!')
80 defaults['num_colls'] = max([defaults['num_colls'], len(params[key])])
82 self.comm = description['sweeper_params'].get('comm', None)
83 if self.comm:
84 from mpi4py import MPI
86 self.MPI_SUM = MPI.SUM
88 return {**defaults, **super().setup(controller, params, description, **kwargs)}
90 def matmul(self, A, b):
91 """
92 Matrix vector multiplication, possibly MPI parallel.
93 The parallel implementation performs a reduce operation in every row of the matrix. While communicating the
94 entire vector once could reduce the number of communications, this way we never need to store the entire vector
95 on any specific rank.
97 Args:
98 A (2d np.ndarray): Matrix
99 b (list): Vector
101 Returns:
102 List: Axb
103 """
104 if self.comm:
105 res = [A[i, 0] * b[0] if b[i] is not None else None for i in range(A.shape[0])]
106 buf = b[0] * 0.0
107 for i in range(1, A.shape[1]):
108 self.comm.Reduce(A[i, self.comm.rank + 1] * b[self.comm.rank + 1], buf, op=self.MPI_SUM, root=i - 1)
109 if i == self.comm.rank + 1:
110 res[i] += buf
111 return res
112 else:
113 return A @ b
115 def switch_sweeper(self, S):
116 """
117 Update to the next sweeper in line.
119 Args:
120 S (pySDC.Step.step): The current step
122 Returns:
123 None
124 """
126 # generate dictionaries with the new parameters
127 new_params_sweeper = {key: self.params.get(key)[S.status.active_coll] for key in self.params.vary_keys_sweeper}
128 sweeper_params = self.params.sweeper_params.copy()
129 update_params_sweeper = {**sweeper_params, **new_params_sweeper}
131 new_params_level = {key: self.params.get(key)[S.status.active_coll] for key in self.params.vary_keys_level}
133 # update sweeper for all levels
134 for L in S.levels:
135 P = L.prob
137 # store solution of current level which will be interpolated to new level
138 u_old = [me.flatten() if me is not None else me for me in L.u]
139 nodes_old = L.sweep.coll.nodes.copy()
141 # change sweeper
142 L.sweep.__init__(update_params_sweeper, L)
144 # reset level to tell it the new structure of the solution
145 L.params.__dict__.update(new_params_level)
146 L.reset_level(reset_status=False)
148 # interpolate solution of old collocation problem to new one
149 nodes_new = L.sweep.coll.nodes.copy()
150 interpolator = LagrangeApproximation(points=np.append(0, nodes_old))
152 u_inter = self.matmul(interpolator.getInterpolationMatrix(np.append(0, nodes_new)), u_old)
154 # assign the interpolated values to the nodes in the level
155 for i in range(0, len(u_inter)):
156 if u_inter[i] is not None:
157 me = P.dtype_u(P.init)
158 me[:] = np.reshape(u_inter[i], P.init[0])
159 L.u[i] = me
161 # reevaluate rhs
162 for i in range(L.sweep.coll.num_nodes + 1):
163 if L.u[i] is not None:
164 L.f[i] = L.prob.eval_f(L.u[i], L.time)
166 # log the new parameters
167 self.log(f'Switching to collocation {S.status.active_coll + 1} of {self.params.num_colls}', S, level=20)
168 msg = 'New quadrature:'
169 for key in list(sweeper_params.keys()) + list(new_params_level.keys()):
170 if key in self.params.vary_keys_sweeper:
171 msg += f'\n--> {key}: {update_params_sweeper[key]}'
172 elif key in self.params.vary_keys_level:
173 msg += f'\n--> {key}: {new_params_level[key]}'
174 else:
175 msg += f'\n {key}: {update_params_sweeper[key]}'
176 self.log(msg, S)
178 def setup_status_variables(self, controller, **kwargs):
179 """
180 Add an index for which collocation method to use.
182 This lives on the step rather than on this object: every step walks through the collocation
183 methods on its own, and there is one instance of this class per controller, which is a whole
184 block of steps when running with more than one process.
186 Args:
187 controller (pySDC.Controller.controller): The controller
189 Returns:
190 None
191 """
192 self.add_status_variable_to_step('active_coll', 0)
194 def reset_status_variables(self, controller, **kwargs):
195 """
196 Reset the status variables between time steps.
198 Args:
199 controller (pySDC.Controller.controller): The controller
201 Returns:
202 None
203 """
204 self.set_step_status_variable('active_coll', 0)
206 def post_iteration_processing(self, controller, S, **kwargs):
207 """
208 Switch to the next collocation method if the current one is converged.
210 Args:
211 controller (pySDC.Controller.controller): The controller
212 S (pySDC.Step.step): The current step
214 Returns:
215 None
216 """
217 if (S.status.active_coll < self.params.num_colls - 1) and S.status.done:
218 S.status.active_coll += 1
219 S.status.done = False
220 self.switch_sweeper(S)
222 def post_spread_processing(self, controller, S, **kwargs):
223 """
224 Overwrite the sweeper parameters with the first collocation parameters set up here.
226 Args:
227 controller (pySDC.Controller.controller): The controller
228 S (pySDC.Step.step): The current step
230 Returns:
231 None
232 """
233 self.switch_sweeper(S)
235 def check_parameters(self, controller, params, description, **kwargs):
236 """
237 Check if we allow the scheme to solve the collocation problems to convergence.
239 Args:
240 controller (pySDC.Controller): The controller
241 params (dict): The params passed for this specific convergence controller
242 description (dict): The description object used to instantiate the controller
244 Returns:
245 bool: Whether the parameters are compatible
246 str: The error message
247 """
248 if description["level_params"].get("restol", -1.0) <= 1e-16:
249 return (
250 False,
251 "Switching the collocation problems requires solving them to some tolerance that can be reached. Please set attainable `restol` in the level params",
252 )
254 return True, ""