HyPar  1.0
Finite-Difference Hyperbolic-Parabolic PDE Solver on Cartesian Grids
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
Interp1PrimFifthOrderUpwind.c
Go to the documentation of this file.
1 
6 #include <stdlib.h>
7 #include <basic.h>
8 #include <arrayfunctions.h>
9 #include <mathfunctions.h>
10 #include <interpolation.h>
11 #include <mpivars.h>
12 #include <hypar.h>
13 
14 #ifdef with_omp
15 #include <omp.h>
16 #endif
17 
18 #undef _MINIMUM_GHOSTS_
19 
23 #define _MINIMUM_GHOSTS_ 3
24 
69  double *fI,
70  double *fC,
71  double *u,
72  double *x,
73  int upw,
74  int dir,
75  void *s,
76  void *m,
77  int uflag
78  )
79 {
80  HyPar *solver = (HyPar*) s;
81 
82  int ghosts = solver->ghosts;
83  int ndims = solver->ndims;
84  int nvars = solver->nvars;
85  int *dim = solver->dim_local;
86  int *stride= solver->stride_with_ghosts;
87 
88  /* define some constants */
89  static const double one_by_thirty = 1.0/30.0,
90  thirteen_by_sixty = 13.0/60.0,
91  fortyseven_by_sixty = 47.0/60.0,
92  twentyseven_by_sixty = 27.0/60.0,
93  one_by_twenty = 1.0/20.0;
94 
95  /* create index and bounds for the outer loop, i.e., to loop over all 1D lines along
96  dimension "dir" */
97  int indexC[ndims], indexI[ndims], index_outer[ndims], bounds_outer[ndims], bounds_inter[ndims];
98  _ArrayCopy1D_(dim,bounds_outer,ndims); bounds_outer[dir] = 1;
99  _ArrayCopy1D_(dim,bounds_inter,ndims); bounds_inter[dir] += 1;
100  int N_outer; _ArrayProduct1D_(bounds_outer,ndims,N_outer);
101 
102  int i;
103 #pragma omp parallel for schedule(auto) default(shared) private(i,index_outer,indexC,indexI)
104  for (i=0; i<N_outer; i++) {
105  _ArrayIndexnD_(ndims,i,bounds_outer,index_outer,0);
106  _ArrayCopy1D_(index_outer,indexC,ndims);
107  _ArrayCopy1D_(index_outer,indexI,ndims);
108  for (indexI[dir] = 0; indexI[dir] < dim[dir]+1; indexI[dir]++) {
109  int qm1,qm2,qm3,qp1,qp2,p;
110  _ArrayIndex1D_(ndims,bounds_inter,indexI,0,p);
111  if (upw > 0) {
112  indexC[dir] = indexI[dir]-1; _ArrayIndex1D_(ndims,dim,indexC,ghosts,qm1);
113  qm3 = qm1 - 2*stride[dir];
114  qm2 = qm1 - stride[dir];
115  qp1 = qm1 + stride[dir];
116  qp2 = qm1 + 2*stride[dir];
117  } else {
118  indexC[dir] = indexI[dir] ; _ArrayIndex1D_(ndims,dim,indexC,ghosts,qm1);
119  qm3 = qm1 + 2*stride[dir];
120  qm2 = qm1 + stride[dir];
121  qp1 = qm1 - stride[dir];
122  qp2 = qm1 - 2*stride[dir];
123  }
124 
125  /* Defining stencil points */
126  double *fm3, *fm2, *fm1, *fp1, *fp2;
127  fm3 = (fC+qm3*nvars);
128  fm2 = (fC+qm2*nvars);
129  fm1 = (fC+qm1*nvars);
130  fp1 = (fC+qp1*nvars);
131  fp2 = (fC+qp2*nvars);
132 
133  int v;
134  for (v=0; v<nvars; v++) {
135  (fI+p*nvars)[v] = one_by_thirty * fm3[v]
136  - thirteen_by_sixty * fm2[v]
137  + fortyseven_by_sixty * fm1[v]
138  + twentyseven_by_sixty * fp1[v]
139  - one_by_twenty * fp2[v];
140  }
141  }
142  }
143 
144  return(0);
145 }
Definitions for the functions computing the interpolated value of the primitive at the cell interface...
#define _ArrayIndexnD_(N, index, imax, i, ghost)
int * dim_local
Definition: hypar.h:37
MPI related function definitions.
int Interp1PrimFifthOrderUpwind(double *, double *, double *, double *, int, int, void *, void *, int)
5th order upwind reconstruction (component-wise) on a uniform grid
int ghosts
Definition: hypar.h:52
#define _ArrayIndex1D_(N, imax, i, ghost, index)
Contains function definitions for common mathematical functions.
int * stride_with_ghosts
Definition: hypar.h:414
#define _ArrayCopy1D_(x, y, size)
int nvars
Definition: hypar.h:29
Contains structure definition for hypar.
Some basic definitions and macros.
int ndims
Definition: hypar.h:26
Contains macros and function definitions for common array operations.
#define _ArrayProduct1D_(x, size, p)
Structure containing all solver-specific variables and functions.
Definition: hypar.h:23