HyPar  1.0
Finite-Difference Hyperbolic-Parabolic PDE Solver on Cartesian Grids
NavierStokes3DPreStep.c File Reference

Pre-step function for 3D Navier Stokes equations. More...

#include <arrayfunctions.h>
#include <mathfunctions.h>
#include <matmult_native.h>
#include <physicalmodels/navierstokes3d.h>
#include <hypar.h>

Go to the source code of this file.

Functions

int NavierStokes3DPreStep (double *u, void *s, void *m, double waqt)
 

Detailed Description

Pre-step function for 3D Navier Stokes equations.

Author
Debojyoti Ghosh

Definition in file NavierStokes3DPreStep.c.

Function Documentation

◆ NavierStokes3DPreStep()

int NavierStokes3DPreStep ( double *  u,
void *  s,
void *  m,
double  waqt 
)

Pre-step function for the 3D Navier Stokes equations: This function is called at the beginning of each time step.

  • The solution at the beginning of the step is copied into NavierStokes3D::solution for the linearized flux partitioning.
  • The linearized "fast" Jacobian (representing the acoustic modes) NavierStokes3D::fast_jac is computed as:

    \begin{equation} A_f\left({\bf u}\right) = R\left({\bf u}\right)\Lambda_f\left({\bf u}\right)L\left({\bf u}\right) \end{equation}

    where \(R\) and \(L\) are the matrices of right and left eigenvectors, and,

    \begin{equation} \Lambda_f = diag\left[0,0,0,u+c,u-c \right] \end{equation}

    with \(c\) being the speed of sound.

    Reference:
  • Ghosh, D., Constantinescu, E. M., Semi-Implicit Time Integration of Atmospheric Flows with Characteristic-Based Flux Partitioning, SIAM Journal on Scientific Computing, 38 (3), 2016, A1848-A1875, http://dx.doi.org/10.1137/15M1044369.
Parameters
uSolution vector
sSolver object of type HyPar
mMPI object of type MPIVariables
waqtCurrent simulation time

Definition at line 32 of file NavierStokes3DPreStep.c.

38 {
39  HyPar *solver = (HyPar*) s;
40  NavierStokes3D *param = (NavierStokes3D*) solver->physics;
41  int *dim = solver->dim_local;
42  int ghosts = solver->ghosts, dir, p;
43  double *A;
44 
45  static const int ndims = _MODEL_NDIMS_;
46  static const int JacSize = _MODEL_NVARS_*_MODEL_NVARS_;
47  static int index[_MODEL_NDIMS_], bounds[_MODEL_NDIMS_], offset[_MODEL_NDIMS_];
50 
51  /* set bounds for array index to include ghost points */
52  _ArrayAddCopy1D_(dim,(2*ghosts),bounds,ndims);
53  /* set offset such that index is compatible with ghost point arrangement */
54  _ArraySetValue_(offset,ndims,-ghosts);
55  /* copy the solution to act as a reference for linearization */
57 
58  int done = 0; _ArraySetValue_(index,ndims,0);
59  while (!done) {
60  _ArrayIndex1DWO_(ndims,dim,index,offset,ghosts,p);
61  int q = _MODEL_NVARS_*p;
62 
63  dir = _XDIR_;
64  A = (param->fast_jac + _MODEL_NDIMS_*JacSize*p + dir*JacSize);
65  /* get the eigenvalues, and left & right eigenvectors */
69  /* remove the entropy modes (corresponding to eigenvalues u) */
70  D[0] = D[12] = D[18] = 0.0;
71  /* assemble the Jacobian */
72  MatMult5(_MODEL_NVARS_,DL,D,L );
73  MatMult5(_MODEL_NVARS_,A ,R,DL);
74 
75  dir = _YDIR_;
76  A = (param->fast_jac + _MODEL_NDIMS_*JacSize*p + dir*JacSize);
77  /* get the eigenvalues, and left & right eigenvectors */
81  /* remove the entropy modes (corresponding to eigenvalues v) */
82  D[0] = D[6] = D[18] = 0.0;
83  /* assemble the Jacobian */
84  MatMult5(_MODEL_NVARS_,DL,D,L );
85  MatMult5(_MODEL_NVARS_,A ,R,DL);
86 
87  dir = _ZDIR_;
88  A = (param->fast_jac + _MODEL_NDIMS_*JacSize*p + dir*JacSize);
89  /* get the eigenvalues, and left & right eigenvectors */
93  /* remove the entropy modes (corresponding to eigenvalues v) */
94  D[0] = D[6] = D[12] = 0.0;
95  /* assemble the Jacobian */
96  MatMult5(_MODEL_NVARS_,DL,D,L );
97  MatMult5(_MODEL_NVARS_,A ,R,DL);
98 
99  _ArrayIncrementIndex_(ndims,bounds,index,done);
100  }
101  return(0);
102 }
#define _NavierStokes3DRightEigenvectors_(u, stride, R, ga, dir)
#define _NavierStokes3DLeftEigenvectors_(u, stride, L, ga, dir)
#define _ZDIR_
Structure containing variables and parameters specific to the 3D Navier Stokes equations. This structure contains the physical parameters, variables, and function pointers specific to the 3D Navier-Stokes equations.
#define _ArrayAddCopy1D_(x, a, y, size)
Structure containing all solver-specific variables and functions.
Definition: hypar.h:23
#define _YDIR_
Definition: euler2d.h:41
#define _MODEL_NDIMS_
Definition: euler1d.h:56
#define _NavierStokes3DEigenvalues_(u, stride, D, gamma, dir)
#define _ArrayIndex1DWO_(N, imax, i, offset, ghost, index)
static const int _NavierStokes3D_stride_
#define _ArraySetValue_(x, size, value)
int * dim_local
Definition: hypar.h:37
#define _ArrayIncrementIndex_(N, imax, i, done)
void * physics
Definition: hypar.h:266
#define MatMult5(N, A, X, Y)
#define _XDIR_
Definition: euler1d.h:75
int ghosts
Definition: hypar.h:52
#define _MODEL_NVARS_
Definition: euler1d.h:58
int npoints_local_wghosts
Definition: hypar.h:42
#define _ArrayCopy1D_(x, y, size)