Coverage for pySDC / projects / StroemungsRaum / run_convection_diffusion_equation_FEniCS.py: 72%

1import json 

2from pathlib import Path 

3import dolfin as df 

4 

5from pySDC.projects.StroemungsRaum.problem_classes.ConvectionDiffusion_2D_FEniCS import fenics_ConvDiff2D_mass 

6from pySDC.implementations.sweeper_classes.imex_1st_order_mass import imex_1st_order_mass 

7from pySDC.implementations.controller_classes.controller_nonMPI import controller_nonMPI 

8from pySDC.helpers.stats_helper import get_sorted 

9 

10from pySDC.implementations.hooks.log_solution import LogSolution 

11 

12 

13def setup(t0=None): 

14 """ 

15 Helper routine to set up parameters 

16 

17 Args: 

18 t0: float, 

19 initial time 

20 

21 Returns: 

22 description: dict, 

23 pySDC description dictionary containing problem and method parameters. 

24 controller_params: dict, 

25 Parameters for the pySDC controller. 

26 """ 

27 # time step size 

28 dt = 0.005 

29 

30 # initialize level parameters 

31 level_params = dict() 

32 level_params['restol'] = 1e-12 

33 level_params['dt'] = dt 

34 

35 # initialize step parameters 

36 step_params = dict() 

37 step_params['maxiter'] = 20 

38 

39 # initialize sweeper parameters 

40 sweeper_params = dict() 

41 sweeper_params['quad_type'] = 'RADAU-RIGHT' 

42 sweeper_params['num_nodes'] = 2 

43 

44 problem_params = dict() 

45 problem_params['nu'] = 1e-2 

46 problem_params['t0'] = t0 # ugly, but necessary to set up ProblemClass 

47 problem_params['c_nvars'] = 64 

48 problem_params['family'] = 'CG' 

49 problem_params['order'] = 2 

50 problem_params['sigma'] = 0.05 

51 

52 # initialize controller parameters 

53 controller_params = dict() 

54 controller_params['logger_level'] = 20 

55 controller_params['hook_class'] = LogSolution 

56 

57 description = dict() 

58 

59 description['problem_class'] = fenics_ConvDiff2D_mass 

60 description['sweeper_class'] = imex_1st_order_mass 

61 

62 description['problem_params'] = problem_params 

63 description['sweeper_params'] = sweeper_params 

64 description['level_params'] = level_params 

65 description['step_params'] = step_params 

66 

67 return description, controller_params 

68 

69 

70def run_simulation(description, controller_params, Tend): 

71 """ 

72 Run the time integration for the 2D convection–diffusion benchmark problem. 

73 

74 Args: 

75 description: dict, 

76 pySDC problem and method description. 

77 controller_params: dict, 

78 Parameters for the pySDC controller. 

79 Tend: float, 

80 Final simulation time. 

81 

82 Returns: 

83 P: problem instance, 

84 Problem instance containing the final solution and other problem-related information. 

85 stats: dict, 

86 collected runtime statistics, 

87 rel_err: float, 

88 relative final-time error. 

89 """ 

90 # get initial time from description 

91 t0 = description['problem_params']['t0'] 

92 

93 controller = controller_nonMPI(num_procs=1, controller_params=controller_params, description=description) 

94 

95 # get the initial value 

96 P = controller.MS[0].levels[0].prob 

97 uinit = P.u_exact(t0) 

98 

99 # get exact solution at final time for error calculation 

100 uex = P.u_exact(Tend) 

101 

102 # run the simulation, which returns the final solution and collected statistics 

103 uend, stats = controller.run(u0=uinit, t0=t0, Tend=Tend) 

104 

105 # compute relative error at final time 

106 rel_err = abs(uex - uend) / abs(uex) 

107 

108 return P, stats, rel_err 

109 

110 

111def run_postprocessing(description, problem, stats): 

112 """ 

113 Postprocess and store simulation results for visualization and analysis. 

114 

115 Args: 

116 description: dict, 

117 pySDC description containing problem parameters. 

118 problem: Problem instance, 

119 Problem instance containing the final solution and other problem-related information. 

120 stats: dict, 

121 collected runtime statistics, 

122 

123 Returns: None 

124 """ 

125 # Get the data directory 

126 import os 

127 

128 path = f"{os.path.dirname(__file__)}/data/convection_diffusion/" 

129 

130 # If it does not exist, create the 'data' directory at the specified path, including any necessary parent directories 

131 Path(path).mkdir(parents=True, exist_ok=True) 

132 

133 # get the solution at every time step, sorted by time 

134 Solutions = get_sorted(stats, type='u', sortby='time') 

135 

136 # Save parameters 

137 parameters = description['problem_params'] 

138 parameters.update(description['level_params']) 

139 parameters['Tend'] = Solutions[-1][0] 

140 

141 with open(path + "convection_diffusion_FEniCS_parameters.json", 'w') as f: 

142 json.dump(parameters, f) 

143 

144 # Create XDMF file for visualization output 

145 xdmffile_u = df.XDMFFile(path + "convection_diffusion_FEniCS_solutions.xdmf") 

146 

147 for i in range(len(Solutions)): 

148 time = Solutions[i][0] 

149 # 

150 un = Solutions[i][1] 

151 ux = problem.u_exact(time) 

152 # 

153 xdmffile_u.write_checkpoint(un.values, "un", time, df.XDMFFile.Encoding.HDF5, True) 

154 xdmffile_u.write_checkpoint(ux.values, "ux", time, df.XDMFFile.Encoding.HDF5, True) 

155 # 

156 xdmffile_u.close() 

157 

158 return None 

159 

160 

161if __name__ == "__main__": 

162 

163 t0 = 0.0 

164 Tend = 0.1 

165 

166 # run the setup to get description and controller parameters 

167 description, controller_params = setup(t0=t0) 

168 

169 # run the simulation and get the problem, stats and relative error 

170 problem, stats, rel_err = run_simulation(description, controller_params, Tend) 

171 print('The relative error at time ', Tend, 'is ', rel_err) 

172 

173 # run postprocessing to save parameters and solution for visualization 

174 run_postprocessing(description, problem, stats)