Coverage for pySDC/projects/AsympConv/PFASST_conv_tests.py: 100%
156 statements
« prev ^ index » next coverage.py v7.5.0, created at 2024-04-29 09:02 +0000
1import os
2import pickle
4import matplotlib
5import numpy as np
7matplotlib.use('Agg')
8import matplotlib.pyplot as plt
10from pySDC.implementations.problem_classes.HeatEquation_ND_FD import heatNd_unforced
11from pySDC.implementations.problem_classes.AdvectionEquation_ND_FD import advectionNd
12from pySDC.implementations.sweeper_classes.generic_implicit import generic_implicit
13from pySDC.implementations.transfer_classes.TransferMesh import mesh_to_mesh
14from pySDC.implementations.controller_classes.controller_nonMPI import controller_nonMPI
16from pySDC.helpers.stats_helper import get_sorted
19def main():
20 """
21 Main driver running diffusion and advection tests
22 """
23 nsweeps = 3
24 run_diffusion(nsweeps=nsweeps)
25 run_advection(nsweeps=nsweeps)
26 plot_results(nsweeps=nsweeps)
28 nsweeps = 2
29 run_diffusion(nsweeps=nsweeps)
30 run_advection(nsweeps=nsweeps)
31 plot_results(nsweeps=nsweeps)
34def run_diffusion(nsweeps):
35 """
36 A simple test program to test PFASST convergence for the heat equation with random initial data
38 Args:
39 nsweeps: number of fine sweeps to perform
40 """
42 # initialize level parameters
43 level_params = dict()
44 level_params['restol'] = 1e-08
45 level_params['dt'] = 0.25
46 level_params['nsweeps'] = [nsweeps, 1]
48 # initialize sweeper parameters
49 sweeper_params = dict()
50 sweeper_params['quad_type'] = 'RADAU-RIGHT'
51 sweeper_params['num_nodes'] = [3]
52 sweeper_params['QI'] = ['LU']
53 sweeper_params['initial_guess'] = 'zero'
55 # initialize problem parameters
56 problem_params = dict()
57 problem_params['freq'] = 1 # frequency for the test value
58 problem_params['nvars'] = [127, 63] # number of degrees of freedom for each level
59 problem_params['bc'] = 'dirichlet-zero' # boundary conditions
61 # initialize step parameters
62 step_params = dict()
63 step_params['maxiter'] = 50
65 # initialize space transfer parameters
66 space_transfer_params = dict()
67 space_transfer_params['rorder'] = 2
68 space_transfer_params['iorder'] = 2
69 space_transfer_params['periodic'] = False
71 # initialize controller parameters
72 controller_params = dict()
73 controller_params['logger_level'] = 30
75 # fill description dictionary for easy step instantiation
76 description = dict()
77 description['problem_class'] = heatNd_unforced # pass problem class
78 description['sweeper_class'] = generic_implicit # pass sweeper (see part B)
79 description['sweeper_params'] = sweeper_params # pass sweeper parameters
80 description['level_params'] = level_params # pass level parameters
81 description['step_params'] = step_params # pass step parameters
82 description['space_transfer_class'] = mesh_to_mesh # pass spatial transfer class
83 description['space_transfer_params'] = space_transfer_params # pass paramters for spatial transfer
85 # set time parameters
86 t0 = 0.0
87 Tend = 4 * level_params['dt']
89 # set up number of parallel time-steps to run PFASST with
90 num_proc = 4
92 results = dict()
94 for i in range(-3, 10):
95 ratio = level_params['dt'] / (1.0 / (problem_params['nvars'][0] + 1)) ** 2
97 problem_params['nu'] = 10.0**i / ratio # diffusion coefficient
98 description['problem_params'] = problem_params # pass problem parameters
100 out = 'Working on c = %6.4e' % problem_params['nu']
101 print(out)
102 cfl = ratio * problem_params['nu']
103 out = ' CFL number: %4.2e' % cfl
104 print(out)
106 # instantiate controller
107 controller = controller_nonMPI(num_procs=num_proc, controller_params=controller_params, description=description)
109 # get initial values on finest level
110 P = controller.MS[0].levels[0].prob
111 uinit = P.u_exact(t0)
113 # call main function to get things done...
114 uend, stats = controller.run(u0=uinit, t0=t0, Tend=Tend)
116 # filter statistics by type (number of iterations)
117 iter_counts = get_sorted(stats, type='niter', sortby='time')
119 niters = np.array([item[1] for item in iter_counts])
121 out = ' Mean number of iterations: %4.2f' % np.mean(niters)
122 print(out)
124 if nsweeps == 3 and (i == -3 or i == 9):
125 assert np.mean(niters) <= 2, 'ERROR: too much iterations for diffusive asymptotics, got %s' % np.mean(
126 niters
127 )
129 results[cfl] = np.mean(niters)
131 fname = 'data/results_conv_diffusion_NS' + str(nsweeps) + '.pkl'
132 file = open(fname, 'wb')
133 pickle.dump(results, file)
134 file.close()
136 assert os.path.isfile(fname), 'ERROR: pickle did not create file'
139def run_advection(nsweeps):
140 """
141 A simple test program to test PFASST convergence for the periodic advection equation
143 Args:
144 nsweeps: number of fine sweeps to perform
145 """
147 # initialize level parameters
148 level_params = dict()
149 level_params['restol'] = 1e-08
150 level_params['dt'] = 0.25
151 level_params['nsweeps'] = [nsweeps, 1]
153 # initialize sweeper parameters
154 sweeper_params = dict()
155 sweeper_params['quad_type'] = 'RADAU-RIGHT'
156 sweeper_params['num_nodes'] = [3]
157 sweeper_params['QI'] = ['LU'] # For the IMEX sweeper, the LU-trick can be activated for the implicit part
158 sweeper_params['initial_guess'] = 'zero'
160 # initialize problem parameters
161 problem_params = dict()
162 problem_params['freq'] = 64 # frequency for the test value
163 problem_params['nvars'] = [128, 64] # number of degrees of freedom for each level
164 problem_params['order'] = 2
165 problem_params['stencil_type'] = 'center'
166 problem_params['bc'] = 'periodic' # boundary conditions
168 # initialize step parameters
169 step_params = dict()
170 step_params['maxiter'] = 50
172 # initialize space transfer parameters
173 space_transfer_params = dict()
174 space_transfer_params['rorder'] = 2
175 space_transfer_params['iorder'] = 2
176 space_transfer_params['periodic'] = True
178 # initialize controller parameters
179 controller_params = dict()
180 controller_params['logger_level'] = 30
182 # fill description dictionary for easy step instantiation
183 description = dict()
184 description['problem_class'] = advectionNd # pass problem class
185 description['sweeper_class'] = generic_implicit # pass sweeper (see part B)
186 description['sweeper_params'] = sweeper_params # pass sweeper parameters
187 description['level_params'] = level_params # pass level parameters
188 description['step_params'] = step_params # pass step parameters
189 description['space_transfer_class'] = mesh_to_mesh # pass spatial transfer class
190 description['space_transfer_params'] = space_transfer_params # pass paramters for spatial transfer
192 # set time parameters
193 t0 = 0.0
194 Tend = 4 * level_params['dt']
196 # set up number of parallel time-steps to run PFASST with
197 num_proc = 4
199 results = dict()
201 for i in range(-3, 10):
202 ratio = level_params['dt'] / (1.0 / (problem_params['nvars'][0] + 1))
204 problem_params['c'] = 10.0**i / ratio # diffusion coefficient
205 description['problem_params'] = problem_params # pass problem parameters
207 out = 'Working on nu = %6.4e' % problem_params['c']
208 print(out)
209 cfl = ratio * problem_params['c']
210 out = ' CFL number: %4.2e' % cfl
211 print(out)
213 # instantiate controller
214 controller = controller_nonMPI(num_procs=num_proc, controller_params=controller_params, description=description)
216 # get initial values on finest level
217 P = controller.MS[0].levels[0].prob
218 uinit = P.u_exact(t0)
220 # call main function to get things done...
221 uend, stats = controller.run(u0=uinit, t0=t0, Tend=Tend)
223 # filter statistics by type (number of iterations)
224 iter_counts = get_sorted(stats, type='niter', sortby='time')
226 niters = np.array([item[1] for item in iter_counts])
228 out = ' Mean number of iterations: %4.2f' % np.mean(niters)
229 print(out)
231 if nsweeps == 3 and (i == -3 or i == 9):
232 assert np.mean(niters) <= 2, 'ERROR: too much iterations for advective asymptotics, got %s' % np.mean(
233 niters
234 )
236 results[cfl] = np.mean(niters)
238 fname = 'data/results_conv_advection_NS' + str(nsweeps) + '.pkl'
239 file = open(fname, 'wb')
240 pickle.dump(results, file)
241 file.close()
243 assert os.path.isfile(fname), 'ERROR: pickle did not create file'
246def plot_results(nsweeps):
247 """
248 Plotting routine for iteration counts
250 Args:
251 nsweeps: number of fine sweeps used
252 """
254 fname = 'data/results_conv_diffusion_NS' + str(nsweeps) + '.pkl'
255 file = open(fname, 'rb')
256 results_diff = pickle.load(file)
257 file.close()
259 fname = 'data/results_conv_advection_NS' + str(nsweeps) + '.pkl'
260 file = open(fname, 'rb')
261 results_adv = pickle.load(file)
262 file.close()
264 xvalues_diff = sorted(results_diff.keys())
265 niter_diff = []
266 for x in xvalues_diff:
267 niter_diff.append(results_diff[x])
269 xvalues_adv = sorted(results_adv.keys())
270 niter_adv = []
271 for x in xvalues_adv:
272 niter_adv.append(results_adv[x])
274 # set up plotting parameters
275 params = {
276 'legend.fontsize': 20,
277 'figure.figsize': (12, 8),
278 'axes.labelsize': 20,
279 'axes.titlesize': 20,
280 'xtick.labelsize': 16,
281 'ytick.labelsize': 16,
282 'lines.linewidth': 3,
283 }
284 plt.rcParams.update(params)
286 # set up figure
287 plt.figure()
288 plt.xlabel(r'$\mu$')
289 plt.ylabel('no. of iterations')
290 plt.xlim(min(xvalues_diff + xvalues_adv) / 10.0, max(xvalues_diff + xvalues_adv) * 10.0)
291 plt.ylim(min(niter_diff + niter_adv) - 1, max(niter_diff + niter_adv) + 1)
292 plt.grid()
294 # plot
295 plt.semilogx(xvalues_diff, niter_diff, 'r-', marker='s', markersize=10, label='diffusion')
296 plt.semilogx(xvalues_adv, niter_adv, 'b-', marker='o', markersize=10, label='advection')
298 plt.legend(loc=1, ncol=1, numpoints=1)
300 # plt.show()
301 # save plot, beautify
302 fname = 'data/conv_test_niter_NS' + str(nsweeps) + '.png'
303 plt.savefig(fname, bbox_inches='tight')
305 assert os.path.isfile(fname), 'ERROR: plotting did not create file'
308if __name__ == "__main__":
309 main()