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

Evaluate the parabolic term using a conservative spatial 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 ParabolicFunctionCons1Stage (double *par, double *u, void *s, void *m, double t)
 

Detailed Description

Evaluate the parabolic term using a conservative spatial discretization.

Author
Debojyoti Ghosh

Definition in file ParabolicFunctionCons1Stage.c.

Function Documentation

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

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

\begin{equation} {\bf P}\left({\bf u}\right) = \sum_{d=0}^{D-1} \frac {\partial^2 {\bf g}_d\left(\bf u\right)} {\partial x_d^2}, \end{equation}

which is discretized at a grid point as:

\begin{equation} \left.{\bf P}\left({\bf u}\right)\right|_j = \sum_{d=0}^{D-1} \frac { \hat{\bf g}_{j+1/2} - \hat{\bf g}_{j-1/2} } {\Delta x_d^2}, \end{equation}

where \(d\) is the spatial dimension index, \(D\) is the total number of spatial dimensions (HyPar::ndims), and \(j\) is the grid index along \(d\). \(\hat{\bf g}_d\) is the numerical approximation to the second primitive of \({\bf g}_d\left({\bf u}\right)\), computed using HyPar::InterpolateInterfacesPar.

Note: this form of the parabolic term does not allow for cross-derivatives.

To use this form of the parabolic term:

Reference: Liu, Y., Shu, C.-W., Zhang, M., "High order finite difference WENO schemes for nonlinear degenerate parabolic equations", SIAM J. Sci. Comput., 33 (2), 2011, pp. 939-965, http://dx.doi.org/10.1137/100791002

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 34 of file ParabolicFunctionCons1Stage.c.

41 {
42  HyPar *solver = (HyPar*) s;
43  MPIVariables *mpi = (MPIVariables*) m;
44  double *FluxI = solver->fluxI; /* interface flux array */
45  double *Func = solver->fluxC; /* diffusion function array */
46  int d, v, i, done;
48 
49  int ndims = solver->ndims;
50  int nvars = solver->nvars;
51  int ghosts = solver->ghosts;
52  int *dim = solver->dim_local;
53  double *dxinv = solver->dxinv;
54  int size = solver->npoints_local_wghosts;
55 
56  int index[ndims], index1[ndims], index2[ndims], dim_interface[ndims];
57 
58  _ArraySetValue_(par,size*nvars,0.0);
59  if (!solver->GFunction) return(0); /* zero parabolic term */
60  solver->count_par++;
61 
62  int offset = 0;
63  for (d = 0; d < ndims; d++) {
64  _ArrayCopy1D_(dim,dim_interface,ndims); dim_interface[d]++;
65  int size_cellcenter = 1; for (i = 0; i < ndims; i++) size_cellcenter *= (dim[i] + 2*ghosts);
66  int size_interface = 1; for (i = 0; i < ndims; i++) size_interface *= dim_interface[i];
67 
68  /* evaluate cell-centered flux */
69  IERR solver->GFunction(Func,u,d,solver,t); CHECKERR(ierr);
70  /* compute interface fluxes */
71  IERR solver->InterpolateInterfacesPar(FluxI,Func,d,solver,mpi); CHECKERR(ierr);
72 
73  /* calculate the second derivative */
74  done = 0; _ArraySetValue_(index,ndims,0);
75  int p, p1, p2;
76  while (!done) {
77  _ArrayCopy1D_(index,index1,ndims);
78  _ArrayCopy1D_(index,index2,ndims); index2[d]++;
79  _ArrayIndex1D_(ndims,dim ,index ,ghosts,p);
80  _ArrayIndex1D_(ndims,dim_interface,index1,0 ,p1);
81  _ArrayIndex1D_(ndims,dim_interface,index2,0 ,p2);
82  for (v=0; v<nvars; v++)
83  par[nvars*p+v] += ((dxinv[offset+ghosts+index[d]] * dxinv[offset+ghosts+index[d]])
84  * (FluxI[nvars*p2+v] - FluxI[nvars*p1+v]));
85  _ArrayIncrementIndex_(ndims,dim,index,done);
86  }
87 
88  offset += dim[d] + 2*ghosts;
89  }
90 
91  if (solver->flag_ib) _ArrayBlockMultiply_(par,solver->iblank,size,nvars);
92  return(0);
93 }
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)
int flag_ib
Definition: hypar.h:441
int * dim_local
Definition: hypar.h:37
int ghosts
Definition: hypar.h:52
#define _ArrayIndex1D_(N, imax, i, ghost, index)
int count_par
Definition: hypar.h:418
double * fluxI
Definition: hypar.h:136
int(* InterpolateInterfacesPar)(double *, double *, int, void *, void *)
Definition: hypar.h:239
#define _ArrayBlockMultiply_(x, a, n, bs)
#define _ArrayCopy1D_(x, y, size)
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
Structure of MPI-related variables.
Structure containing all solver-specific variables and functions.
Definition: hypar.h:23
int(* GFunction)(double *, double *, int, void *, double)
Definition: hypar.h:310
#define _DECLARE_IERR_
Definition: basic.h:17