5 int beg,
int end,
Grid *grid)
45 double *inv_dx, *inv_dy, *inv_dz;
46 double *inv_dxi, *inv_dyi, *inv_dzi;
47 double dl1, dl2, dl3, theta, r_1, s_1;
61 #if GEOMETRY == SPHERICAL
66 dl2 = dx2 = inv_dy[j]; ,
67 dl3 = dx3 = inv_dz[
k];
69 for (i = beg; i <= end; i++){
71 #if GEOMETRY == POLAR && DIMENSIONS >= 2
73 #elif GEOMETRY == SPHERICAL
79 gradT[i][0] = (T[k][j][i+1] - T[k][j][i])*dl1; ,
80 gradT[
i][1] = 0.25*( T[
k][j+1][
i] + T[
k][j+1][i+1]
81 - T[
k][j-1][
i] - T[
k][j-1][i+1])*dl2; ,
82 gradT[
i][2] = 0.25*( T[k+1][
j][
i] + T[k+1][
j][i+1]
83 - T[k-1][
j][
i] - T[k-1][
j][i+1])*dl3;
91 dl1 = dx1 = inv_dx[i]; ,
93 dl3 = dx3 = inv_dz[
k];
95 for (j = beg; j <= end; j++){
99 #elif GEOMETRY == SPHERICAL
103 dl3 = dx3*r_1/sin(theta);)
106 gradT[j][0] = 0.25*( T[k][j][i+1] + T[k][j+1][i+1]
107 - T[k][j][i-1] - T[k][j+1][i-1])*dl1; ,
108 gradT[
j][1] = (T[
k][j+1][
i] - T[
k][
j][
i])*dl2; ,
109 gradT[
j][2] = 0.25*( T[k+1][
j][
i] + T[k+1][j+1][
i]
110 - T[k-1][
j][
i] - T[k-1][j+1][
i])*dl3;
118 #if GEOMETRY == POLAR || GEOMETRY == SPHERICAL
122 #if GEOMETRY == SPHERICAL
124 s_1 = 1.0/sin(theta);
127 for (k = beg; k <= end; k++){
129 #if GEOMETRY == SPHERICAL
132 gradT[
k][0] = 0.25*( T[
k][
j][i+1] + T[k+1][
j][i+1]
133 - T[
k][
j][i-1] - T[k+1][
j][i-1])*dl1;
134 gradT[
k][1] = 0.25*( T[
k][j+1][
i] + T[k+1][j+1][
i]
135 - T[
k][j-1][
i] - T[k+1][j-1][
i])*dl2;
136 gradT[
k][2] = (T[k+1][
j][
i] - T[
k][
j][
i])*dl3;
int g_i
x1 grid index when sweeping along the x2 or x3 direction.
int g_dir
Specifies the current sweep or direction of integration.
int g_j
x2 grid index when sweeping along the x1 or x3 direction.
int g_k
x3 grid index when sweeping along the x1 or x2 direction.
D_EXPAND(tot/[n]=(double) grid[IDIR].np_int_glob;, tot/[n]=(double) grid[JDIR].np_int_glob;, tot/[n]=(double) grid[KDIR].np_int_glob;)
void GetGradient(double ***T, double **gradT, int beg, int end, Grid *grid)
double * inv_dxi
inverse of the distance between the center of two cells, inv_dxi = .