Coverage for pySDC/projects/parallelSDC/preconditioner_playground.py: 97%
142 statements
« prev ^ index » next coverage.py v7.5.0, created at 2024-04-29 09:02 +0000
1import os
2import pickle
3from collections import namedtuple
5import numpy as np
7import pySDC.helpers.plot_helper as plt_helper
8from pySDC.helpers.stats_helper import get_sorted
10from pySDC.implementations.controller_classes.controller_nonMPI import controller_nonMPI
11from pySDC.implementations.problem_classes.AdvectionEquation_ND_FD import advectionNd
12from pySDC.implementations.problem_classes.GeneralizedFisher_1D_FD_implicit import generalized_fisher
13from pySDC.implementations.problem_classes.HeatEquation_ND_FD import heatNd_unforced
14from pySDC.implementations.problem_classes.Van_der_Pol_implicit import vanderpol
15from pySDC.implementations.sweeper_classes.generic_implicit import generic_implicit
17ID = namedtuple('ID', ['setup', 'qd_type', 'param'])
20def main():
21 # initialize level parameters (part I)
22 level_params = dict()
23 level_params['restol'] = 1e-08
25 # initialize sweeper parameters (part I)
26 sweeper_params = dict()
27 sweeper_params['quad_type'] = 'RADAU-RIGHT'
28 sweeper_params['num_nodes'] = 3
30 # initialize step parameters
31 step_params = dict()
32 step_params['maxiter'] = 100
34 # initialize controller parameters
35 controller_params = dict()
36 controller_params['logger_level'] = 30
38 # set up list of Q-delta types and setups
39 qd_list = ['LU', 'IE', 'IEpar', 'Qpar', 'MIN', 'MIN3', 'MIN_GT']
40 setup_list = [
41 ('heat', 63, [10.0**i for i in range(-3, 3)]),
42 ('advection', 64, [10.0**i for i in range(-3, 3)]),
43 ('vanderpol', 2, [0.1 * 2**i for i in range(0, 10)]),
44 ('fisher', 63, [2**i for i in range(-2, 3)]),
45 ]
46 # setup_list = [('fisher', 63, [2 * i for i in range(1, 6)])]
48 # pre-fill results with lists of setups
49 results = dict()
50 for setup, nvars, param_list in setup_list:
51 results[setup] = (nvars, param_list)
53 # loop over all Q-delta matrix types
54 for qd_type in qd_list:
55 # assign implicit Q-delta matrix
56 sweeper_params['QI'] = qd_type
58 # loop over all setups
59 for setup, nvars, param_list in setup_list:
60 # initialize problem parameters (part I)
61 problem_params = dict()
62 if setup != 'vanderpol':
63 problem_params['nvars'] = nvars # number of degrees of freedom for each level
65 # loop over all parameters
66 for param in param_list:
67 # fill description for the controller
68 description = dict()
69 description['sweeper_class'] = generic_implicit # pass sweeper
70 description['sweeper_params'] = sweeper_params # pass sweeper parameters
71 description['step_params'] = step_params # pass step parameters
73 print('working on: %s - %s - %s' % (qd_type, setup, param))
75 # decide which setup to take
76 if setup == 'heat':
77 problem_params['nu'] = param
78 problem_params['freq'] = 2
79 problem_params['bc'] = 'dirichlet-zero' # boundary conditions
81 level_params['dt'] = 0.1
83 description['problem_class'] = heatNd_unforced
84 description['problem_params'] = problem_params
85 description['level_params'] = level_params # pass level parameters
87 elif setup == 'advection':
88 problem_params['c'] = param
89 problem_params['order'] = 2
90 problem_params['freq'] = 2
91 problem_params['stencil_type'] = 'center' # boundary conditions
92 problem_params['bc'] = 'periodic' # boundary conditions
94 level_params['dt'] = 0.1
96 description['problem_class'] = advectionNd
97 description['problem_params'] = problem_params
98 description['level_params'] = level_params # pass level parameters
100 elif setup == 'vanderpol':
101 problem_params['newton_tol'] = 1e-09
102 problem_params['newton_maxiter'] = 20
103 problem_params['mu'] = param
104 problem_params['u0'] = np.array([2.0, 0])
106 level_params['dt'] = 0.1
108 description['problem_class'] = vanderpol
109 description['problem_params'] = problem_params
110 description['level_params'] = level_params
112 elif setup == 'fisher':
113 problem_params['nu'] = 1
114 problem_params['lambda0'] = param
115 problem_params['newton_maxiter'] = 20
116 problem_params['newton_tol'] = 1e-10
117 problem_params['interval'] = (-5, 5)
119 level_params['dt'] = 0.01
121 description['problem_class'] = generalized_fisher
122 description['problem_params'] = problem_params
123 description['level_params'] = level_params
125 else:
126 print('Setup not implemented..', setup)
127 exit()
129 # instantiate controller
130 controller = controller_nonMPI(
131 num_procs=1, controller_params=controller_params, description=description
132 )
134 # get initial values on finest level
135 P = controller.MS[0].levels[0].prob
136 uinit = P.u_exact(0)
138 # call main function to get things done...
139 uend, stats = controller.run(u0=uinit, t0=0, Tend=level_params['dt'])
141 # filter statistics by type (number of iterations)
142 iter_counts = get_sorted(stats, type='niter', sortby='time')
144 # just one time-step, grep number of iteration and store
145 niter = iter_counts[0][1]
146 id = ID(setup=setup, qd_type=qd_type, param=param)
147 results[id] = niter
149 assert len(results) == (6 + 6 + 10 + 5) * 7 + 4, 'ERROR: did not get all results, got %s' % len(results)
151 # write out for later visualization
152 file = open('data/parallelSDC_iterations_precond.pkl', 'wb')
153 pickle.dump(results, file)
155 assert os.path.isfile('data/parallelSDC_iterations_precond.pkl'), 'ERROR: pickle did not create file'
158def plot_iterations():
159 """
160 Helper routine to plot iteration counts
161 """
163 file = open('data/parallelSDC_iterations_precond.pkl', 'rb')
164 results = pickle.load(file)
166 # find the lists/header required for plotting
167 qd_type_list = []
168 setup_list = []
169 for key in results.keys():
170 if isinstance(key, ID):
171 if key.qd_type not in qd_type_list:
172 qd_type_list.append(key.qd_type)
173 elif isinstance(key, str):
174 setup_list.append(key)
175 print('Found these type of preconditioners:', qd_type_list)
176 print('Found these setups:', setup_list)
178 assert len(qd_type_list) == 7, 'ERROR did not find five preconditioners, got %s' % qd_type_list
179 assert len(setup_list) == 4, 'ERROR: did not find three setup, got %s' % setup_list
181 qd_type_list = ['LU', 'IE', 'IEpar', 'Qpar', 'MIN', 'MIN3', 'MIN_GT']
182 marker_list = [None, None, 's', 'o', '^', 'd', 'x']
183 color_list = ['k', 'k', 'r', 'g', 'b', 'c', 'm']
185 plt_helper.setup_mpl()
187 # loop over setups and Q-delta types: one figure per setup, all Qds in one plot
188 for setup in setup_list:
189 plt_helper.newfig(textwidth=238.96, scale=0.89)
191 for qd_type, marker, color in zip(qd_type_list, marker_list, color_list):
192 niter = np.zeros(len(results[setup][1]))
193 for key in results.keys():
194 if isinstance(key, ID):
195 if key.setup == setup and key.qd_type == qd_type:
196 xvalue = results[setup][1].index(key.param)
197 niter[xvalue] = results[key]
198 if qd_type == 'LU':
199 ls = '--'
200 lw = 0.5
201 elif qd_type == 'IE':
202 ls = '-.'
203 lw = 0.5
204 else:
205 ls = '-'
206 lw = 1
207 plt_helper.plt.semilogx(
208 results[setup][1],
209 niter,
210 label=qd_type,
211 lw=lw,
212 linestyle=ls,
213 color=color,
214 marker=marker,
215 markeredgecolor='k',
216 )
218 if setup == 'heat':
219 xlabel = r'$\nu$'
220 elif setup == 'advection':
221 xlabel = r'$c$'
222 elif setup == 'fisher':
223 xlabel = r'$\lambda_0$'
224 elif setup == 'vanderpol':
225 xlabel = r'$\mu$'
226 else:
227 print('Setup not implemented..', setup)
228 exit()
230 plt_helper.plt.ylim([0, 60])
231 plt_helper.plt.legend(loc=2, ncol=1)
232 plt_helper.plt.ylabel('number of iterations')
233 plt_helper.plt.xlabel(xlabel)
234 plt_helper.plt.grid()
236 # save plot as PDF and PGF
237 fname = 'data/parallelSDC_preconditioner_' + setup
238 plt_helper.savefig(fname)
240 assert os.path.isfile(fname + '.pdf'), 'ERROR: plotting did not create PDF file'
241 # assert os.path.isfile(fname + '.pgf'), 'ERROR: plotting did not create PGF file'
242 assert os.path.isfile(fname + '.png'), 'ERROR: plotting did not create PNG file'
245if __name__ == "__main__":
246 main()
247 plot_iterations()