HyPar  1.0
Finite-Difference Hyperbolic-Parabolic PDE Solver on Cartesian Grids
main.cpp File Reference
#include <stdio.h>
#include <stdlib.h>
#include <sys/time.h>
#include <string>
#include <petscinterface.h>
#include <Python.h>
#include <numpy/arrayobject.h>
#include <mpivars_cpp.h>
#include <simulation_library.h>

Go to the source code of this file.

Macros

#define NPY_NO_DEPRECATED_API   NPY_1_7_API_VERSION
 

Functions

static void initializePython (int)
 
static void initializePythonPlotting (int)
 
int main (int argc, char **argv)
 Main driver. More...
 

Variables

static const char help [] = "HyPar - A finite-difference algorithm for solving hyperbolic-parabolic PDEs"
 

Macro Definition Documentation

◆ NPY_NO_DEPRECATED_API

#define NPY_NO_DEPRECATED_API   NPY_1_7_API_VERSION

Definition at line 161 of file main.cpp.

Function Documentation

◆ initializePython()

void initializePython ( int  a_rank)
static

Initialize Python interpreter

Parameters
a_rankMPI rank

Definition at line 441 of file main.cpp.

442 {
443  Py_Initialize();
444  if (!a_rank) printf("Initialized Python.\n");
445  PyRun_SimpleString("import os");
446  PyRun_SimpleString("hypar_dir = os.environ.get('HYPAR_DIR')");
447 #ifndef serial
448  MPI_Barrier(MPI_COMM_WORLD);
449 #endif
450  return;
451 }

◆ initializePythonPlotting()

void initializePythonPlotting ( int  a_rank)
static

Initialize Python plotting stuff

Parameters
a_rankMPI rank

Definition at line 454 of file main.cpp.

455 {
456  /* load plotting tools and scripts */
457  PyRun_SimpleString("hypar_dir_plt_py = hypar_dir + '/src/PlottingFunctions'");
458  PyRun_SimpleString("import sys");
459  PyRun_SimpleString("sys.path.append(hypar_dir_plt_py)");
460  if (!a_rank) {
461  PyRun_SimpleString
462  ("print('Added plotting script directory (%s) to Python path.' % hypar_dir_plt_py)");
463  }
464 #ifndef serial
465  MPI_Barrier(MPI_COMM_WORLD);
466 #endif
467  return;
468 }

◆ main()

int main ( int  argc,
char **  argv 
)

Main driver.

The main driver function that calls the initialization, solving, and cleaning up functions.

Definition at line 181 of file main.cpp.

