Coverage for pySDC / projects / RayleighBenard / RBC3D_configs.py: 77%

1from pySDC.projects.GPU.configs.base_config import Config 

2 

3 

4def get_config(args): 

5 name = args['config'] 

6 if name == 'RBC3D': 

7 return RayleighBenard3DRegular(args) 

8 elif name in globals().keys(): 

9 return globals()[name](args) 

10 else: 

11 raise NotImplementedError(f'There is no configuration called {name!r}!') 

12 

13 

14class RayleighBenard3DRegular(Config): 

15 sweeper_type = 'IMEX' 

16 Tend = 50 

17 gamma = 1 

18 res_ratio = 1 

19 dealiasing = 3.0 / 2.0 

20 

21 def get_sweeper(self, useMPI): 

22 if useMPI: 

23 from pySDC.projects.RayleighBenard.sweepers import imex_1st_order_MPI_fixed_k as sweeper 

24 else: 

25 from pySDC.projects.RayleighBenard.sweepers import imex_1st_order_diagonal_serial as sweeper 

26 return sweeper 

27 

28 def get_file_name(self): 

29 res = self.args['res'] 

30 return f'{self.base_path}/data/{type(self).__name__}-res{res}.pySDC' 

31 

32 def get_LogToFile(self, *args, **kwargs): 

33 if self.comms[1].rank > 0: 

34 return None 

35 import numpy as np 

36 from pySDC.implementations.hooks.log_solution import LogToFile 

37 

38 LogToFile.filename = self.get_file_name() 

39 LogToFile.time_increment = 5e-1 

40 # LogToFile.allow_overwriting = True 

41 

42 return LogToFile 

43 

44 def get_controller_params(self, *args, **kwargs): 

45 from pySDC.implementations.hooks.log_step_size import LogStepSize 

46 

47 controller_params = super().get_controller_params(*args, **kwargs) 

48 controller_params['hook_class'] += [LogStepSize] 

49 return controller_params 

50 

51 def get_description(self, *args, MPIsweeper=False, res=-1, **kwargs): 

52 from pySDC.implementations.problem_classes.RayleighBenard3D import ( 

53 RayleighBenard3D, 

54 ) 

55 from pySDC.implementations.problem_classes.generic_spectral import ( 

56 compute_residual_DAE, 

57 compute_residual_DAE_MPI, 

58 ) 

59 from pySDC.implementations.convergence_controller_classes.step_size_limiter import StepSizeSlopeLimiter 

60 from pySDC.implementations.convergence_controller_classes.crash import StopAtNan 

61 

62 desc = super().get_description(*args, MPIsweeper=MPIsweeper, **kwargs) 

63 

64 if MPIsweeper: 

65 desc['sweeper_class'].compute_residual = compute_residual_DAE_MPI 

66 else: 

67 desc['sweeper_class'].compute_residual = compute_residual_DAE 

68 

69 desc['level_params']['dt'] = 0.01 

70 desc['level_params']['restol'] = 1e-7 

71 

72 desc['convergence_controllers'][StepSizeSlopeLimiter] = {'dt_rel_min_slope': 0.1} 

73 desc['convergence_controllers'][StopAtNan] = {} 

74 

75 desc['sweeper_params']['quad_type'] = 'RADAU-RIGHT' 

76 desc['sweeper_params']['num_nodes'] = 2 

77 desc['sweeper_params']['QI'] = 'MIN-SR-S' 

78 desc['sweeper_params']['QE'] = 'PIC' 

79 

80 res = 64 if res == -1 else res 

81 desc['problem_params']['Rayleigh'] = 1e8 

82 desc['problem_params']['nx'] = self.res_ratio * res 

83 desc['problem_params']['ny'] = self.res_ratio * res 

84 desc['problem_params']['nz'] = res 

85 desc['problem_params']['Lx'] = self.gamma 

86 desc['problem_params']['Ly'] = self.gamma 

87 desc['problem_params']['Lz'] = 1 

88 desc['problem_params']['heterogeneous'] = True 

89 desc['problem_params']['dealiasing'] = self.dealiasing 

90 

91 desc['step_params']['maxiter'] = 3 

92 

93 desc['problem_class'] = RayleighBenard3D 

94 

95 return desc 

96 

97 def get_initial_condition(self, P, *args, restart_idx=0, **kwargs): 

98 

99 if restart_idx == 0: 

100 u0 = P.u_exact(t=0, seed=P.comm.rank, noise_level=1e-3) 

101 u0_with_pressure = P.solve_system(u0, 1e-9, u0) 

