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

Slip-wall boundary conditions. More...

#include <stdlib.h>
#include <basic.h>
#include <arrayfunctions.h>
#include <boundaryconditions.h>
#include <physicalmodels/euler1d.h>
#include <physicalmodels/euler2d.h>
#include <physicalmodels/navierstokes3d.h>

Go to the source code of this file.

Functions

int BCSlipWallU (void *b, void *m, int ndims, int nvars, int *size, int ghosts, double *phi, double waqt)
 

Detailed Description

Slip-wall boundary conditions.

Author
Debojyoti Ghosh

Definition in file BCSlipWall.c.

Function Documentation

int BCSlipWallU ( void *  b,
void *  m,
int  ndims,
int  nvars,
int *  size,
int  ghosts,
double *  phi,
double  waqt 
)

Applies the slip-wall boundary condition: This is specific to the two and three dimensional Euler and Navier-Stokes systems (Euler2D, NavierStokes2D, NavierStokes3D). It is used for simulating inviscid walls or symmetric boundaries. The pressure, density, and tangential velocity at the ghost points are extrapolated from the interior, while the normal velocity at the ghost points is set such that the interpolated value at the boundary face is equal to the specified wall velocity.

Parameters
bBoundary object of type DomainBoundary
mMPI object of type MPIVariables
ndimsNumber of spatial dimensions
nvarsNumber of variables/DoFs per grid point
sizeInteger array with the number of grid points in each spatial dimension
ghostsNumber of ghost points
phiThe solution array on which to apply the boundary condition
waqtCurrent solution time

Definition at line 22 of file BCSlipWall.c.

