HyPar  1.0
Finite-Difference Hyperbolic-Parabolic PDE Solver on Cartesian Grids
PetscPostStage.cpp
Go to the documentation of this file.
1 
5 #ifdef with_petsc
6 
7 #include <stdio.h>
8 #include <basic.h>
9 #include <simulation_object.h>
10 #include <petscinterface.h>
11 
12 #undef __FUNCT__
13 #define __FUNCT__ "PetscPostStage"
14 
16 PetscErrorCode PetscPostStage( TS ts,
17  PetscReal stagetime,
18  PetscInt stageindex,
19  Vec *Y )
21 {
22  PETScContext* context(nullptr);
23 
24  PetscFunctionBegin;
25 
26  TSGetApplicationContext(ts,&context);
27  if (!context) {
28  fprintf(stderr,"Error in PetscPreTimeStep: Null context!\n");
29  return(1);
30  }
31  SimulationObject* sim = (SimulationObject*) context->simobj;
32  int nsims = context->nsims;
33 
34  TSType time_scheme;
35  TSGetType(ts,&time_scheme);
36 
37  TSGetTimeStep(ts,&(context->dt));
38  context->stage_index = stageindex;
39  if (context->stage_times.size() == stageindex) {
40  context->stage_times.push_back(stagetime/context->dt);
41  }
42 
43  for (int ns = 0; ns < nsims; ns++) {
44 
45  HyPar* solver = &(sim[ns].solver);
46  MPIVariables* mpi = &(sim[ns].mpi);
47 
48  /* get solution */
49  TransferVecFromPETSc(solver->u,Y[stageindex],context,ns,context->offsets[ns]);
50 
51  /* apply immersed boundaries */
52  solver->ApplyIBConditions(solver,mpi,solver->u,stagetime);
53 
54  /* If using a non-linear scheme with ARKIMEX methods,
55  compute the non-linear finite-difference operator */
56  if (!strcmp(time_scheme,TSARKIMEX)) {
57  solver->NonlinearInterp( solver->u,
58  solver,
59  mpi,
60  (double)stagetime,
61  solver->FFunction );
62  }
63 
64  /* Call any physics-specific post-stage function, if available */
65  if (solver->PostStage) {
66  solver->PostStage(solver->u,solver,mpi,stagetime);
67  }
68 
69  TransferVecToPETSc(solver->u,Y[stageindex],context,ns,context->offsets[ns]);
70 
71  }
72 
73  PetscFunctionReturn(0);
74 }
75 
76 #endif
int TransferVecToPETSc(const double *const, Vec, void *, const int, const int)
Structure defining a simulation.
std::vector< double > stage_times
int(* NonlinearInterp)(double *, void *, void *, double, int(*)(double *, double *, int, void *, double))
Definition: hypar.h:228
double * u
Definition: hypar.h:116
Some basic definitions and macros.
Simulation object.
int(* ApplyIBConditions)(void *, void *, double *, double)
Definition: hypar.h:217
int(* FFunction)(double *, double *, int, void *, double)
Definition: hypar.h:276
Structure containing all solver-specific variables and functions.
Definition: hypar.h:23
PetscErrorCode PetscPostStage(TS ts, PetscReal stagetime, PetscInt stageindex, Vec *Y)
Structure containing the variables for time-integration with PETSc.
int(* PostStage)(double *, void *, void *, double)
Definition: hypar.h:336
Structure of MPI-related variables.
int TransferVecFromPETSc(double *const, const Vec, void *, const int, const int)