182 {
183  int ierr = 0, d, n;
184  struct timeval main_start, solve_start;
185  struct timeval main_end , solve_end ;
186 #ifdef with_petsc
187  PetscBool use_petscts;
188 #endif
189  int use_petsc = 0;
190 
191 #ifdef serial
192  int world = 0;
193  int rank = 0;
194  int nproc = 1;
195  printf("HyPar - Serial Version\n");
196 #else
197  MPI_Comm world;
198  int rank, nproc;
199  MPI_Init(&argc,&argv);
200  MPI_Comm_dup(MPI_COMM_WORLD, &world);
201  MPI_Comm_rank(MPI_COMM_WORLD,&rank );
202  MPI_Comm_size(MPI_COMM_WORLD,&nproc);
203  if (!rank) printf("HyPar - Parallel (MPI) version with %d processes\n",nproc);
204 #endif
205 
206 #ifdef with_petsc
207  PetscInitialize(&argc,&argv,(char*)0,help);
208  if (!rank) printf("Compiled with PETSc time integration.\n");
209 #endif
210 
211 #ifdef with_python
212  initializePython(rank);
214 #endif
215 
216  gettimeofday(&main_start,NULL);
217 
218  int sim_type = -1;
219  Simulation *sim = NULL;
220 
221  if (!rank) {
222 
223  std::string ensemble_sim_fname(_ENSEMBLE_SIM_INP_FNAME_);
224  std::string sparsegrids_sim_fname(_SPARSEGRIDS_SIM_INP_FNAME_);
225 
226  FILE *f_ensemble_sim = fopen(ensemble_sim_fname.c_str(), "r");
227  FILE *f_sparsegrids_sim = fopen(sparsegrids_sim_fname.c_str(), "r");
228 
229  if (f_ensemble_sim && f_sparsegrids_sim) {
230 
231  fprintf(stderr,"Error: Cannot have both %s and %s input files.\n",
233  fprintf(stderr, "Remove one or both of them depending on the kind of simulation you want to run.\n");
234  fclose(f_ensemble_sim);
235  fclose(f_sparsegrids_sim);
236 
237  } else if (f_ensemble_sim) {
238 
239  sim_type = _SIM_TYPE_ENSEMBLE_;
240  fclose(f_ensemble_sim);
241 
242  } else if (f_sparsegrids_sim) {
243 
244  sim_type = _SIM_TYPE_SPARSE_GRIDS_;
245  fclose(f_sparsegrids_sim);
246 
247  } else {
248 
249  sim_type = _SIM_TYPE_SINGLE_;
250 
251  }
252 
253  }
254 
255 #ifndef serial
256  MPI_Bcast(&sim_type, 1, MPI_INT, 0, MPI_COMM_WORLD);
257 #endif
258 
259  if (sim_type == _SIM_TYPE_SINGLE_) {
260  sim = new SingleSimulation;
261  } else if (sim_type == _SIM_TYPE_ENSEMBLE_) {
262  if (!rank) printf("-- Ensemble Simulation --\n");
263  sim = new EnsembleSimulation;
264  } else if (sim_type == _SIM_TYPE_SPARSE_GRIDS_) {
265  if (!rank) printf("-- Sparse Grids Simulation --\n");
266  sim = new SparseGridsSimulation;
267  } else {
268  fprintf(stderr, "ERROR: invalid sim_type (%d) on rank %d.\n",
269  sim_type, rank);
270  }
271 
272  if (sim == NULL) {
273  fprintf(stderr, "ERROR: unable to create sim on rank %d.\n",
274  rank );
275  return 1;
276  }
277 
278  /* Allocate simulation objects */
279  ierr = sim->define(rank, nproc);
280  if (!sim->isDefined()) {
281  printf("Error: Simulation::define() failed on rank %d\n",
282  rank);
283  return 1;
284  }
285  if (ierr) {
286  printf("Error: Simulation::define() returned with status %d on process %d.\n",
287  ierr, rank);
288  return(ierr);
289  }
290 
291 #ifndef serial
292  ierr = sim->mpiCommDup();
293 #endif
294 
295 #ifdef with_petsc
296  use_petscts = PETSC_FALSE; /* default value */
297  ierr = PetscOptionsGetBool( nullptr,nullptr,
298  "-use-petscts",
299  &use_petscts,
300  nullptr); CHKERRQ(ierr);
301  if (use_petscts == PETSC_TRUE) use_petsc = 1;
302  sim->usePetscTS(use_petscts);
303 #endif
304 
305  /* Read Inputs */
306  ierr = sim->ReadInputs();
307  if (ierr) {
308  printf("Error: Simulation::ReadInputs() returned with status %d on process %d.\n",ierr,rank);
309  return(ierr);
310  }
311 
312  /* Initialize and allocate arrays */
313  ierr = sim->Initialize();
314  if (ierr) {
315  printf("Error: Simulation::Initialize() returned with status %d on process %d.\n",ierr,rank);
316  return(ierr);
317  }
318 
319  /* read and set grid & initial solution */
320  ierr = sim->InitialSolution();
321  if (ierr) {
322  printf("Error: Simulation::InitialSolution() returned with status %d on process %d.\n",ierr,rank);
323  return(ierr);
324  }
325 
326  /* Initialize domain boundaries */
327  ierr = sim->InitializeBoundaries();
328  if (ierr) {
329  printf("Error: Simulation::InitializeBoundaries() returned with status %d on process %d.\n",ierr,rank);
330  return(ierr);
331  }
332 
333  /* Initialize immersed boundaries */
334  ierr = sim->InitializeImmersedBoundaries();
335  if (ierr) {
336  printf("Error: Simulation::InitializeImmersedBoundaries() returned with status %d on process %d.\n",ierr,rank);
337  return(ierr);
338  }
339 
340  /* Initialize solvers */
341  ierr = sim->InitializeSolvers();
342  if (ierr) {
343  printf("Error: Simulation::InitializeSolvers() returned with status %d on process %d.\n",ierr,rank);
344  return(ierr);
345  }
346 
347  /* Initialize physics */
348  ierr = sim->InitializePhysics();
349  if (ierr) {
350  printf("Error: Simulation::InitializePhysics() returned with status %d on process %d.\n",ierr,rank);
351  return(ierr);
352  }
353 
354  /* Initialize physics data */
355  ierr = sim->InitializePhysicsData();
356  if (ierr) {
357  printf("Error: Simulation::InitializePhysicsData() returned with status %d on process %d.\n",ierr,rank);
358  return(ierr);
359  }
360 
361  /* Wrap up initializations */
362  ierr = sim->InitializationWrapup();
363  if (ierr) {
364  printf("Error: Simulation::InitializationWrapup() returned with status %d on process %d.\n",ierr,rank);
365  return(ierr);
366  }
367 
368  /* Initializations complete */
369 
370  /* Run the solver */
371 #ifndef serial
372  MPI_Barrier(MPI_COMM_WORLD);
373 #endif
374  gettimeofday(&solve_start,NULL);
375 #ifdef with_petsc
376  if (use_petsc == 1) {
377  /* Use PETSc time-integration */
378  ierr = sim->SolvePETSc();
379  if (ierr) {
380  printf("Error: Simulation::SolvePETSc() returned with status %d on process %d.\n",ierr,rank);
381  return(ierr);
382  }
383  } else {
384  /* Use native time-integration */
385  ierr = sim->Solve();
386  if (ierr) {
387  printf("Error: Simulation::Solve() returned with status %d on process %d.\n",ierr,rank);
388  return(ierr);
389  }
390  }
391 #else
392  /* Use native time-integration */
393  ierr = sim->Solve();
394  if (ierr) {
395  printf("Error: Simulation::Solve() returned with status %d on process %d.\n",ierr,rank);
396  return(ierr);
397  }
398 #endif
399  gettimeofday(&solve_end,NULL);
400 #ifndef serial
401  MPI_Barrier(MPI_COMM_WORLD);
402 #endif
403  gettimeofday(&main_end,NULL);
404 
405  /* calculate solver and total runtimes */
406  long long walltime;
407  walltime = ( (main_end.tv_sec * 1000000 + main_end.tv_usec )
408  - (main_start.tv_sec * 1000000 + main_start.tv_usec));
409  double main_runtime = (double) walltime / 1000000.0;
410  ierr = MPIMax_double(&main_runtime,&main_runtime,1,&world); if(ierr) return(ierr);
411  walltime = ( (solve_end.tv_sec * 1000000 + solve_end.tv_usec )
412  - (solve_start.tv_sec * 1000000 + solve_start.tv_usec));
413  double solver_runtime = (double) walltime / 1000000.0;
414  ierr = MPIMax_double(&solver_runtime,&solver_runtime,1,&world); if(ierr) return(ierr);
415 
416  /* Write errors and other data */
417  sim->WriteErrors(solver_runtime, main_runtime);
418 
419  /* Cleaning up */
420  delete sim;
421  if (!rank) printf("Finished.\n");
422 
423 #ifdef with_python
424  Py_Finalize();
425 #endif
426 
427 #ifdef with_petsc
428  PetscFinalize();
429 #endif
430 
431 #ifndef serial
432  MPI_Comm_free(&world);
433  MPI_Finalize();
434 #endif
435 
436  return(0);
437 }
virtual int SolvePETSc()=0
Class describing ensemble simulations.
virtual int InitialSolution()=0
virtual int define(int, int)=0
virtual int InitializationWrapup()
Definition: simulation.h:83
virtual int mpiCommDup()=0
Class describing sparse grids simulations.
virtual int Initialize()=0
virtual int InitializeSolvers()=0
static const char help[]
Definition: main.cpp:169
virtual int InitializeImmersedBoundaries()=0
#define _SIM_TYPE_SINGLE_
virtual int ReadInputs()=0
static void initializePythonPlotting(int)
Definition: main.cpp:454
virtual void WriteErrors(double, double)=0
virtual bool isDefined() const =0
virtual int InitializeBoundaries()=0
int MPIMax_double(double *, double *, int, void *)
Definition: MPIMax.c:38
virtual int Solve()=0
#define _SIM_TYPE_SPARSE_GRIDS_
virtual int InitializePhysicsData()=0
#define _SPARSEGRIDS_SIM_INP_FNAME_
static void initializePython(int)
Definition: main.cpp:441
#define _SIM_TYPE_ENSEMBLE_
Base class for a simulation.
Definition: simulation.h:46
Class describing a single simulation.
#define _ENSEMBLE_SIM_INP_FNAME_
virtual void usePetscTS(PetscBool)=0
virtual int InitializePhysics()=0

Variable Documentation

◆ help

const char help[] = "HyPar - A finite-difference algorithm for solving hyperbolic-parabolic PDEs"
static

Definition at line 169 of file main.cpp.