HyPar  1.0
Finite-Difference Hyperbolic-Parabolic PDE Solver on Cartesian Grids
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
HyperbolicFunction.c
Go to the documentation of this file.
1 
6 #include <time.h>
7 
8 #include <stdlib.h>
9 #include <basic.h>
10 #include <arrayfunctions.h>
11 #include <mpivars.h>
12 #include <hypar.h>
13 
14 static int ReconstructHyperbolic (double*,double*,double*,double*,int,void*,void*,double,int,
15  int(*)(double*,double*,double*,double*,double*,double*,int,void*,double));
16 static int DefaultUpwinding (double*,double*,double*,double*,double*,double*,int,void*,double);
17 
32  double *hyp,
33  double *u,
34  void *s,
35  void *m,
36  double t,
37  int LimFlag,
41  int(*FluxFunction)(double*,double*,int,void*,double),
43  int(*UpwindFunction)(double*,double*,double*,double*,double*,double*,int,void*,double)
44  )
45 {
46  HyPar *solver = (HyPar*) s;
47  MPIVariables *mpi = (MPIVariables*) m;
48  int d, v, i, done;
49  double *FluxI = solver->fluxI; /* interface flux */
50  double *FluxC = solver->fluxC; /* cell centered flux */
52 
53  int ndims = solver->ndims;
54  int nvars = solver->nvars;
55  int ghosts = solver->ghosts;
56  int *dim = solver->dim_local;
57  int size = solver->npoints_local_wghosts;
58  double *x = solver->x;
59  double *dxinv = solver->dxinv;
60  int index[ndims], index1[ndims], index2[ndims], dim_interface[ndims];
61 
62  LimFlag = (LimFlag && solver->flag_nonlinearinterp && solver->SetInterpLimiterVar);
63 
64  _ArraySetValue_(hyp,size*nvars,0.0);
65  _ArraySetValue_(solver->StageBoundaryIntegral,2*ndims*nvars,0.0);
66  if (!FluxFunction) return(0); /* zero hyperbolic term */
67  solver->count_hyp++;
68 
69 #if defined(CPU_STAT)
70  double cpu_time = 0.0;
71  clock_t cpu_start, cpu_end;
72 #endif
73 
74  int offset = 0;
75  for (d = 0; d < ndims; d++) {
76  _ArrayCopy1D_(dim,dim_interface,ndims); dim_interface[d]++;
77  int size_cellcenter = 1; for (i = 0; i < ndims; i++) size_cellcenter *= (dim[i] + 2*ghosts);
78  int size_interface = 1; for (i = 0; i < ndims; i++) size_interface *= dim_interface[i];
79 
80  /* evaluate cell-centered flux */
81  IERR FluxFunction(FluxC,u,d,solver,t); CHECKERR(ierr);
82  /* compute interface fluxes */
83  IERR ReconstructHyperbolic(FluxI,FluxC,u,x+offset,d,solver,mpi,t,LimFlag,UpwindFunction);
84  CHECKERR(ierr);
85 
86  /* calculate the first derivative */
87  done = 0; _ArraySetValue_(index,ndims,0);
88  int p, p1, p2;
89 
90 #if defined(CPU_STAT)
91  cpu_start = clock();
92 #endif
93 
94  while (!done) {
95  _ArrayCopy1D_(index,index1,ndims);
96  _ArrayCopy1D_(index,index2,ndims); index2[d]++;
97  _ArrayIndex1D_(ndims,dim ,index ,ghosts,p);
98  _ArrayIndex1D_(ndims,dim_interface,index1,0 ,p1);
99  _ArrayIndex1D_(ndims,dim_interface,index2,0 ,p2);
100  for (v=0; v<nvars; v++) hyp[nvars*p+v] += dxinv[offset+ghosts+index[d]]
101  * (FluxI[nvars*p2+v]-FluxI[nvars*p1+v]);
102  /* boundary flux integral */
103  if (index[d] == 0)
104  for (v=0; v<nvars; v++) solver->StageBoundaryIntegral[(2*d+0)*nvars+v] -= FluxI[nvars*p1+v];
105  if (index[d] == dim[d]-1)
106  for (v=0; v<nvars; v++) solver->StageBoundaryIntegral[(2*d+1)*nvars+v] += FluxI[nvars*p2+v];
107 
108  _ArrayIncrementIndex_(ndims,dim,index,done);
109  }
110 
111 #if defined(CPU_STAT)
112  cpu_end = clock();
113  cpu_time += (double)(cpu_end - cpu_start) / CLOCKS_PER_SEC;
114 #endif
115 
116  offset += dim[d] + 2*ghosts;
117  }
118 
119 #if defined(CPU_STAT)
120  printf("HyperbolicFunction CPU time = %8.6lf\n", cpu_time);
121 #endif
122 
123  if (solver->flag_ib) _ArrayBlockMultiply_(hyp,solver->iblank,size,nvars);
124 
125  return(0);
126 }
127 
168  double *fluxI,
171  double *fluxC,
173  double *u,
174  double *x,
175  int dir,
176  void *s,
177  void *m,
178  double t,
179  int LimFlag,
183  int(*UpwindFunction)(double*,double*,double*,double*,double*,double*,int,void*,double)
184  )
185 {
186  HyPar *solver = (HyPar*) s;
187  MPIVariables *mpi = (MPIVariables*) m;
189 
190  double *uC = NULL;
191  double *uL = solver->uL;
192  double *uR = solver->uR;
193  double *fluxL = solver->fL;
194  double *fluxR = solver->fR;
195 
196  /*
197  precalculate the non-linear interpolation coefficients if required
198  else reuse the weights previously calculated
199  */
200 
201  if (LimFlag) IERR solver->SetInterpLimiterVar(fluxC,u,x,dir,solver,mpi);
202 
203  /* if defined, calculate the modified u-function to be used for upwinding
204  e.g.: used in well-balanced schemes for Euler/Navier-Stokes with gravity
205  otherwise, just copy u to uC */
206  if (solver->UFunction) {
207  uC = solver->uC;
208  IERR solver->UFunction(uC,u,dir,solver,mpi,t); CHECKERR(ierr);
209  } else uC = u;
210 
211  /* Interpolation -> to calculate left and right-biased interface flux and state variable*/
212  IERR solver->InterpolateInterfacesHyp(uL ,uC ,u,x, 1,dir,solver,mpi,1); CHECKERR(ierr);
213  IERR solver->InterpolateInterfacesHyp(uR ,uC ,u,x,-1,dir,solver,mpi,1); CHECKERR(ierr);
214  IERR solver->InterpolateInterfacesHyp(fluxL,fluxC,u,x, 1,dir,solver,mpi,0); CHECKERR(ierr);
215  IERR solver->InterpolateInterfacesHyp(fluxR,fluxC,u,x,-1,dir,solver,mpi,0); CHECKERR(ierr);
216 
217  /* Upwind -> to calculate the final interface flux */
218  if (UpwindFunction) { IERR UpwindFunction (fluxI,fluxL,fluxR,uL ,uR ,u ,dir,solver,t); CHECKERR(ierr); }
219  else { IERR DefaultUpwinding (fluxI,fluxL,fluxR,NULL,NULL,NULL,dir,solver,t); CHECKERR(ierr); }
220 
221  return(0);
222 }
223 
227  double *fI,
228  double *fL,
229  double *fR,
230  double *uL,
231  double *uR,
232  double *u,
233  int dir,
234  void *s,
235  double t
236  )
237 {
238  HyPar *solver = (HyPar*) s;
239  int done;
240 
241  int *dim = solver->dim_local;
242  int ndims = solver->ndims;
243  int nvars = solver->nvars;
244 
245  int bounds_outer[ndims], bounds_inter[ndims];
246  _ArrayCopy1D_(dim,bounds_outer,ndims); bounds_outer[dir] = 1;
247  _ArrayCopy1D_(dim,bounds_inter,ndims); bounds_inter[dir] += 1;
248 
249  done = 0; int index_outer[ndims], index_inter[ndims];
250  _ArraySetValue_(index_outer,ndims,0);
251  while (!done) {
252  _ArrayCopy1D_(index_outer,index_inter,ndims);
253  for (index_inter[dir] = 0; index_inter[dir] < bounds_inter[dir]; index_inter[dir]++) {
254  int p; _ArrayIndex1D_(ndims,bounds_inter,index_inter,0,p);
255  int v; for (v=0; v<nvars; v++) fI[nvars*p+v] = 0.5 * (fL[nvars*p+v]+fR[nvars*p+v]);
256  }
257  _ArrayIncrementIndex_(ndims,bounds_outer,index_outer,done);
258  }
259 
260  return(0);
261 }
int npoints_local_wghosts
Definition: hypar.h:42
#define _ArraySetValue_(x, size, value)
double * uR
Definition: hypar.h:139
double * iblank
Definition: hypar.h:436
double * uL
Definition: hypar.h:139
static int DefaultUpwinding(double *, double *, double *, double *, double *, double *, int, void *, double)
#define _ArrayIncrementIndex_(N, imax, i, done)
int HyperbolicFunction(double *hyp, double *u, void *s, void *m, double t, int LimFlag, int(*FluxFunction)(double *, double *, int, void *, double), int(*UpwindFunction)(double *, double *, double *, double *, double *, double *, int, void *, double))
double * fL
Definition: hypar.h:139
int flag_ib
Definition: hypar.h:441
int * dim_local
Definition: hypar.h:37
MPI related function definitions.
int ghosts
Definition: hypar.h:52
#define _ArrayIndex1D_(N, imax, i, ghost, index)
int(* SetInterpLimiterVar)(double *, double *, double *, int, void *, void *)
Definition: hypar.h:234
static int ReconstructHyperbolic(double *, double *, double *, double *, int, void *, void *, double, int, int(*)(double *, double *, double *, double *, double *, double *, int, void *, double))
double * fluxI
Definition: hypar.h:136
double * uC
Definition: hypar.h:131
#define _ArrayBlockMultiply_(x, a, n, bs)
#define _ArrayCopy1D_(x, y, size)
double * fluxC
Definition: hypar.h:128
int nvars
Definition: hypar.h:29
int(* InterpolateInterfacesHyp)(double *, double *, double *, double *, int, int, void *, void *, int)
Definition: hypar.h:224
#define CHECKERR(ierr)
Definition: basic.h:18
Contains structure definition for hypar.
double * StageBoundaryIntegral
Definition: hypar.h:382
int(* UFunction)(double *, double *, int, void *, void *, double)
Definition: hypar.h:321
int count_hyp
Definition: hypar.h:418
double * dxinv
Definition: hypar.h:110
Some basic definitions and macros.
int ndims
Definition: hypar.h:26
Contains macros and function definitions for common array operations.
#define IERR
Definition: basic.h:16
Structure of MPI-related variables.
double * x
Definition: hypar.h:107
int flag_nonlinearinterp
Definition: hypar.h:411
Structure containing all solver-specific variables and functions.
Definition: hypar.h:23
double * fR
Definition: hypar.h:139
#define _DECLARE_IERR_
Definition: basic.h:17