HyPar  1.0
Finite-Difference Hyperbolic-Parabolic PDE Solver on Cartesian Grids
ShallowWater2DUpwind.c
Go to the documentation of this file.
1 
6 #include <stdlib.h>
7 #include <math.h>
8 #include <basic.h>
9 #include <arrayfunctions.h>
10 #include <mathfunctions.h>
11 #include <matmult_native.h>
12 #include <physicalmodels/shallowwater2d.h>
13 #include <hypar.h>
14 
31 int ShallowWater2DUpwindRoe(
32  double *fI,
33  double *fL,
34  double *fR,
35  double *uL,
36  double *uR,
37  double *u,
38  int dir,
39  void *s,
40  double t
41  )
42 {
43  HyPar *solver = (HyPar*) s;
44  ShallowWater2D *param = (ShallowWater2D*) solver->physics;
45  int done,k;
46 
47  int ndims = solver->ndims;
48  int ghosts= solver->ghosts;
49  int *dim = solver->dim_local;
50 
51  int index_outer[ndims], index_inter[ndims], bounds_outer[ndims], bounds_inter[ndims];
52  _ArrayCopy1D_(dim,bounds_outer,ndims); bounds_outer[dir] = 1;
53  _ArrayCopy1D_(dim,bounds_inter,ndims); bounds_inter[dir] += 1;
54  static double R[_MODEL_NVARS_*_MODEL_NVARS_], D[_MODEL_NVARS_*_MODEL_NVARS_], L[_MODEL_NVARS_*_MODEL_NVARS_],
55  DL[_MODEL_NVARS_*_MODEL_NVARS_], modA[_MODEL_NVARS_*_MODEL_NVARS_];
56 
57  done = 0; _ArraySetValue_(index_outer,ndims,0);
58  while (!done) {
59  _ArrayCopy1D_(index_outer,index_inter,ndims);
60  for (index_inter[dir] = 0; index_inter[dir] < bounds_inter[dir]; index_inter[dir]++) {
61  int p; _ArrayIndex1D_(ndims,bounds_inter,index_inter,0,p);
62  int indexL[ndims]; _ArrayCopy1D_(index_inter,indexL,ndims); indexL[dir]--;
63  int indexR[ndims]; _ArrayCopy1D_(index_inter,indexR,ndims);
64  int pL; _ArrayIndex1D_(ndims,dim,indexL,ghosts,pL);
65  int pR; _ArrayIndex1D_(ndims,dim,indexR,ghosts,pR);
66  double udiff[_MODEL_NVARS_], uavg[_MODEL_NVARS_],udiss[_MODEL_NVARS_];
67 
68  /* Roe's upwinding scheme */
69 
70  udiff[0] = 0.5 * (uR[_MODEL_NVARS_*p+0] - uL[_MODEL_NVARS_*p+0]);
71  udiff[1] = 0.5 * (uR[_MODEL_NVARS_*p+1] - uL[_MODEL_NVARS_*p+1]);
72  udiff[2] = 0.5 * (uR[_MODEL_NVARS_*p+2] - uL[_MODEL_NVARS_*p+2]);
73 
74  _ShallowWater2DRoeAverage_ (uavg,(u+_MODEL_NVARS_*pL),(u+_MODEL_NVARS_*pR),param);
75  _ShallowWater2DEigenvalues_ (uavg,D,param,dir);
76  _ShallowWater2DLeftEigenvectors_ (uavg,L,param,dir);
77  _ShallowWater2DRightEigenvectors_ (uavg,R,param,dir);
78 
79  k = 0; D[k] = absolute(D[k]);
80  k = 4; D[k] = absolute(D[k]);
81  k = 8; D[k] = absolute(D[k]);
82 
83  MatMult3(_MODEL_NVARS_,DL,D,L);
84  MatMult3(_MODEL_NVARS_,modA,R,DL);
85 
86  udiss[0] = modA[0*_MODEL_NVARS_+0]*udiff[0] + modA[0*_MODEL_NVARS_+1]*udiff[1] + modA[0*_MODEL_NVARS_+2]*udiff[2];
87  udiss[1] = modA[1*_MODEL_NVARS_+0]*udiff[0] + modA[1*_MODEL_NVARS_+1]*udiff[1] + modA[1*_MODEL_NVARS_+2]*udiff[2];
88  udiss[2] = modA[2*_MODEL_NVARS_+0]*udiff[0] + modA[2*_MODEL_NVARS_+1]*udiff[1] + modA[2*_MODEL_NVARS_+2]*udiff[2];
89 
90  fI[_MODEL_NVARS_*p+0] = 0.5 * (fL[_MODEL_NVARS_*p+0]+fR[_MODEL_NVARS_*p+0]) - udiss[0];
91  fI[_MODEL_NVARS_*p+1] = 0.5 * (fL[_MODEL_NVARS_*p+1]+fR[_MODEL_NVARS_*p+1]) - udiss[1];
92  fI[_MODEL_NVARS_*p+2] = 0.5 * (fL[_MODEL_NVARS_*p+2]+fR[_MODEL_NVARS_*p+2]) - udiss[2];
93  }
94  _ArrayIncrementIndex_(ndims,bounds_outer,index_outer,done);
95  }
96 
97  return(0);
98 }
99 
118 int ShallowWater2DUpwindLLF(
119  double *fI,
120  double *fL,
121  double *fR,
122  double *uL,
123  double *uR,
124  double *u,
125  int dir,
126  void *s,
127  double t
128  )
129 {
130  HyPar *solver = (HyPar*) s;
131  ShallowWater2D *param = (ShallowWater2D*) solver->physics;
132  int done,k;
133 
134  int ndims = solver->ndims;
135  int *dim = solver->dim_local;
136  int ghosts = solver->ghosts;
137 
138  int index_outer[ndims], index_inter[ndims], bounds_outer[ndims], bounds_inter[ndims];
139  _ArrayCopy1D_(dim,bounds_outer,ndims); bounds_outer[dir] = 1;
140  _ArrayCopy1D_(dim,bounds_inter,ndims); bounds_inter[dir] += 1;
141  static double R[_MODEL_NVARS_*_MODEL_NVARS_], D[_MODEL_NVARS_*_MODEL_NVARS_], L[_MODEL_NVARS_*_MODEL_NVARS_];
142 
143  done = 0; _ArraySetValue_(index_outer,ndims,0);
144  while (!done) {
145  _ArrayCopy1D_(index_outer,index_inter,ndims);
146  for (index_inter[dir] = 0; index_inter[dir] < bounds_inter[dir]; index_inter[dir]++) {
147  int p; _ArrayIndex1D_(ndims,bounds_inter,index_inter,0,p);
148  int indexL[ndims]; _ArrayCopy1D_(index_inter,indexL,ndims); indexL[dir]--;
149  int indexR[ndims]; _ArrayCopy1D_(index_inter,indexR,ndims);
150  int pL; _ArrayIndex1D_(ndims,dim,indexL,ghosts,pL);
151  int pR; _ArrayIndex1D_(ndims,dim,indexR,ghosts,pR);
152  double uavg[_MODEL_NVARS_], fcL[_MODEL_NVARS_], fcR[_MODEL_NVARS_],
153  ucL[_MODEL_NVARS_], ucR[_MODEL_NVARS_], fc[_MODEL_NVARS_];
154 
155  /* Local Lax-Friedrich upwinding scheme */
156 
157  _ShallowWater2DRoeAverage_ (uavg,(u+_MODEL_NVARS_*pL),(u+_MODEL_NVARS_*pR),param);
158  _ShallowWater2DEigenvalues_ (uavg,D,param,dir);
159  _ShallowWater2DLeftEigenvectors_ (uavg,L,param,dir);
160  _ShallowWater2DRightEigenvectors_(uavg,R,param,dir);
161 
162  /* calculate characteristic fluxes and variables */
163  MatVecMult3(_MODEL_NVARS_,ucL,L,(uL+_MODEL_NVARS_*p));
164  MatVecMult3(_MODEL_NVARS_,ucR,L,(uR+_MODEL_NVARS_*p));
165  MatVecMult3(_MODEL_NVARS_,fcL,L,(fL+_MODEL_NVARS_*p));
166  MatVecMult3(_MODEL_NVARS_,fcR,L,(fR+_MODEL_NVARS_*p));
167 
168  for (k = 0; k < _MODEL_NVARS_; k++) {
169  double eigL,eigC,eigR;
170  _ShallowWater2DEigenvalues_((u+_MODEL_NVARS_*pL),D,param,dir);
171  eigL = D[k*_MODEL_NVARS_+k];
172  _ShallowWater2DEigenvalues_((u+_MODEL_NVARS_*pR),D,param,dir);
173  eigR = D[k*_MODEL_NVARS_+k];
174  _ShallowWater2DEigenvalues_(uavg,D,param,dir);
175  eigC = D[k*_MODEL_NVARS_+k];
176 
177  double alpha = max3(absolute(eigL),absolute(eigC),absolute(eigR));
178  fc[k] = 0.5 * (fcL[k] + fcR[k] + alpha * (ucL[k]-ucR[k]));
179  }
180 
181  /* calculate the interface flux from the characteristic flux */
182  MatVecMult3(_MODEL_NVARS_,(fI+_MODEL_NVARS_*p),R,fc);
183  }
184  _ArrayIncrementIndex_(ndims,bounds_outer,index_outer,done);
185  }
186 
187  return(0);
188 }
189 
max3
#define max3(a, b, c)
Definition: math_ops.h:27
absolute
#define absolute(a)
Definition: math_ops.h:32
mathfunctions.h
Contains function definitions for common mathematical functions.
basic.h
Some basic definitions and macros.
_ShallowWater2DRoeAverage_
#define _ShallowWater2DRoeAverage_(uavg, uL, uR, p)
Definition: shallowwater2d.h:108
shallowwater2d.h
2D Shallow Water Equations
HyPar::ndims
int ndims
Definition: hypar.h:26
ShallowWater2D
Structure containing variables and parameters specific to the 2D Shallow Water equations. This structure contains the physical parameters, variables, and function pointers specific to the 2D ShallowWater equations.
Definition: shallowwater2d.h:240
HyPar
Structure containing all solver-specific variables and functions.
Definition: hypar.h:23
_ShallowWater2DLeftEigenvectors_
#define _ShallowWater2DLeftEigenvectors_(u, L, p, dir)
Definition: shallowwater2d.h:154
hypar.h
Contains structure definition for hypar.
MatVecMult3
#define MatVecMult3(N, y, A, x)
Definition: matmult_native.h:137
_ArrayIndex1D_
#define _ArrayIndex1D_(N, imax, i, ghost, index)
Definition: arrayfunctions.h:40
_ShallowWater2DEigenvalues_
#define _ShallowWater2DEigenvalues_(u, D, p, dir)
Definition: shallowwater2d.h:129
ShallowWater2DUpwindRoe
int ShallowWater2DUpwindRoe(double *fI, double *fL, double *fR, double *uL, double *uR, double *u, int dir, void *s, double t)
Definition: ShallowWater2DUpwind.c:31
_ArraySetValue_
#define _ArraySetValue_(x, size, value)
Definition: arrayfunctions.h:169
HyPar::dim_local
int * dim_local
Definition: hypar.h:37
_ArrayIncrementIndex_
#define _ArrayIncrementIndex_(N, imax, i, done)
Definition: arrayfunctions.h:126
HyPar::physics
void * physics
Definition: hypar.h:266
HyPar::ghosts
int ghosts
Definition: hypar.h:52
matmult_native.h
Contains macros and function definitions for common matrix multiplication.
MatMult3
#define MatMult3(N, A, X, Y)
Definition: matmult_native.h:43
ShallowWater2DUpwindLLF
int ShallowWater2DUpwindLLF(double *fI, double *fL, double *fR, double *uL, double *uR, double *u, int dir, void *s, double t)
Definition: ShallowWater2DUpwind.c:118
_MODEL_NVARS_
#define _MODEL_NVARS_
Definition: euler1d.h:58
_ShallowWater2DRightEigenvectors_
#define _ShallowWater2DRightEigenvectors_(u, R, p, dir)
Definition: shallowwater2d.h:197
_ArrayCopy1D_
#define _ArrayCopy1D_(x, y, size)
Definition: arrayfunctions.h:319
arrayfunctions.h
Contains macros and function definitions for common array operations.