102 P.cached_factorizations.pop(1e-9) 

103 return u0_with_pressure, 0 

104 else: 

105 from pySDC.helpers.fieldsIO import FieldsIO 

106 

107 P.setUpFieldsIO() 

108 outfile = FieldsIO.fromFile(self.get_file_name()) 

109 

110 t0, solution = outfile.readField(restart_idx) 

111 solution = solution[: P.spectral.ncomponents, ...] 

112 

113 u0 = P.u_init 

114 

115 if P.spectral_space: 

116 u0[...] = P.transform(solution) 

117 else: 

118 u0[...] = solution 

119 

120 return u0, t0 

121 

122 def prepare_caches(self, prob): 

123 """ 

124 Cache the fft objects, which are expensive to create on GPU because graphs have to be initialized. 

125 """ 

126 prob.eval_f(prob.u_init) 

127 

128 def prepare_for_benchmark(self): 

129 def _pass(*args, **kwargs): 

130 pass 

131 

132 self.get_LogToFile = _pass 

133 

134 def prepare_caches_for_benchmark(self, prob, controller): 

135 _rhs = prob.u_init 

136 sweeper = controller.MS[0].levels[0].sweep 

137 _dt = sweeper.level.dt 

138 

139 hooks = controller.hooks 

140 

141 # mute controller 

142 type(controller).hooks = [] 

143 controller.run(_rhs, 0, _dt) 

144 

145 # unmute controller 

146 controller.hooks = hooks 

147 

148 try: 

149 import cuda as cp 

150 

151 cp.cuda.get_current_stream().synchronize() 

152 except ModuleNotFoundError: 

153 pass 

154 from mpi4py import MPI 

155 

156 MPI.COMM_WORLD.Barrier() 

157 controller.logger.critical('Set up caches for benchmarking') 

158 

159 def prepare_description_for_benchmark(self, description, controller_params): 

160 from pySDC.projects.RayleighBenard.benchmarks.print_timings_hook import PrintCPUTimings, PrintGPUTimings 

161 from mpi4py import MPI 

162 

163 self.Tend = 5 * description['level_params']['dt'] 

164 

165 controller_params['logger_level'] = 40 

166 controller_params['hook_class'] += [PrintCPUTimings] 

167 if description['problem_params']['useGPU']: 

168 controller_params['hook_class'] += [PrintGPUTimings] 

169 

170 description['problem_params']['max_cached_factorizations'] = 99 

171 

172 time_rank = 0 

173 if 'comm' in description['sweeper_params'].keys(): 

174 time_rank = description['sweeper_params']['comm'].rank 

175 for i in range(MPI.COMM_WORLD.size): 

176 if MPI.COMM_WORLD.rank == i: 

177 print( 

178 f'Global rank {MPI.COMM_WORLD.rank} is {time_rank} in time and {description["problem_params"]["comm"].rank} in space', 

179 flush=True, 

180 ) 

181 MPI.COMM_WORLD.barrier() 

182 

183 

184class RBC3Dverification(RayleighBenard3DRegular): 

185 converged = 0 

186 dt = 1e-2 

187 ic_config = { 

188 'config': None, 

189 'res': -1, 

190 'dt': -1, 

191 } 

192 res = None 

193 Ra = None 

194 Tend = 100 

195 res_ratio = 4 

196 gamma = 4 

197 

198 def get_file_name(self): 

199 res = self.args['res'] 

200 dt = self.args['dt'] 

201 return f'{self.base_path}/data/{type(self).__name__}-res{res}-dt{dt:.0e}.pySDC' 

202 

203 def get_description(self, *args, res=-1, dt=-1, **kwargs): 

204 desc = super().get_description(*args, **kwargs) 

205 desc['level_params']['nsweeps'] = 4 

206 desc['level_params']['restol'] = -1 

207 desc['step_params']['maxiter'] = 1 

208 desc['sweeper_params']['QI'] = 'MIN-SR-S' 

209 desc['sweeper_params']['skip_residual_computation'] = ('IT_CHECK', 'IT_DOWN', 'IT_UP', 'IT_FINE', 'IT_COARSE') 

210 desc['sweeper_params']['num_nodes'] = 4 

211 Ra = int(type(self).__name__[-3]) * 10 ** int(type(self).__name__[-1]) 

212 desc['problem_params']['Rayleigh'] = Ra 

213 desc['problem_params']['Prandtl'] = 0.7 

214 

215 _res = self.res if res == -1 else res 

216 desc['problem_params']['nx'] = _res * self.res_ratio 

