|
| 1 | +import numpy as np |
| 2 | + |
| 3 | +from pySDC.implementations.problem_classes.HeatEquation_1D_FD_forced import heat1d_forced |
| 4 | +from pySDC.implementations.datatype_classes.mesh import mesh, rhs_imex_mesh |
| 5 | +from pySDC.implementations.collocation_classes.gauss_radau_right import CollGaussRadau_Right |
| 6 | +from pySDC.implementations.sweeper_classes.imex_1st_order import imex_1st_order |
| 7 | +from pySDC.implementations.transfer_classes.TransferMesh import mesh_to_mesh |
| 8 | +from pySDC.implementations.controller_classes.allinclusive_multigrid_nonMPI import allinclusive_multigrid_nonMPI |
| 9 | + |
| 10 | +from pySDC.helpers.stats_helper import filter_stats, sort_stats |
| 11 | + |
| 12 | + |
| 13 | +def main(): |
| 14 | + """ |
| 15 | + A simple test program to do PFASST runs for the heat equation |
| 16 | + """ |
| 17 | + |
| 18 | + # initialize level parameters |
| 19 | + level_params = dict() |
| 20 | + level_params['restol'] = 1E-10 |
| 21 | + level_params['dt'] = 0.25 |
| 22 | + |
| 23 | + # initialize sweeper parameters |
| 24 | + sweeper_params = dict() |
| 25 | + sweeper_params['collocation_class'] = CollGaussRadau_Right |
| 26 | + sweeper_params['num_nodes'] = 3 |
| 27 | + sweeper_params['QI'] = 'LU' # For the IMEX sweeper, the LU-trick can be activated for the implicit part |
| 28 | + |
| 29 | + # initialize problem parameters |
| 30 | + problem_params = dict() |
| 31 | + problem_params['nu'] = 0.1 # diffusion coefficient |
| 32 | + problem_params['freq'] = 8 # frequency for the test value |
| 33 | + problem_params['nvars'] = 511 # number of degrees of freedom for each level |
| 34 | + |
| 35 | + # initialize step parameters |
| 36 | + step_params = dict() |
| 37 | + step_params['maxiter'] = 50 |
| 38 | + |
| 39 | + # initialize space transfer parameters |
| 40 | + space_transfer_params = dict() |
| 41 | + space_transfer_params['rorder'] = 2 |
| 42 | + space_transfer_params['iorder'] = 2 |
| 43 | + |
| 44 | + # initialize controller parameters |
| 45 | + controller_params = dict() |
| 46 | + controller_params['logger_level'] = 20 |
| 47 | + |
| 48 | + # fill description dictionary for easy step instantiation |
| 49 | + description = dict() |
| 50 | + description['problem_class'] = heat1d_forced # pass problem class |
| 51 | + description['problem_params'] = problem_params # pass problem parameters |
| 52 | + description['dtype_u'] = mesh # pass data type for u |
| 53 | + description['dtype_f'] = rhs_imex_mesh # pass data type for f |
| 54 | + description['sweeper_class'] = imex_1st_order # pass sweeper |
| 55 | + description['sweeper_params'] = sweeper_params # pass sweeper parameters |
| 56 | + description['level_params'] = level_params # pass level parameters |
| 57 | + description['step_params'] = step_params # pass step parameters |
| 58 | + description['space_transfer_class'] = mesh_to_mesh # pass spatial transfer class |
| 59 | + description['space_transfer_params'] = space_transfer_params # pass paramters for spatial transfer |
| 60 | + |
| 61 | + # set time parameters |
| 62 | + t0 = 0.0 |
| 63 | + Tend = 4.0 |
| 64 | + |
| 65 | + # instantiate controller |
| 66 | + controller = allinclusive_multigrid_nonMPI(num_procs=16, controller_params=controller_params, |
| 67 | + description=description) |
| 68 | + |
| 69 | + # get initial values on finest level |
| 70 | + P = controller.MS[0].levels[0].prob |
| 71 | + uinit = P.u_exact(t0) |
| 72 | + |
| 73 | + # call main function to get things done... |
| 74 | + uend, stats = controller.run(u0=uinit, t0=t0, Tend=Tend) |
| 75 | + |
| 76 | + # compute exact solution and compare |
| 77 | + uex = P.u_exact(Tend) |
| 78 | + err = abs(uex - uend) |
| 79 | + |
| 80 | + print('\nError (vs. exact solution) at time %4.2f: %6.4e\n' % (Tend, err)) |
| 81 | + |
| 82 | + filtered_stats = filter_stats(stats, type='niter') |
| 83 | + |
| 84 | + # convert filtered statistics to list of iterations count, sorted by process |
| 85 | + iter_counts = sort_stats(filtered_stats, sortby='time') |
| 86 | + |
| 87 | + # compute and print statistics |
| 88 | + for item in iter_counts: |
| 89 | + print('Number of iterations for time %4.2f: %2i' % item) |
| 90 | + |
| 91 | + niters = np.array([item[1] for item in iter_counts]) |
| 92 | + |
| 93 | + print(' Mean number of iterations: %4.2f' % np.mean(niters)) |
| 94 | + print(' Range of values for number of iterations: %2i ' % np.ptp(niters)) |
| 95 | + print(' Position of max/min number of iterations: %2i -- %2i' % (int(np.argmax(niters)), int(np.argmin(niters)))) |
| 96 | + |
| 97 | + |
| 98 | +if __name__ == "__main__": |
| 99 | + main() |
0 commit comments