Coverage for pySDC/tutorial/step_4/C_SDC_vs_MLSDC.py: 100%

1from pathlib import Path 

2 

3from pySDC.helpers.stats_helper import get_sorted 

4 

5from pySDC.implementations.controller_classes.controller_nonMPI import controller_nonMPI 

6from pySDC.implementations.problem_classes.HeatEquation_ND_FD import heatNd_unforced 

7from pySDC.implementations.sweeper_classes.generic_implicit import generic_implicit 

8from pySDC.implementations.transfer_classes.TransferMesh import mesh_to_mesh 

9 

10 

11def main(): 

12 """ 

13 A simple test program to compare SDC and MLSDC 

14 """ 

15 

16 # initialize level parameters 

17 level_params = dict() 

18 level_params['restol'] = 1e-09 

19 level_params['dt'] = 0.1 

20 

21 # initialize sweeper parameters 

22 sweeper_params_sdc = dict() 

23 sweeper_params_sdc['node_type'] = 'LEGENDRE' 

24 sweeper_params_sdc['quad_type'] = 'RADAU-RIGHT' 

25 sweeper_params_sdc['num_nodes'] = 5 

26 sweeper_params_sdc['QI'] = 'LU' 

27 

28 sweeper_params_mlsdc = dict() 

29 sweeper_params_mlsdc['node_type'] = 'LEGENDRE' 

30 sweeper_params_mlsdc['quad_type'] = 'RADAU-RIGHT' 

31 sweeper_params_mlsdc['num_nodes'] = [5, 3, 2] 

32 sweeper_params_mlsdc['QI'] = 'LU' 

33 

34 # initialize problem parameters 

35 problem_params_sdc = dict() 

36 problem_params_sdc['nu'] = 0.1 # diffusion coefficient 

37 problem_params_sdc['freq'] = 4 # frequency for the test value 

38 problem_params_sdc['nvars'] = 1023 # number of degrees of freedom for each level 

39 problem_params_sdc['bc'] = 'dirichlet-zero' # boundary conditions 

40 

41 problem_params_mlsdc = dict() 

42 problem_params_mlsdc['nu'] = 0.1 # diffusion coefficient 

43 problem_params_mlsdc['freq'] = 4 # frequency for the test value 

44 problem_params_mlsdc['nvars'] = [1023, 511, 255] # number of degrees of freedom for each level 

45 problem_params_mlsdc['bc'] = 'dirichlet-zero' # boundary conditions 

46 

47 # initialize step parameters 

48 step_params = dict() 

49 step_params['maxiter'] = 20 

50 

51 # initialize space transfer parameters 

52 space_transfer_params = dict() 

53 space_transfer_params['rorder'] = 2 

54 space_transfer_params['iorder'] = 6 

55 

56 # initialize controller parameters 

57 controller_params = dict() 

58 controller_params['logger_level'] = 30 

59 

60 # fill description dictionary for SDC 

61 description_sdc = dict() 

62 description_sdc['problem_class'] = heatNd_unforced 

63 description_sdc['problem_params'] = problem_params_sdc 

64 description_sdc['sweeper_class'] = generic_implicit 

65 description_sdc['sweeper_params'] = sweeper_params_sdc 

66 description_sdc['level_params'] = level_params 

67 description_sdc['step_params'] = step_params 

68 

69 # fill description dictionary for MLSDC 

70 description_mlsdc = dict() 

71 description_mlsdc['problem_class'] = heatNd_unforced 

72 description_mlsdc['problem_params'] = problem_params_mlsdc 

73 description_mlsdc['sweeper_class'] = generic_implicit 

74 description_mlsdc['sweeper_params'] = sweeper_params_mlsdc 

75 description_mlsdc['level_params'] = level_params 

76 description_mlsdc['step_params'] = step_params 

77 description_mlsdc['space_transfer_class'] = mesh_to_mesh 

78 description_mlsdc['space_transfer_params'] = space_transfer_params 

79 

80 # instantiate the controller (no controller parameters used here) 

81 controller_sdc = controller_nonMPI(num_procs=1, controller_params=controller_params, description=description_sdc) 

82 controller_mlsdc = controller_nonMPI( 

83 num_procs=1, controller_params=controller_params, description=description_mlsdc 

84 ) 

85 

86 # set time parameters 

87 t0 = 0.0 

88 Tend = 0.1 

89 

90 # get initial values on finest level 

91 P = controller_sdc.MS[0].levels[0].prob 

92 uinit = P.u_exact(t0) 

93 

94 # call main functions to get things done... 

95 uend_sdc, stats_sdc = controller_sdc.run(u0=uinit, t0=t0, Tend=Tend) 

96 uend_mlsdc, stats_mlsdc = controller_mlsdc.run(u0=uinit, t0=t0, Tend=Tend) 

97 

98 # get number of iterations for both 

99 niter_sdc = get_sorted(stats_sdc, type='niter', sortby='time')[0][1] 

100 niter_mlsdc = get_sorted(stats_mlsdc, type='niter', sortby='time')[0][1] 

101 

102 # compute exact solution and compare both 

103 uex = P.u_exact(Tend) 

104 err_sdc = abs(uex - uend_sdc) 

105 err_mlsdc = abs(uex - uend_mlsdc) 

106 diff = abs(uend_mlsdc - uend_sdc) 

107 

108 # print out and check 

109 Path("data").mkdir(parents=True, exist_ok=True) 

110 f = open('data/step_4_C_out.txt', 'a') 

111 out = 'Error SDC and MLSDC: %12.8e -- %12.8e' % (err_sdc, err_mlsdc) 

112 f.write(out + '\n') 

113 print(out) 

114 out = 'Difference SDC vs. MLSDC: %12.8e' % diff 

115 f.write(out + '\n') 

116 print(out) 

117 out = 'Number of iterations SDC and MLSDC: %2i -- %2i' % (niter_sdc, niter_mlsdc) 

118 f.write(out + '\n') 

119 print(out) 

120 

121 assert diff < 6e-10, "ERROR: difference between MLSDC and SDC is higher than expected, got %s" % diff 

122 assert niter_sdc - niter_mlsdc <= 6, "ERROR: MLSDC required more iterations than expected, got %s" % niter_mlsdc 

123 

124 

125if __name__ == "__main__": 

126 main()