Compute interpolation coefficients for PPM in every coordinate system.
31 double i1, i2, rp, dr, den;
48 for (d = 0; d <
DIMENSIONS; d++) uniform_grid[d] = 1;
52 for (i = 0; i < grid[d].
np_tot-1; i++){
53 if (fabs(grid[d].
dx[i+1]/grid[d].
dx[i]-1.0) > 1.e-9){
66 #if PPM_VERSION == PPM3
68 #elif PPM_VERSION == PPM4
70 #elif PPM_VERSION == PPM5
85 end = grid[d].
np_tot - 1 - iR;
98 #if GEOMETRY == CYLINDRICAL || GEOMETRY == POLAR
99 if (d ==
IDIR)
for (i = beg; i <=
end; i++){
105 #if PPM_VERSION == PPM3
106 den = 12.0*
POLY_3(1.0, -1.0, -3.0, 2.0, i1);
109 wp[
i][-1] =
POLY_3( -3.0, 2.0, 6.0, -4.0, i1)/den;
110 wp[
i][ 0] =
POLY_3( 11.0, -13.0, -28.0, 20.0, i1)/den;
111 wp[
i][ 1] =
POLY_3( 4.0, -1.0, -14.0, 8.0, i1)/den;
115 wm[
i][-1] =
POLY_3( 3.0, -5.0, -10.0, 8.0, i1)/den;
116 wm[
i][ 0] =
POLY_3( 10.0, -9.0, -32.0, 20.0, i1)/den;
117 wm[
i][ 1] =
POLY_3( -1.0, 2.0, 6.0, -4.0, i1)/den;
119 #elif PPM_VERSION == PPM4
120 den = 24.0*
POLY_2(4.0, -15.0, 5.0, i2);
123 wp[
i][-1] =
POLY_4(-12.0, -1.0, 30.0, -1.0, -10.0, i1)/den;
124 wp[
i][ 0] =
POLY_4( 60.0, -27.0, -210.0, 13.0, 70.0, i1)/den;
125 wp[
i][ 1] =
POLY_4( 60.0, 27.0, -210.0, -13.0, 70.0, i1)/den;
126 wp[
i][ 2] =
POLY_4(-12.0, 1.0, 30.0, 1.0, -10.0, i1)/den;
128 #elif PPM_VERSION == PPM5
129 den = 120.0*(2.0*i1-1.0)*
POLY_4(12.0, 16.0, -13.0, -6.0, 3.0, i1);
131 wp[
i][-2] =
POLY_5( -80, 32, 200, -80, -60, 24, i1)/den;
132 wp[
i][-1] =
POLY_5( 492, -193, -1230, 535, 384, -156, i1)/den;
133 wp[
i][ 0] =
POLY_5(-1276, 1157, 4090, -2075, -1332, 564, i1)/den;
134 wp[
i][ 1] =
POLY_5( -684, 27, 2610, -885, -888, 324, i1)/den;
135 wp[
i][ 2] =
POLY_5( 108, -63, -270, 105, 96, -36, i1)/den;
137 wm[
i][-2] =
POLY_5( 60, -84, -261, 129, 84, -36, i1)/den;
138 wm[
i][-1] =
POLY_5( -504, 660, 2133, -1197, -732, 324, i1)/den;
139 wm[
i][ 0] =
POLY_5(-1128, 604, 4487, -1763, -1488, 564, i1)/den;
140 wm[
i][ 1] =
POLY_5( 168, -292, -1119, 511, 396, -156, i1)/den;
141 wm[
i][ 2] =
POLY_5( -36, 72, 160, -80, -60, 24, i1)/den;
150 #if GEOMETRY == SPHERICAL
151 if (d ==
IDIR)
for (i = beg; i <=
end; i++){
157 #if PPM_VERSION == PPM3
158 den = 18.0*
POLY_6(4.0, -6.0, -9.0, 20.0, 15.0, -30.0, 10.0, i1);
161 wp[
i][-1] = -
POLY_2(7.0, -9.0, 3.0, i1)
162 *
POLY_2(3.0, -9.0, 10.0, i2)/den;
163 wp[
i][ 0] =
POLY_2(1.0, -3.0, 3.0, i1)
164 *
POLY_4(69.0, 96.0, -63.0, -90.0, 50.0, i1)/den;
165 wp[
i][ 1] = 2.0*
POLY_2( 1.0, 3.0, 3.0, i1)
166 *
POLY_4(12.0, -48.0, 72.0, -45.0, 10.0, i1)/den;
170 wm[
i][-1] = 2.0*
POLY_2(7.0, -9.0, 3.0, i1)
171 *
POLY_4(1.0, -1.0, -3.0, 5.0, 10.0, i1)/den;
172 wm[
i][ 0] =
POLY_2(1.0, -3.0, 3.0, i1)
173 *
POLY_4(62.0, 100.0, -33.0, -110.0, 50.0, i1)/den;
174 wm[
i][ 1] = -
POLY_2(1.0, 3.0, 3.0, i1)
175 *
POLY_4(4.0, -22.0, 51.0, -40.0, 10.0, i1)/den;
177 #elif PPM_VERSION == PPM4
178 den = 36.0*
POLY_4(16.0, -60.0, 150.0, -85.0, 15.0, i2);
181 wp[
i][-1] = -
POLY_2( 7, -9, 3, i1)/den
182 *
POLY_6(12, 16, -30, -48.0, 23, 48, 15, i1);
183 wp[
i][ 0] =
POLY_2( 1, -3, 3, i1)/den
184 *
POLY_6(372, 1008.0, 510, -720, -487, 144, 105, i1);
185 wp[
i][ 1] =
POLY_2( 1, 3, 3.0, i1)/den
186 *
POLY_6(372, -1008, 510, 720, -487, -144, 105, i1);
187 wp[
i][ 2] = -
POLY_2( 7, 9, 3, i1)/den
188 *
POLY_6(12, -16, -30, 48, 23, -48, 15, i1);
189 #elif PPM_VERSION == PPM5
190 den =
POLY_10( 48.0, -48.0, -164.0, 200.0, 390.0,
191 -399.0, -161.0, 210.0, 0.0, -35.0, 7.0, i1);
194 wp[
i][-2] = 2.0*
POLY_2(19., -15., 3., i1)/den
195 *
POLY_4(16., -60., 94., -45., 7., i2);
197 wp[
i][-1] = -
POLY_2(7.0, -9.0, 3.0, i1)/den;
198 wp[
i][-1] *=
POLY_8(508., 240., -1740., -795., 2417.,
199 930., -780., -175., 91., i1);
201 wp[
i][0] =
POLY_2(1.0, -3.0, 3.0, i1)/den;
202 wp[
i][0] *=
POLY_8(8132., 15120., -5700., -20325., 3863.,
203 8670., -1800., -1225., 329., i1);
205 wp[
i][1] =
POLY_2(1.0, 3.0, 3.0, i1)/den;
206 wp[
i][1] *=
POLY_8(4212., -15120., 16560., 1275., -11517.,
207 4350., 1620., -1225., 189., i1);
209 wp[
i][2] = -
POLY_2(7.0, 9.0, 3.0, i1)/den;
210 wp[
i][2] *=
POLY_8( 108., -240., -120., 645., -223.,
211 -510., 510., -175., 21., i1);
213 wm[
i][-2] = -
POLY_2( 19., -15., 3., i1)/den;
214 wm[
i][-2] *=
POLY_8( 16., -16., -60., 96., 222.,
215 -51., -127., 7., 21., i1);
217 wm[
i][-1] =
POLY_2( 7., -9., 3., i1)/den;
218 wm[
i][-1] *=
POLY_8( 344., -164., -1350., 1184., 4888.,
219 1071., -1663., -287., 189., i1);
221 wm[
i][0] =
POLY_2(1.0, -3.0, 3.0, i1)/den;
222 wm[
i][0] *=
POLY_8(7064., 15196., -310., -21376., 368.,
223 9431., -1163., -1407., 329., i1);
225 wm[
i][1] = -
POLY_2(1.0, 3.0, 3.0, i1)/den;
226 wm[
i][1] *=
POLY_8( 696., -3516., 6850., -1544., -4388.,
227 2329., 543., -553., 91., i1);
229 wm[
i][2] = 2.0*
POLY_2(7.0, 9.0, 3.0, i1)/den;
230 wm[
i][2] *=
POLY_8( 12., -42., 25., 132., -91.,
231 -122., 151., -56., 7., i1);
233 }
else if (d ==
JDIR) {
237 }
else if (d ==
KDIR){
void PPM_Q6_Coeffs(double *hp, double *hm, int dir, Grid *grid)
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)
int beg
Global start index for the local array.
static void PPM_FindWeights(double **, double **, int, Grid *)
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)
static void PPM_CartCoeff(double **, double **, int, int)
#define POLY_8(a0, a1, a2, a3, a4, a5, a6, a7, a8, x)