32 {
33  DomainBoundary *boundary = (DomainBoundary*) b;
34 
35  int dim = boundary->dim;
36  int face = boundary->face;
37 
38  if (ndims == 1) {
39 
40  /* create a fake physics object */
41  Euler1D physics;
42  double gamma;
43  gamma = physics.gamma = boundary->gamma;
44  double inv_gamma_m1 = 1.0/(gamma-1.0);
45 
46  if (boundary->on_this_proc) {
47  int bounds[ndims], indexb[ndims], indexi[ndims];
48  _ArraySubtract1D_(bounds,boundary->ie,boundary->is,ndims);
49  _ArraySetValue_(indexb,ndims,0);
50  int done = 0;
51  while (!done) {
52  int p1, p2;
53  _ArrayCopy1D_(indexb,indexi,ndims);
54  _ArrayAdd1D_(indexi,indexi,boundary->is,ndims);
55  if (face == 1) indexi[dim] = ghosts-1-indexb[dim];
56  else if (face == -1) indexi[dim] = size[dim]-indexb[dim]-1;
57  else return(1);
58  _ArrayIndex1DWO_(ndims,size,indexb,boundary->is,ghosts,p1);
59  _ArrayIndex1D_(ndims,size,indexi,ghosts,p2);
60 
61  /* flow variables in the interior */
62  double rho, uvel, energy, pressure;
63  double rho_gpt, uvel_gpt, energy_gpt, pressure_gpt;
64  _Euler1DGetFlowVar_((phi+nvars*p2),rho,uvel,energy,pressure,(&physics));
65  /* set the ghost point values */
66  rho_gpt = rho;
67  pressure_gpt = pressure;
68  uvel_gpt = 2.0*boundary->FlowVelocity[_XDIR_] - uvel;
69  energy_gpt = inv_gamma_m1*pressure_gpt
70  + 0.5 * rho_gpt * uvel_gpt*uvel_gpt;
71 
72  phi[nvars*p1+0] = rho_gpt;
73  phi[nvars*p1+1] = rho_gpt * uvel_gpt;
74  phi[nvars*p1+2] = energy_gpt;
75 
76  _ArrayIncrementIndex_(ndims,bounds,indexb,done);
77  }
78  }
79 
80  } else if (ndims == 2) {
81 
82  /* create a fake physics object */
83  Euler2D physics;
84  double gamma;
85  gamma = physics.gamma = boundary->gamma;
86  double inv_gamma_m1 = 1.0/(gamma-1.0);
87 
88  if (boundary->on_this_proc) {
89  int bounds[ndims], indexb[ndims], indexi[ndims];
90  _ArraySubtract1D_(bounds,boundary->ie,boundary->is,ndims);
91  _ArraySetValue_(indexb,ndims,0);
92  int done = 0;
93  while (!done) {
94  int p1, p2;
95  _ArrayCopy1D_(indexb,indexi,ndims);
96  _ArrayAdd1D_(indexi,indexi,boundary->is,ndims);
97  if (face == 1) indexi[dim] = ghosts-1-indexb[dim];
98  else if (face == -1) indexi[dim] = size[dim]-indexb[dim]-1;
99  else return(1);
100  _ArrayIndex1DWO_(ndims,size,indexb,boundary->is,ghosts,p1);
101  _ArrayIndex1D_(ndims,size,indexi,ghosts,p2);
102 
103  /* flow variables in the interior */
104  double rho, uvel, vvel, energy, pressure;
105  double rho_gpt, uvel_gpt, vvel_gpt, energy_gpt, pressure_gpt;
106  _Euler2DGetFlowVar_((phi+nvars*p2),rho,uvel,vvel,energy,pressure,(&physics));
107  /* set the ghost point values */
108  rho_gpt = rho;
109  pressure_gpt = pressure;
110  if (dim == _XDIR_) {
111  uvel_gpt = 2.0*boundary->FlowVelocity[_XDIR_] - uvel;
112  vvel_gpt = vvel;
113  } else if (dim == _YDIR_) {
114  uvel_gpt = uvel;
115  vvel_gpt = 2.0*boundary->FlowVelocity[_YDIR_] - vvel;
116  } else {
117  uvel_gpt = 0.0;
118  vvel_gpt = 0.0;
119  }
120  energy_gpt = inv_gamma_m1*pressure_gpt
121  + 0.5 * rho_gpt * (uvel_gpt*uvel_gpt + vvel_gpt*vvel_gpt);
122 
123  phi[nvars*p1+0] = rho_gpt;
124  phi[nvars*p1+1] = rho_gpt * uvel_gpt;
125  phi[nvars*p1+2] = rho_gpt * vvel_gpt;
126  phi[nvars*p1+3] = energy_gpt;
127 
128  _ArrayIncrementIndex_(ndims,bounds,indexb,done);
129  }
130  }
131 
132  } else if (ndims == 3) {
133 
134  /* create a fake physics object */
135  double gamma;
136  gamma = boundary->gamma;
137  double inv_gamma_m1 = 1.0/(gamma-1.0);
138 
139  if (boundary->on_this_proc) {
140  int bounds[ndims], indexb[ndims], indexi[ndims];
141  _ArraySubtract1D_(bounds,boundary->ie,boundary->is,ndims);
142  _ArraySetValue_(indexb,ndims,0);
143  int done = 0;
144  while (!done) {
145  int p1, p2;
146  _ArrayCopy1D_(indexb,indexi,ndims);
147  _ArrayAdd1D_(indexi,indexi,boundary->is,ndims);
148  if (face == 1) indexi[dim] = ghosts-1-indexb[dim];
149  else if (face == -1) indexi[dim] = size[dim]-indexb[dim]-1;
150  else return(1);
151  _ArrayIndex1DWO_(ndims,size,indexb,boundary->is,ghosts,p1);
152  _ArrayIndex1D_(ndims,size,indexi,ghosts,p2);
153 
154  /* flow variables in the interior */
155  double rho, uvel, vvel, wvel, energy, pressure;
156  double rho_gpt, uvel_gpt, vvel_gpt, wvel_gpt, energy_gpt, pressure_gpt;
157  _NavierStokes3DGetFlowVar_((phi+nvars*p2),_NavierStokes3D_stride_,rho,uvel,vvel,wvel,energy,pressure,gamma);
158  /* set the ghost point values */
159  rho_gpt = rho;
160  pressure_gpt = pressure;
161  if (dim == _XDIR_) {
162  uvel_gpt = 2.0*boundary->FlowVelocity[_XDIR_] - uvel;
163  vvel_gpt = vvel;
164  wvel_gpt = wvel;
165  } else if (dim == _YDIR_) {
166  uvel_gpt = uvel;
167  vvel_gpt = 2.0*boundary->FlowVelocity[_YDIR_] - vvel;
168  wvel_gpt = wvel;
169  } else if (dim == _ZDIR_) {
170  uvel_gpt = uvel;
171  vvel_gpt = vvel;
172  wvel_gpt = 2.0*boundary->FlowVelocity[_ZDIR_] - wvel;
173  } else {
174  uvel_gpt = 0.0;
175  vvel_gpt = 0.0;
176  wvel_gpt = 0.0;
177  }
178  energy_gpt = inv_gamma_m1*pressure_gpt
179  + 0.5 * rho_gpt
180  * (uvel_gpt*uvel_gpt + vvel_gpt*vvel_gpt + wvel_gpt*wvel_gpt);
181 
182  phi[nvars*p1+0] = rho_gpt;
183  phi[nvars*p1+1] = rho_gpt * uvel_gpt;
184  phi[nvars*p1+2] = rho_gpt * vvel_gpt;
185  phi[nvars*p1+3] = rho_gpt * wvel_gpt;
186  phi[nvars*p1+4] = energy_gpt;
187 
188  _ArrayIncrementIndex_(ndims,bounds,indexb,done);
189  }
190  }
191 
192  }
193  return(0);
194 }
#define _YDIR_
Definition: euler2d.h:41
#define _NavierStokes3DGetFlowVar_(u, stride, rho, vx, vy, vz, e, P, gamma)
#define _ArraySetValue_(x, size, value)
#define _Euler1DGetFlowVar_(u, rho, v, e, P, p)
Definition: euler1d.h:81
#define _ArrayIncrementIndex_(N, imax, i, done)
#define _Euler2DGetFlowVar_(u, rho, vx, vy, e, P, p)
Definition: euler2d.h:44
#define _ArraySubtract1D_(x, a, b, size)
#define _ArrayIndex1D_(N, imax, i, ghost, index)
#define _ZDIR_
Structure containing the variables and function pointers defining a boundary.
#define _ArrayIndex1DWO_(N, imax, i, offset, ghost, index)
#define _ArrayCopy1D_(x, y, size)
Structure containing variables and parameters specific to the 1D Euler equations. This structure cont...
Definition: euler1d.h:273
double gamma
Definition: euler1d.h:274
#define _ArrayAdd1D_(x, a, b, size)
#define _XDIR_
Definition: euler1d.h:75
double gamma
Definition: euler2d.h:245
static const int _NavierStokes3D_stride_