Compute interpolation coefficients for PPM in every coordinate system.
66 double i1, i2, rp, dr, den;
84 for (d = 0; d <
DIMENSIONS; d++) uniform_grid[d] = 1;
88 for (i = 0; i < grid[d].
np_tot-1; i++){
89 if (fabs(grid[d].
dx[i+1]/grid[d].
dx[i]-1.0) > 1.e-9){
132 if (!uniform_grid[d]) {
138 end = grid[d].
np_tot - 1 - iR;
151 #if GEOMETRY == CYLINDRICAL || GEOMETRY == POLAR
152 if (d ==
IDIR)
for (i = beg; i <=
end; i++){
159 den = 12.0*
POLY_3(1.0, -1.0, -3.0, 2.0, i1);
162 wp[
i][-1] =
POLY_3( -3.0, 2.0, 6.0, -4.0, i1)/den;
163 wp[
i][ 0] =
POLY_3( 11.0, -13.0, -28.0, 20.0, i1)/den;
164 wp[
i][ 1] =
POLY_3( 4.0, -1.0, -14.0, 8.0, i1)/den;
168 wm[
i][-1] =
POLY_3( 3.0, -5.0, -10.0, 8.0, i1)/den;
169 wm[
i][ 0] =
POLY_3( 10.0, -9.0, -32.0, 20.0, i1)/den;
170 wm[
i][ 1] =
POLY_3( -1.0, 2.0, 6.0, -4.0, i1)/den;
173 den = 24.0*
POLY_2(4.0, -15.0, 5.0, i2);
176 wp[
i][-1] =
POLY_4(-12.0, -1.0, 30.0, -1.0, -10.0, i1)/den;
177 wp[
i][ 0] =
POLY_4( 60.0, -27.0, -210.0, 13.0, 70.0, i1)/den;
178 wp[
i][ 1] =
POLY_4( 60.0, 27.0, -210.0, -13.0, 70.0, i1)/den;
179 wp[
i][ 2] =
POLY_4(-12.0, 1.0, 30.0, 1.0, -10.0, i1)/den;
182 den = 120.0*(2.0*i1-1.0)*
POLY_4(12.0, 16.0, -13.0, -6.0, 3.0, i1);
184 wp[
i][-2] =
POLY_5( -80, 32, 200, -80, -60, 24, i1)/den;
185 wp[
i][-1] =
POLY_5( 492, -193, -1230, 535, 384, -156, i1)/den;
186 wp[
i][ 0] =
POLY_5(-1276, 1157, 4090, -2075, -1332, 564, i1)/den;
187 wp[
i][ 1] =
POLY_5( -684, 27, 2610, -885, -888, 324, i1)/den;
188 wp[
i][ 2] =
POLY_5( 108, -63, -270, 105, 96, -36, i1)/den;
190 wm[
i][-2] =
POLY_5( 60, -84, -261, 129, 84, -36, i1)/den;
191 wm[
i][-1] =
POLY_5( -504, 660, 2133, -1197, -732, 324, i1)/den;
192 wm[
i][ 0] =
POLY_5(-1128, 604, 4487, -1763, -1488, 564, i1)/den;
193 wm[
i][ 1] =
POLY_5( 168, -292, -1119, 511, 396, -156, i1)/den;
194 wm[
i][ 2] =
POLY_5( -36, 72, 160, -80, -60, 24, i1)/den;
203 #if GEOMETRY == SPHERICAL
204 if (d ==
IDIR)
for (i = beg; i <=
end; i++){
211 den = 18.0*
POLY_6(4.0, -6.0, -9.0, 20.0, 15.0, -30.0, 10.0, i1);
214 wp[
i][-1] = -
POLY_2(7.0, -9.0, 3.0, i1)
215 *
POLY_2(3.0, -9.0, 10.0, i2)/den;
216 wp[
i][ 0] =
POLY_2(1.0, -3.0, 3.0, i1)
217 *
POLY_4(69.0, 96.0, -63.0, -90.0, 50.0, i1)/den;
218 wp[
i][ 1] = 2.0*
POLY_2( 1.0, 3.0, 3.0, i1)
219 *
POLY_4(12.0, -48.0, 72.0, -45.0, 10.0, i1)/den;
223 wm[
i][-1] = 2.0*
POLY_2(7.0, -9.0, 3.0, i1)
224 *
POLY_4(1.0, -1.0, -3.0, 5.0, 10.0, i1)/den;
225 wm[
i][ 0] =
POLY_2(1.0, -3.0, 3.0, i1)
226 *
POLY_4(62.0, 100.0, -33.0, -110.0, 50.0, i1)/den;
227 wm[
i][ 1] = -
POLY_2(1.0, 3.0, 3.0, i1)
228 *
POLY_4(4.0, -22.0, 51.0, -40.0, 10.0, i1)/den;
231 den = 36.0*
POLY_4(16.0, -60.0, 150.0, -85.0, 15.0, i2);
234 wp[
i][-1] = -
POLY_2( 7, -9, 3, i1)/den
235 *
POLY_6(12, 16, -30, -48.0, 23, 48, 15, i1);
236 wp[
i][ 0] =
POLY_2( 1, -3, 3, i1)/den
237 *
POLY_6(372, 1008.0, 510, -720, -487, 144, 105, i1);
238 wp[
i][ 1] =
POLY_2( 1, 3, 3.0, i1)/den
239 *
POLY_6(372, -1008, 510, 720, -487, -144, 105, i1);
240 wp[
i][ 2] = -
POLY_2( 7, 9, 3, i1)/den
241 *
POLY_6(12, -16, -30, 48, 23, -48, 15, i1);
243 den =
POLY_10( 48.0, -48.0, -164.0, 200.0, 390.0,
244 -399.0, -161.0, 210.0, 0.0, -35.0, 7.0, i1);
247 wp[
i][-2] = 2.0*
POLY_2(19., -15., 3., i1)/den
248 *
POLY_4(16., -60., 94., -45., 7., i2);
250 wp[
i][-1] = -
POLY_2(7.0, -9.0, 3.0, i1)/den;
251 wp[
i][-1] *=
POLY_8(508., 240., -1740., -795., 2417.,
252 930., -780., -175., 91., i1);
254 wp[
i][0] =
POLY_2(1.0, -3.0, 3.0, i1)/den;
255 wp[
i][0] *=
POLY_8(8132., 15120., -5700., -20325., 3863.,
256 8670., -1800., -1225., 329., i1);
258 wp[
i][1] =
POLY_2(1.0, 3.0, 3.0, i1)/den;
259 wp[
i][1] *=
POLY_8(4212., -15120., 16560., 1275., -11517.,
260 4350., 1620., -1225., 189., i1);
262 wp[
i][2] = -
POLY_2(7.0, 9.0, 3.0, i1)/den;
263 wp[
i][2] *=
POLY_8( 108., -240., -120., 645., -223.,
264 -510., 510., -175., 21., i1);
266 wm[
i][-2] = -
POLY_2( 19., -15., 3., i1)/den;
267 wm[
i][-2] *=
POLY_8( 16., -16., -60., 96., 222.,
268 -51., -127., 7., 21., i1);
270 wm[
i][-1] =
POLY_2( 7., -9., 3., i1)/den;
271 wm[
i][-1] *=
POLY_8( 344., -164., -1350., 1184., 4888.,
272 1071., -1663., -287., 189., i1);
274 wm[
i][0] =
POLY_2(1.0, -3.0, 3.0, i1)/den;
275 wm[
i][0] *=
POLY_8(7064., 15196., -310., -21376., 368.,
276 9431., -1163., -1407., 329., i1);
278 wm[
i][1] = -
POLY_2(1.0, 3.0, 3.0, i1)/den;
279 wm[
i][1] *=
POLY_8( 696., -3516., 6850., -1544., -4388.,
280 2329., 543., -553., 91., i1);
282 wm[
i][2] = 2.0*
POLY_2(7.0, 9.0, 3.0, i1)/den;
283 wm[
i][2] *=
POLY_8( 12., -42., 25., 132., -91.,
284 -122., 151., -56., 7., i1);
286 }
else if (d ==
JDIR) {
290 }
else if (d ==
KDIR){
int end
Global end index for the local array.
#define POLY_6(a0, a1, a2, a3, a4, a5, a6, x)
double *** ArrayBox(long int nrl, long int nrh, long int ncl, long int nch, long int ndl, long int ndh)
#define POLY_5(a0, a1, a2, a3, a4, a5, x)
static void PPM_FindWeights(double **, double **, int, Grid *)
int beg
Global start index for the local array.
void PPM_Q6_Coeffs(double *hp, double *hm, int dir, Grid *grid)
static void PPM_CartCoeff(double **, double **, int, int)
long int NMAX_POINT
Maximum number of points among the three directions, boundaries excluded.
#define POLY_3(a0, a1, a2, a3, x)
#define POLY_2(a0, a1, a2, x)
#define POLY_4(a0, a1, a2, a3, a4, x)
#define ARRAY_2D(nx, ny, type)
int np_tot
Total number of points in the local domain (boundaries included).
#define POLY_10(a0, a1, a2, a3, a4, a5, a6, a7, a8, a9, a10, x)
#define POLY_8(a0, a1, a2, a3, a4, a5, a6, a7, a8, x)