HyPar  1.0
Finite-Difference Hyperbolic-Parabolic PDE Solver on Cartesian Grids
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
ParabolicFunctionNC2Stage.c File Reference

Evaluate the parabolic term using a 2-stage finite difference discretization. More...

#include <stdlib.h>
#include <basic.h>
#include <arrayfunctions.h>
#include <mpivars.h>
#include <hypar.h>

Go to the source code of this file.

Functions

int ParabolicFunctionNC2Stage (double *par, double *u, void *s, void *m, double t)
 

Detailed Description

Evaluate the parabolic term using a 2-stage finite difference discretization.

Author
Debojyoti Ghosh

Definition in file ParabolicFunctionNC2Stage.c.

Function Documentation

int ParabolicFunctionNC2Stage ( double *  par,
double *  u,
void *  s,
void *  m,
double  t 
)

Evaluate the parabolic term using a "1.5"-stage finite-difference spatial discretization: The parabolic term is assumed to be of the form:

\begin{equation} {\bf P}\left({\bf u}\right) = \sum_{d1=0}^{D-1}\sum_{d2=0}^{D-1} \frac {\partial^2 h_{d1,d2}\left(\bf u\right)} {\partial x_{d1} \partial x_{d2}}, \end{equation}

where \(d1\) and \(d2\) are spatial dimension indices, and \(D\) is the total number of spatial dimensions (HyPar::ndims). This term is discretized at a grid point as:

\begin{equation} \left.{\bf P}\left({\bf u}\right)\right|_j = \sum_{d1=0}^{D-1} \sum_{d2=0}^{D-1} \frac { \mathcal{D}_{d1}\mathcal{D}_{d2} \left[ {\bf h}_{d1,d2} \right] } {\Delta x_{d1} \Delta x_{d2}}, \end{equation}

where \(\mathcal{D}\) denotes the finite-difference approximation to the first derivative. Each of the first derivative approximations are \(\mathcal{D}_{d1}\) and \(\mathcal{D}_{d2}\) are computed separately, and thus the cross-derivative is evaluated in two steps using HyPar::FirstDerivativePar.

Notes:

  • This form of the parabolic term does allow for cross-derivatives ( \( d1 \ne d2 \)).
  • A \(n\)-th order central approximation to the second derivative can be expressed as a conjugation of two \((n-1)\)-th order approximations to the first derivative, one forward and one backward. Computing it this way avoids odd-even decoupling. Thus, where possible HyPar::FirstDerivativePar should point to the function computing \((n-1)\)-th order first derivative where \(n\) is the desired order. Currently, this is implemented only for \(n=2\). For other values of \(n\), the first derivative is also computed with a \(n\)-th order approximation.

To use this form of the parabolic term:

See Also
ParabolicFunctionNC1_5Stage()
Parameters
pararray to hold the computed parabolic term
usolution
sSolver object of type HyPar
mMPI object of type MPIVariables
tCurrent simulation time

Definition at line 40 of file ParabolicFunctionNC2Stage.c.

47 {
48  HyPar *solver = (HyPar*) s;
49  MPIVariables *mpi = (MPIVariables*) m;
50  double *Func = solver->fluxC;
51  double *Deriv1 = solver->Deriv1;
52  double *Deriv2 = solver->Deriv2;
53  int d, d1, d2, v, p, done;
54  double dxinv1, dxinv2;
56 
57  int ndims = solver->ndims;
58  int nvars = solver->nvars;
59  int ghosts = solver->ghosts;
60  int *dim = solver->dim_local;
61  double *dxinv = solver->dxinv;
62  int size = solver->npoints_local_wghosts;
63 
64  printf("HFunction is defined? = %p\n", solver->HFunction);
65 
66  if (!solver->HFunction) return(0); /* zero parabolic terms */
67  solver->count_par++;
68 
69  int index[ndims];
70  _ArraySetValue_(par,size*nvars,0.0);
71 
72  for (d1 = 0; d1 < ndims; d1++) {
73  for (d2 = 0; d2 < ndims; d2++) {
74 
75  /* calculate the diffusion function */
76  IERR solver->HFunction(Func,u,d1,d2,solver,t); CHECKERR(ierr);
77  IERR solver->FirstDerivativePar(Deriv1,Func ,d1, 1,solver,mpi); CHECKERR(ierr);
78  IERR MPIExchangeBoundariesnD(ndims,nvars,dim,ghosts,mpi,Deriv1); CHECKERR(ierr);
79  IERR solver->FirstDerivativePar(Deriv2,Deriv1,d2,-1,solver,mpi); CHECKERR(ierr);
80 
81  /* calculate the final term - second derivative of the diffusion function */
82  done = 0; _ArraySetValue_(index,ndims,0);
83  while (!done) {
84  _ArrayIndex1D_(ndims,dim,index,ghosts,p);
85  _GetCoordinate_(d1,index[d1],dim,ghosts,dxinv,dxinv1);
86  _GetCoordinate_(d2,index[d2],dim,ghosts,dxinv,dxinv2);
87  for (v=0; v<nvars; v++) par[nvars*p+v] += (dxinv1*dxinv2 * Deriv2[nvars*p+v]);
88  _ArrayIncrementIndex_(ndims,dim,index,done);
89  }
90 
91  }
92  }
93 
94  if (solver->flag_ib) _ArrayBlockMultiply_(par,solver->iblank,size,nvars);
95  return(0);
96 }
int npoints_local_wghosts
Definition: hypar.h:42
#define _ArraySetValue_(x, size, value)
double * iblank
Definition: hypar.h:436
#define _ArrayIncrementIndex_(N, imax, i, done)
double * Deriv1
Definition: hypar.h:151
int(* HFunction)(double *, double *, int, int, void *, double)
Definition: hypar.h:313
double * Deriv2
Definition: hypar.h:151
int flag_ib
Definition: hypar.h:441
int * dim_local
Definition: hypar.h:37
#define _GetCoordinate_(dir, i, dim, ghosts, x, coord)
Definition: basic.h:31
int ghosts
Definition: hypar.h:52
#define _ArrayIndex1D_(N, imax, i, ghost, index)
int count_par
Definition: hypar.h:418
#define _ArrayBlockMultiply_(x, a, n, bs)
int(* FirstDerivativePar)(double *, double *, int, int, void *, void *)
Definition: hypar.h:243
double * fluxC
Definition: hypar.h:128
int nvars
Definition: hypar.h:29
#define CHECKERR(ierr)
Definition: basic.h:18
double * dxinv
Definition: hypar.h:110
int ndims
Definition: hypar.h:26
#define IERR
Definition: basic.h:16
int MPIExchangeBoundariesnD(int, int, int *, int, void *, double *)
Structure of MPI-related variables.
Structure containing all solver-specific variables and functions.
Definition: hypar.h:23
#define _DECLARE_IERR_
Definition: basic.h:17