Coverage for pySDC / projects / StroemungsRaum / run_heat_equation_FEniCS.py: 73%

1import json 

2from pathlib import Path 

3import dolfin as df 

4 

5from pySDC.projects.StroemungsRaum.problem_classes.HeatEquation_2D_FEniCS import fenics_heat2D_mass 

6from pySDC.implementations.controller_classes.controller_nonMPI import controller_nonMPI 

7from pySDC.implementations.sweeper_classes.imex_1st_order_mass import imex_1st_order_mass 

8 

9from pySDC.helpers.stats_helper import get_sorted 

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.2 

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'] = 0.1 

46 problem_params['t0'] = t0 

47 problem_params['c_nvars'] = 64 

48 problem_params['family'] = 'CG' 

49 problem_params['order'] = 2 

50 problem_params['c'] = 0.0 

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_heat2D_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 heat equation 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 # quickly generate block of steps 

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

95 

96 # get initial values on finest level 

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

98 uinit = P.u_exact(t0) 

99 

100 # get exact solution at final time for error calculation 

101 uex = P.u_exact(Tend) 

102 

103 # call main function to get things done... 

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

105 

106 # compute relative error at final time 

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

108 

109 return P, stats, rel_err 

110 

111 

112def run_postprocessing(description, problem, stats, Tend): 

113 """ 

114 Postprocess and store simulation results for visualization and analysis. 

115 

116 Args: 

117 description: dict, 

118 pySDC description containing problem parameters. 

119 problem: Problem instance, 

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

121 stats: dict, 

122 collected runtime statistics, 

123 Tend: float, 

124 Final simulation time. 

125 

126 Returns: None 

127 """ 

128 # Get the data directory 

129 import os 

130 

131 path = f"{os.path.dirname(__file__)}/data/heat_equation/" 

132 

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

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

135 

136 # Save parameters 

137 parameters = description['problem_params'] 

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

139 parameters['Tend'] = Tend 

140 json.dump(parameters, open(path + "heat_equation_FEniCS_parameters.json", 'w')) 

141 

142 # Create XDMF file for visualization output 

143 xdmffile_u = df.XDMFFile(path + "heat_equation_FEniCS_Temperature.xdmf") 

144 

145 # Get the solution at every time step, sorted by time 

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

147 

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

149 time = Solutions[i][0] 

150 # 

151 un = Solutions[i][1] 

152 ux = problem.u_exact(time) 

153 # 

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

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

156 # 

157 xdmffile_u.close() 

158 

159 return None 

160 

161 

162if __name__ == "__main__": 

163 

164 t0 = 0.0 

165 Tend = 1.0 

166 

167 # run the setup to get description and controller parameters 

168 description, controller_params = setup(t0=t0) 

169 

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

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

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

173 

174 # run postprocessing to save parameters and solution for visualization 

175 run_postprocessing(description, problem, stats, Tend)