217 desc['problem_params']['ny'] = _res * self.res_ratio 

218 desc['problem_params']['nz'] = _res 

219 

220 _dt = self.dt if dt == -1 else dt 

221 desc['level_params']['dt'] = _dt 

222 

223 desc['problem_params']['Lx'] = float(self.gamma) 

224 desc['problem_params']['Ly'] = float(self.gamma) 

225 desc['problem_params']['Lz'] = 1.0 

226 return desc 

227 

228 def get_initial_condition(self, P, *args, restart_idx=0, **kwargs): 

229 if self.ic_config['config'] is None or restart_idx != 0: 

230 return super().get_initial_condition(P, *args, restart_idx=restart_idx, **kwargs) 

231 

232 # read initial conditions 

233 from pySDC.helpers.fieldsIO import FieldsIO 

234 

235 ic_config = self.ic_config['config']( 

236 args={**self.args, 'res': self.ic_config['res'], 'dt': self.ic_config['dt']} 

237 ) 

238 desc = ic_config.get_description(res=self.ic_config['res'], dt=self.ic_config['dt']) 

239 ic_nx = desc['problem_params']['nx'] 

240 ic_ny = desc['problem_params']['ny'] 

241 ic_nz = desc['problem_params']['nz'] 

242 

243 _P = type(P)(nx=ic_nx, ny=ic_ny, nz=ic_nz, comm=P.comm, useGPU=P.useGPU) 

244 _P.setUpFieldsIO() 

245 filename = ic_config.get_file_name() 

246 ic_file = FieldsIO.fromFile(filename) 

247 t0, ics = ic_file.readField(-1) 

248 ics = ics[: P.spectral.ncomponents, ...] 

249 P.logger.info(f'Loaded initial conditions from {filename!r} at t={t0}.') 

250 

251 # interpolate the initial conditions using padded transforms 

252 padding = (P.nx / ic_nx, P.ny / ic_ny, P.nz / ic_nz) 

253 P.logger.info(f'Interpolating initial conditions from {ic_nx}x{ic_ny}x{ic_nz} to {P.nx}x{P.ny}x{P.nz}') 

254 

255 ics = _P.xp.array(ics) 

256 _ics_hat = _P.transform(ics) 

257 ics_interpolated = _P.itransform(_ics_hat, padding=padding) 

258 

259 self.get_LogToFile() 

260 

261 P.setUpFieldsIO() 

262 if P.spectral_space: 

263 u0_hat = P.u_init_forward 

264 u0_hat[...] = P.transform(ics_interpolated) 

265 return u0_hat, 0 

266 else: 

267 return ics_interpolated, 0 

268 

269 

270class RBC3DM2K3(RBC3Dverification): 

271 

272 def get_description(self, *args, **kwargs): 

273 desc = super().get_description(*args, **kwargs) 

274 desc['level_params']['nsweeps'] = 3 

275 desc['sweeper_params']['num_nodes'] = 2 

276 return desc 

277 

278 

279class RBC3DM2K2(RBC3Dverification): 

280 

281 def get_description(self, *args, **kwargs): 

282 desc = super().get_description(*args, **kwargs) 

283 desc['level_params']['nsweeps'] = 2 

284 desc['sweeper_params']['num_nodes'] = 2 

285 return desc 

286 

287 

288class RBC3DM3K4(RBC3Dverification): 

289 

290 def get_description(self, *args, **kwargs): 

291 desc = super().get_description(*args, **kwargs) 

292 desc['level_params']['nsweeps'] = 4 

293 desc['sweeper_params']['num_nodes'] = 3 

294 return desc 

295 

296 

297class RBC3DM4K4(RBC3Dverification): 

298 

299 def get_description(self, *args, **kwargs): 

300 desc = super().get_description(*args, **kwargs) 

301 desc['level_params']['nsweeps'] = 4 

302 desc['sweeper_params']['num_nodes'] = 4 

303 return desc 

304 

305 

306class RBC3DverificationRK(RBC3Dverification): 

307 

308 def get_description(self, *args, res=-1, dt=-1, **kwargs): 

309 from pySDC.implementations.sweeper_classes.Runge_Kutta import ARK3 

310 

311 desc = super().get_description(*args, res=res, dt=dt, **kwargs) 

312 desc['level_params']['nsweeps'] = 1 

313 desc['level_params']['restol'] = -1 

314 desc['step_params']['maxiter'] = 1 

315 desc['sweeper_params']['skip_residual_computation'] = ('IT_CHECK', 'IT_DOWN', 'IT_UP', 'IT_FINE', 'IT_COARSE') 

316 desc['sweeper_params'].pop('QI') 

317 desc['sweeper_params'].pop('num_nodes') 

318 desc['sweeper_class'] = ARK3 

319 return desc 

320 

321 

322class RBC3DverificationEuler(RBC3DverificationRK): 

323 

324 def get_description(self, *args, res=-1, dt=-1, **kwargs): 

325 from pySDC.implementations.sweeper_classes.Runge_Kutta import IMEXEulerStifflyAccurate 

326 

327 desc = super().get_description(*args, res=res, dt=dt, **kwargs) 

328 desc['sweeper_class'] = IMEXEulerStifflyAccurate 

329 return desc 

330 

331 

332# --- Ra 1e5 --- 

333class RBC3DG4R4SDC22Ra1e5(RBC3DM2K2): 

334 Tend = 200 

335 dt = 4e-2 

336 res = 32 

337 converged = 50 

338 

339 

340class RBC3DG4R4SDC23Ra1e5(RBC3DM2K3): 

341 Tend = 200 

342 dt = 6e-2 

343 res = 32 

344 converged = 50 

345 

346 

347class RBC3DG4R4SDC34Ra1e5(RBC3DM3K4): 

348 Tend = 200 

349 dt = 6e-2 

350 res = 32 

351 converged = 50 

352 

353 

354class RBC3DG4R4SDC44Ra1e5(RBC3DM4K4): 

355 Tend = 200 

356 dt = 6e-2 

357 res = 32 

358 converged = 50 

359 

360 

361class RBC3DG4R4RKRa1e5(RBC3DverificationRK): 

362 Tend = 200 

363 dt = 5e-2 

364 res = 32 

365 converged = 50 

366 

367 

368class RBC3DG4R4EulerRa1e5(RBC3DverificationEuler): 

369 Tend = 200 

370 dt = 2e-2 

371 res = 32 

372 converged = 50 

373 

374 

375# --- Ra 1e6 --- 

376class RBC3DG4R4SDC44Ra1e6(RBC3DM4K4): 

377 Tend = 75 

378 dt = 1e-2 

379 res = 64 

380 converged = 22 

381 ic_config = {'config': RBC3DG4R4SDC34Ra1e5, 'res': 32, 'dt': 0.02} 

382 

383 

384class RBC3DG4R4SDC23Ra1e6(RBC3DM2K3): 

385 Tend = 75 

386 dt = 1e-2 

387 res = 64 

388 converged = 22 

389 ic_config = {'config': RBC3DG4R4SDC34Ra1e5, 'res': 32, 'dt': 0.02} 

390 

391 

392class RBC3DG4R4RKRa1e6(RBC3DverificationRK): 

393 Tend = 75 

394 dt = 1e-2 

395 res = 64 

396 ic_config = {'config': RBC3DG4R4SDC34Ra1e5, 'res': 32, 'dt': 0.02} 

397 converged = 22 

398 

399 

400class RBC3DG4R4EulerRa1e6(RBC3DverificationEuler): 

401 Tend = 75 

402 dt = 5e-3 

403 res = 64 

404 ic_config = {'config': RBC3DG4R4SDC34Ra1e5, 'res': 32, 'dt': 0.02} 

405 converged = 22 

406 

407 

408# --- Ra 1e7 --- 

409class RBC3DG4R4SDC23Ra1e7(RBC3DM2K3): 

410 Tend = 45 

411 dt = 5e-3 

412 res = 128 

413 converged = 25 

414 ic_config = {'config': RBC3DG4R4SDC23Ra1e6, 'res': 64, 'dt': 0.01} 

415 

416 

417class RBC3DG4R4SDC44Ra1e7(RBC3DM4K4): 

418 Tend = 45 

419 dt = 5e-3 

420 res = 128 

421 converged = 25 

422 ic_config = {'config': RBC3DG4R4SDC23Ra1e6, 'res': 64, 'dt': 0.01} 

423 

424 

425class RBC3DG4R4EulerRa1e7(RBC3DverificationEuler): 

426 Tend = 45 

427 dt = 1e-3 

428 res = 128 

429 converged = 25 

430 ic_config = {'config': RBC3DG4R4SDC23Ra1e6, 'res': 64, 'dt': 0.01} 

431 

432 

433class RBC3DG4R4RKRa1e7(RBC3DverificationRK): 

434 Tend = 45 

435 dt = 4e-3 

436 res = 128 

437 converged = 25 

438 ic_config = {'config': RBC3DG4R4SDC23Ra1e6, 'res': 64, 'dt': 0.01}