MFIX  2016-1
calc_s_ddot_s.f
Go to the documentation of this file.
1 !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
2 ! C
3 ! Module name: CALC_S_DDOT_S(IJK1, IJK2, FCELL, COM, M, DEL_DOT_U, C
4 ! S_DDOT_S, S_dd) C
5 ! C
6 ! Purpose: Calculate del.u, S:S and S_xx, S_yy or S_zz at the C
7 ! boundary for use in frictional boundary condition C
8 ! C
9 ! C
10 ! Author: Anuj Srivastava, Princeton University Date: 4-APR-98 C
11 ! Reviewer: Date: C
12 ! C
13 ! C
14 ! Literature/Document References: C
15 ! C
16 ! Variables referenced: C
17 ! Variables modified: C
18 ! C
19 ! Local variables: C
20 ! C
21 !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
22 !
23  SUBROUTINE calc_s_ddot_s(IJK1,IJK2,FCELL,COM,M,DEL_DOT_U,S_DDOT_S,S_DD)
24 !...Translated by Pacific-Sierra Research VAST-90 2.06G5 12:17:31 12/09/98
25 !...Switches: -xf
26 !
27 !-----------------------------------------------
28 ! M o d u l e s
29 !-----------------------------------------------
30  USE param
31  USE param1
32  USE constant
33  USE fldvar
34  USE geometry
35  USE indices
36  USE compar
37  USE fun_avg
38  USE functions
39  IMPLICIT NONE
40 !-----------------------------------------------
41 ! G l o b a l P a r a m e t e r s
42 !-----------------------------------------------
43 !-----------------------------------------------
44 ! D u m m y A r g u m e n t s
45 !-----------------------------------------------
46 !
47 ! IJK indices for wall cell and fluid cell
48  INTEGER IJK1, IJK2
49 !
50 ! Other indices
51  INTEGER IPJK2, IPJMK2, IPJPK2, IPJKM2, IPJKP2
52  INTEGER IJPK2, IJPKM2, IJPKP2, IMJPK2
53  INTEGER IJKP2, IJMKP2, IMJKP2
54 
55 !
56 ! The location (e,w,n...) of fluid cell
57  CHARACTER FCELL
58 !
59 ! Velocity component (U, V, W)
60  CHARACTER COM
61 !
62 ! Solids phase index
63  INTEGER M
64 !
65 ! del.u
66  DOUBLE PRECISION DEL_DOT_U
67 !
68 ! S:S
69  DOUBLE PRECISION S_DDOT_S
70 
71 ! S_dd (dd is the relevant direction x,y or z)
72  DOUBLE PRECISION S_dd
73 !
74 ! The location where D (rate of strain tensor) is calculated in
75 ! located at the center of 4
76 ! cells - 2 fluid cells and 2 wall cells. Coordinates of this are i,j,k
77 
78 ! U_s at the north (i, j+1/2, k)
79  DOUBLE PRECISION U_s_N
80 !
81 ! U_s at the south (i, j-1/2, k)
82  DOUBLE PRECISION U_s_S
83 
84 ! U_s at the east (i+1/2, j, k)
85  DOUBLE PRECISION U_s_E
86 !
87 ! U_s at the west (i-1/2, j, k)
88  DOUBLE PRECISION U_s_W
89 !
90 ! U_s at the top (i, j, k+1/2)
91  DOUBLE PRECISION U_s_T
92 !
93 ! U_s at the bottom (i, j, k-1/2)
94  DOUBLE PRECISION U_s_B
95 !
96 ! U_s at the center (i, j, k)
97 ! Calculated for Cylindrical coordinates only.
98  DOUBLE PRECISION U_s_C
99 !
100 ! V_s at the north (i, j+1/2, k)
101  DOUBLE PRECISION V_s_N
102 !
103 ! V_s at the south (i, j-1/2, k)
104  DOUBLE PRECISION V_s_S
105 !
106 ! V_s at the east (i+1/2, j, k)
107  DOUBLE PRECISION V_s_E
108 !
109 ! V_s at the west (i-1/2, j, k)
110  DOUBLE PRECISION V_s_W
111 
112 ! V_s at the top (i, j, k+1/2)
113  DOUBLE PRECISION V_s_T
114 !
115 ! V_s at the bottom (i, j, k-1/2)
116  DOUBLE PRECISION V_s_B
117 
118 ! W_s at the north (i, j+1/2, k)
119  DOUBLE PRECISION W_s_N
120 !
121 ! W_s at the south (i, j-1/2, k)
122  DOUBLE PRECISION W_s_S
123 !
124 ! W_s at the east (i+1/2, j, k)
125  DOUBLE PRECISION W_s_E
126 !
127 ! W_s at the west (1-1/2, j, k)
128  DOUBLE PRECISION W_s_W
129 !
130 ! W_s at the top (i, j, k+1/2)
131  DOUBLE PRECISION W_s_T
132 !
133 ! W_s at the bottom (i, j, k-1/2)
134  DOUBLE PRECISION W_s_B
135 !
136 ! W_s at the center (i, j, k).
137 ! Calculated for Cylindrical coordinates only.
138  DOUBLE PRECISION W_s_C
139 !-----------------------------------------------
140 
141  SELECT CASE (trim(com))
142  CASE ('U')
143  SELECT CASE (trim(fcell))
144  CASE ('N')
145  ipjk2 = ip_of(ijk2)
146  ipjmk2 = jm_of(ipjk2)
147 !
148  u_s_n = u_s(ijk2,m)
149 !
150  u_s_s = u_s(ijk1,m)
151 !
152  u_s_e = avg_y(avg_x_e(u_s(ijk1,m),u_s(ipjmk2,m),i_of(ipjmk2)),&
153  avg_x_e(u_s(ijk2,m),u_s(ipjk2,m),i_of(ipjk2)),j_of(ipjmk2))
154 !
155  u_s_w = avg_y(avg_x_e(u_s(im_of(ijk1),m),u_s(ijk1,m),i_of(ijk1)),&
156  avg_x_e(u_s(im_of(ijk2),m),u_s(ijk2,m),i_of(ijk2)),j_of(ijk1))
157 !
158  u_s_t = avg_y(avg_z(u_s(ijk1,m),u_s(kp_of(ijk1),m),k_of(ijk1)),&
159  avg_z(u_s(ijk2,m),u_s(kp_of(ijk2),m),k_of(ijk2)),j_of(ijk1))
160 !
161  u_s_b = avg_y(avg_z(u_s(km_of(ijk1),m),u_s(ijk1,m),k_of(km_of(&
162  ijk1))),avg_z(u_s(km_of(ijk2),m),u_s(ijk2,m),k_of(km_of(ijk2))&
163  ),j_of(ijk1))
164 !
165  v_s_n = avg_x(avg_y_n(v_s(ijk1,m),v_s(ijk2,m)),avg_y_n(v_s(&
166  ipjmk2,m),v_s(ipjk2,m)),i_of(ijk2))
167 !
168  v_s_s = avg_x(avg_y_n(v_s(jm_of(ijk1),m),v_s(ijk1,m)),avg_y_n(&
169  v_s(jm_of(ipjmk2),m),v_s(ipjmk2,m)),i_of(ijk1))
170 !
171  v_s_e = zero
172 !
173  v_s_w = zero
174 !
175  v_s_t = zero
176 !
177  v_s_b = zero
178 !
179  w_s_n = avg_x(avg_z_t(w_s(km_of(ijk2),m),w_s(ijk2,m)),avg_z_t(&
180  w_s(km_of(ipjk2),m),w_s(ipjk2,m)),i_of(ijk2))
181 !
182  w_s_s = avg_x(avg_z_t(w_s(km_of(ijk1),m),w_s(ijk1,m)),avg_z_t(&
183  w_s(km_of(ipjmk2),m),w_s(ipjmk2,m)),i_of(ijk1))
184 !
185  w_s_e = avg_y(avg_z_t(w_s(km_of(ipjmk2),m),w_s(ipjmk2,m)),avg_z_t&
186  (w_s(km_of(ipjk2),m),w_s(ipjk2,m)),j_of(ipjmk2))
187 !
188  w_s_w = avg_y(avg_z_t(w_s(km_of(ijk1),m),w_s(ijk1,m)),avg_z_t(&
189  w_s(km_of(ijk2),m),w_s(ijk2,m)),j_of(ijk1))
190 !
191  w_s_t = avg_x(avg_y(w_s(ijk1,m),w_s(ijk2,m),j_of(ijk1)),avg_y(&
192  w_s(ipjmk2,m),w_s(ipjk2,m),j_of(ipjmk2)),i_of(ijk1))
193 !
194  w_s_b = avg_x(avg_y(w_s(km_of(ijk1),m),w_s(km_of(ijk2),m),j_of(&
195  km_of(ijk1))),avg_y(w_s(km_of(ipjmk2),m),w_s(km_of(ipjk2),m),&
196  j_of(km_of(ipjmk2))),i_of(ijk1))
197 !
198  IF (cylindrical) THEN
199  u_s_c = avg_y(u_s(ijk1,m),u_s(ijk2,m),j_of(ijk1))
200  w_s_c = avg_x(w_s_w,w_s_e,i_of(ijk1))
201  ELSE
202  u_s_c = zero
203  w_s_c = zero
204  ENDIF
205 !
206  CALL sddots (ijk1, fcell, 'Y', u_s_n, u_s_s, u_s_e, u_s_w, u_s_t, &
207  u_s_b, u_s_c, v_s_n, v_s_s, v_s_e, v_s_w, v_s_t, v_s_b, w_s_n, &
208  w_s_s, w_s_e, w_s_w, w_s_t, w_s_b, w_s_c, del_dot_u, s_ddot_s, &
209  s_dd)
210 !
211 !
212 !
213 !
214  CASE ('S')
215  ipjk2 = ip_of(ijk2)
216  ipjpk2 = jp_of(ipjk2)
217 !
218  u_s_n = u_s(ijk1,m)
219 !
220  u_s_s = u_s(ijk2,m)
221 !
222  u_s_e = avg_y(avg_x_e(u_s(ijk2,m),u_s(ipjk2,m),i_of(ipjk2)),&
223  avg_x_e(u_s(ijk1,m),u_s(ipjpk2,m),i_of(ipjpk2)),j_of(ipjk2))
224 !
225  u_s_w = avg_y(avg_x_e(u_s(im_of(ijk2),m),u_s(ijk2,m),i_of(ijk2)),&
226  avg_x_e(u_s(im_of(ijk1),m),u_s(ijk1,m),i_of(ijk1)),j_of(ijk2))
227 !
228  u_s_t = avg_y(avg_z(u_s(ijk2,m),u_s(kp_of(ijk2),m),k_of(ijk2)),&
229  avg_z(u_s(ijk1,m),u_s(kp_of(ijk1),m),k_of(ijk1)),j_of(ijk2))
230 !
231  u_s_b = avg_y(avg_z(u_s(km_of(ijk2),m),u_s(ijk2,m),k_of(km_of(&
232  ijk2))),avg_z(u_s(km_of(ijk1),m),u_s(ijk1,m),k_of(km_of(ijk1))&
233  ),j_of(ijk2))
234 !
235  v_s_n = avg_x(avg_y_n(v_s(ijk2,m),v_s(ijk1,m)),avg_y_n(v_s(ipjk2&
236  ,m),v_s(ipjpk2,m)),i_of(ijk1))
237 !
238  v_s_s = avg_x(avg_y_n(v_s(jm_of(ijk2),m),v_s(ijk2,m)),avg_y_n(&
239  v_s(jm_of(ipjk2),m),v_s(ipjk2,m)),i_of(ijk2))
240 !
241  v_s_e = zero
242 !
243  v_s_w = zero
244 !
245  v_s_t = zero
246 !
247  v_s_b = zero
248 !
249  w_s_n = avg_x(avg_z_t(w_s(km_of(ijk1),m),w_s(ijk1,m)),avg_z_t(&
250  w_s(km_of(ipjpk2),m),w_s(ipjpk2,m)),i_of(ijk1))
251 !
252  w_s_s = avg_x(avg_z_t(w_s(km_of(ijk2),m),w_s(ijk2,m)),avg_z_t(&
253  w_s(km_of(ipjk2),m),w_s(ipjk2,m)),i_of(ijk2))
254 !
255  w_s_e = avg_y(avg_z_t(w_s(km_of(ipjk2),m),w_s(ipjk2,m)),avg_z_t(&
256  w_s(km_of(ipjpk2),m),w_s(ipjpk2,m)),j_of(ipjk2))
257 !
258  w_s_w = avg_y(avg_z_t(w_s(km_of(ijk2),m),w_s(ijk2,m)),avg_z_t(&
259  w_s(km_of(ijk1),m),w_s(ijk1,m)),j_of(ijk2))
260 !
261  w_s_t = avg_x(avg_y(w_s(ijk2,m),w_s(ijk1,m),j_of(ijk2)),avg_y(&
262  w_s(ipjk2,m),w_s(ipjpk2,m),j_of(ipjk2)),i_of(ijk2))
263 !
264  w_s_b = avg_x(avg_y(w_s(km_of(ijk2),m),w_s(km_of(ijk1),m),j_of(&
265  km_of(ijk2))),avg_y(w_s(km_of(ipjk2),m),w_s(km_of(ipjpk2),m),&
266  j_of(km_of(ipjk2))),i_of(ijk2))
267 !
268  IF (cylindrical) THEN
269  u_s_c = avg_y(u_s(ijk2,m),u_s(ijk1,m),j_of(ijk2))
270  w_s_c = avg_x(w_s_w,w_s_e,i_of(ijk2))
271  ELSE
272  u_s_c = zero
273  w_s_c = zero
274  ENDIF
275 !
276  CALL sddots (ijk2, fcell, 'Y', u_s_n, u_s_s, u_s_e, u_s_w, u_s_t, &
277  u_s_b, u_s_c, v_s_n, v_s_s, v_s_e, v_s_w, v_s_t, v_s_b, w_s_n, &
278  w_s_s, w_s_e, w_s_w, w_s_t, w_s_b, w_s_c, del_dot_u, s_ddot_s, &
279  s_dd)
280 !
281 !
282  CASE ('T')
283  ipjk2 = ip_of(ijk2)
284  ipjkm2 = km_of(ipjk2)
285 !
286 !
287  u_s_n = avg_z(avg_y(u_s(ijk1,m),u_s(jp_of(ijk1),m),j_of(ijk1)),&
288  avg_y(u_s(ijk2,m),u_s(jp_of(ijk2),m),j_of(ijk2)),k_of(ijk1))
289 !
290  u_s_s = avg_z(avg_y(u_s(jm_of(ijk1),m),u_s(ijk1,m),j_of(jm_of(&
291  ijk1))),avg_y(u_s(jm_of(ijk2),m),u_s(ijk2,m),j_of(jm_of(ijk2))&
292  ),k_of(ijk1))
293 !
294  u_s_e = avg_z(avg_x_e(u_s(ijk1,m),u_s(ipjkm2,m),i_of(ipjkm2)),&
295  avg_x_e(u_s(ijk2,m),u_s(ipjk2,m),i_of(ipjk2)),k_of(ipjkm2))
296 !
297  u_s_w = avg_z(avg_x_e(u_s(im_of(ijk1),m),u_s(ijk1,m),i_of(ijk1)),&
298  avg_x_e(u_s(im_of(ijk2),m),u_s(ijk2,m),i_of(ijk2)),k_of(ijk1))
299 !
300  u_s_t = u_s(ijk2,m)
301 !
302  u_s_b = u_s(ijk1,m)
303 !
304  v_s_n = avg_x(avg_z(v_s(ijk1,m),v_s(ijk2,m),k_of(ijk1)),avg_z(&
305  v_s(ipjkm2,m),v_s(ipjk2,m),k_of(ipjkm2)),i_of(ijk1))
306 !
307  v_s_s = avg_x(avg_z(v_s(jm_of(ijk1),m),v_s(jm_of(ijk2),m),k_of(&
308  jm_of(ijk1))),avg_z(v_s(jm_of(ipjkm2),m),v_s(jm_of(ipjk2),m),&
309  k_of(jm_of(ipjkm2))),i_of(ijk1))
310 !
311  v_s_e = avg_z(avg_y_n(v_s(jm_of(ipjkm2),m),v_s(ipjkm2,m)),avg_y_n&
312  (v_s(jm_of(ipjk2),m),v_s(ipjk2,m)),k_of(ipjkm2))
313 !
314  v_s_w = avg_z(avg_y_n(v_s(jm_of(ijk1),m),v_s(ijk1,m)),avg_y_n(&
315  v_s(jm_of(ijk2),m),v_s(ijk2,m)),k_of(ijk1))
316 !
317  v_s_t = avg_x(avg_y_n(v_s(jm_of(ijk2),m),v_s(ijk2,m)),avg_y_n(&
318  v_s(jm_of(ipjk2),m),v_s(ipjk2,m)),i_of(ijk2))
319 !
320  v_s_b = avg_x(avg_y_n(v_s(jm_of(ijk1),m),v_s(ijk1,m)),avg_y_n(&
321  v_s(jm_of(ipjkm2),m),v_s(ipjkm2,m)),i_of(ijk1))
322 !
323  w_s_n = zero
324 !
325  w_s_s = zero
326 !
327  w_s_e = zero
328 !
329  w_s_w = zero
330 !
331  w_s_t = avg_x(avg_z_t(w_s(ijk1,m),w_s(ijk2,m)),avg_z_t(w_s(&
332  ipjkm2,m),w_s(ipjk2,m)),i_of(ijk2))
333 !
334  w_s_b = avg_x(avg_z_t(w_s(km_of(ijk1),m),w_s(ijk1,m)),avg_z_t(&
335  w_s(km_of(ipjkm2),m),w_s(ipjkm2,m)),i_of(ijk1))
336 !
337  IF (cylindrical) THEN
338  u_s_c = avg_z(u_s(ijk1,m),u_s(ijk2,m),k_of(ijk1))
339  w_s_c = avg_x(w_s(ijk1,m),w_s(ipjkm2,m),i_of(ijk1))
340  ELSE
341  u_s_c = zero
342  w_s_c = zero
343  ENDIF
344 !
345  CALL sddots (ijk1, fcell, 'Y', u_s_n, u_s_s, u_s_e, u_s_w, u_s_t, &
346  u_s_b, u_s_c, v_s_n, v_s_s, v_s_e, v_s_w, v_s_t, v_s_b, w_s_n, &
347  w_s_s, w_s_e, w_s_w, w_s_t, w_s_b, w_s_c, del_dot_u, s_ddot_s, &
348  s_dd)
349 !
350 !
351  CASE ('B')
352  ipjk2 = ip_of(ijk2)
353  ipjkp2 = kp_of(ipjk2)
354 !
355  u_s_n = avg_z(avg_y(u_s(ijk2,m),u_s(jp_of(ijk2),m),j_of(ijk2)),&
356  avg_y(u_s(ijk1,m),u_s(jp_of(ijk1),m),j_of(ijk1)),k_of(ijk2))
357 !
358  u_s_s = avg_z(avg_y(u_s(jm_of(ijk2),m),u_s(ijk2,m),j_of(jm_of(&
359  ijk2))),avg_y(u_s(jm_of(ijk1),m),u_s(ijk1,m),j_of(jm_of(ijk1))&
360  ),k_of(ijk2))
361 !
362  u_s_e = avg_z(avg_x_e(u_s(ijk2,m),u_s(ipjk2,m),i_of(ipjk2)),&
363  avg_x_e(u_s(ijk1,m),u_s(ipjkp2,m),i_of(ipjkp2)),k_of(ipjk2))
364 !
365  u_s_w = avg_z(avg_x_e(u_s(im_of(ijk2),m),u_s(ijk2,m),i_of(ijk2)),&
366  avg_x_e(u_s(im_of(ijk1),m),u_s(ijk1,m),i_of(ijk1)),k_of(ijk2))
367 !
368  u_s_t = u_s(ijk1,m)
369 !
370  u_s_b = u_s(ijk2,m)
371 !
372  v_s_n = avg_x(avg_z(v_s(ijk2,m),v_s(ijk1,m),k_of(ijk2)),avg_z(&
373  v_s(ipjk2,m),v_s(ipjkp2,m),k_of(ipjk2)),i_of(ijk2))
374 !
375  v_s_s = avg_x(avg_z(v_s(jm_of(ijk2),m),v_s(jm_of(ijk1),m),k_of(&
376  jm_of(ijk2))),avg_z(v_s(jm_of(ipjk2),m),v_s(jm_of(ipjkp2),m),&
377  k_of(jm_of(ipjk2))),i_of(ijk2))
378 !
379  v_s_e = avg_z(avg_y_n(v_s(jm_of(ipjk2),m),v_s(ipjk2,m)),avg_y_n(&
380  v_s(jm_of(ipjkp2),m),v_s(ipjkp2,m)),k_of(ipjk2))
381 !
382  v_s_w = avg_z(avg_y_n(v_s(jm_of(ijk2),m),v_s(ijk2,m)),avg_y_n(&
383  v_s(jm_of(ijk1),m),v_s(ijk1,m)),k_of(ijk2))
384 !
385  v_s_t = avg_x(avg_y_n(v_s(jm_of(ijk1),m),v_s(ijk1,m)),avg_y_n(&
386  v_s(jm_of(ipjkp2),m),v_s(ipjkp2,m)),i_of(ijk1))
387 !
388  v_s_b = avg_x(avg_y_n(v_s(jm_of(ijk2),m),v_s(ijk2,m)),avg_y_n(&
389  v_s(jm_of(ipjk2),m),v_s(ipjk2,m)),i_of(ijk2))
390 !
391  w_s_n = zero
392 !
393  w_s_s = zero
394 !
395  w_s_e = zero
396 !
397  w_s_w = zero
398 !
399  w_s_t = avg_x(avg_z_t(w_s(ijk2,m),w_s(ijk1,m)),avg_z_t(w_s(ipjk2&
400  ,m),w_s(ipjkp2,m)),i_of(ijk1))
401 !
402  w_s_b = avg_x(avg_z_t(w_s(km_of(ijk2),m),w_s(ijk2,m)),avg_z_t(&
403  w_s(km_of(ipjk2),m),w_s(ipjk2,m)),i_of(ijk2))
404 !
405  IF (cylindrical) THEN
406  u_s_c = avg_z(u_s(ijk2,m),u_s(ijk1,m),k_of(ijk2))
407  w_s_c = avg_x(w_s(ijk2,m),w_s(ipjk2,m),i_of(ijk2))
408  ELSE
409  u_s_c = zero
410  w_s_c = zero
411  ENDIF
412 !
413  CALL sddots (ijk2, fcell, 'Y', u_s_n, u_s_s, u_s_e, u_s_w, u_s_t, &
414  u_s_b, u_s_c, v_s_n, v_s_s, v_s_e, v_s_w, v_s_t, v_s_b, w_s_n, &
415  w_s_s, w_s_e, w_s_w, w_s_t, w_s_b, w_s_c, del_dot_u, s_ddot_s, &
416  s_dd)
417 !
418  END SELECT
419 !
420  CASE ('V')
421  SELECT CASE (trim(fcell))
422  CASE ('T')
423  ijpk2 = jp_of(ijk2)
424  ijpkm2 = km_of(ijpk2)
425 !
426  u_s_n = avg_z(avg_x_e(u_s(im_of(ijpkm2),m),u_s(ijpkm2,m),i_of(&
427  ijpkm2)),avg_x_e(u_s(im_of(ijpk2),m),u_s(ijpk2,m),i_of(ijpk2))&
428  ,k_of(ijpkm2))
429 !
430  u_s_s = avg_z(avg_x_e(u_s(im_of(ijk1),m),u_s(ijk1,m),i_of(ijk1)),&
431  avg_x_e(u_s(im_of(ijk2),m),u_s(ijk2,m),i_of(ijk2)),k_of(ijk1))
432 !
433  u_s_e = avg_z(avg_y(u_s(ijk1,m),u_s(ijpkm2,m),j_of(ijk1)),avg_y(&
434  u_s(ijk2,m),u_s(ijpk2,m),j_of(ijk2)),k_of(ijk1))
435 !
436  u_s_w = avg_z(avg_y(u_s(im_of(ijk1),m),u_s(im_of(ijpkm2),m),j_of(&
437  im_of(ijk1))),avg_y(u_s(im_of(ijk2),m),u_s(im_of(ijpk2),m),&
438  j_of(im_of(ijk2))),k_of(ijk1))
439 !
440  u_s_t = avg_y(avg_x_e(u_s(im_of(ijk2),m),u_s(ijk2,m),i_of(ijk2)),&
441  avg_x_e(u_s(im_of(ijpk2),m),u_s(ijpk2,m),i_of(ijpk2)),j_of(&
442  ijk2))
443 !
444  u_s_b = avg_y(avg_x_e(u_s(im_of(ijk1),m),u_s(ijk1,m),i_of(ijk1)),&
445  avg_x_e(u_s(im_of(ijpkm2),m),u_s(ijpkm2,m),i_of(ijpkm2)),j_of(&
446  ijk1))
447 !
448  v_s_n = avg_z(avg_y_n(v_s(ijk1,m),v_s(ijpkm2,m)),avg_y_n(v_s(&
449  ijk2,m),v_s(ijpk2,m)),k_of(ijpkm2))
450 !
451  v_s_s = avg_z(avg_y_n(v_s(jm_of(ijk1),m),v_s(ijk1,m)),avg_y_n(&
452  v_s(jm_of(ijk2),m),v_s(ijk2,m)),k_of(ijk1))
453 !
454  v_s_e = avg_z(avg_x(v_s(ijk1,m),v_s(ip_of(ijk1),m),i_of(ijk1)),&
455  avg_x(v_s(ijk2,m),v_s(ip_of(ijk2),m),i_of(ijk2)),k_of(ijk1))
456 !
457  v_s_w = avg_z(avg_x(v_s(im_of(ijk1),m),v_s(ijk1,m),i_of(im_of(&
458  ijk1))),avg_x(v_s(im_of(ijk2),m),v_s(ijk2,m),i_of(im_of(ijk2))&
459  ),k_of(ijk1))
460 !
461  v_s_t = v_s(ijk2,m)
462 !
463  v_s_b = v_s(ijk1,m)
464 !
465  w_s_n = zero
466 !
467  w_s_s = zero
468 !
469  w_s_e = zero
470 !
471  w_s_w = zero
472 !
473  w_s_t = avg_y(avg_z_t(w_s(ijk1,m),w_s(ijk2,m)),avg_z_t(w_s(&
474  ijpkm2,m),w_s(ijpk2,m)),j_of(ijk2))
475 !
476  w_s_b = avg_y(avg_z_t(w_s(km_of(ijk1),m),w_s(ijk1,m)),avg_z_t(&
477  w_s(km_of(ijpkm2),m),w_s(ijpkm2,m)),j_of(ijk1))
478 !
479  IF (cylindrical) THEN
480  u_s_c = avg_z(u_s_b,u_s_t,k_of(ijk1))
481  w_s_c = avg_y(w_s(ijk1,m),w_s(ijpkm2,m),j_of(ijk1))
482  ELSE
483  u_s_c = zero
484  w_s_c = zero
485  ENDIF
486 !
487  CALL sddots (ijk1, fcell, 'N', u_s_n, u_s_s, u_s_e, u_s_w, u_s_t, &
488  u_s_b, u_s_c, v_s_n, v_s_s, v_s_e, v_s_w, v_s_t, v_s_b, w_s_n, &
489  w_s_s, w_s_e, w_s_w, w_s_t, w_s_b, w_s_c, del_dot_u, s_ddot_s, &
490  s_dd)
491 !
492 !
493  CASE ('B')
494  ijpk2 = jp_of(ijk2)
495  ijpkp2 = kp_of(ijpk2)
496 !
497  u_s_n = avg_z(avg_x_e(u_s(im_of(ijpk2),m),u_s(ijpk2,m),i_of(ijpk2&
498  )),avg_x_e(u_s(im_of(ijpkp2),m),u_s(ijpkp2,m),i_of(ijpkp2)),&
499  k_of(ijpk2))
500 !
501  u_s_s = avg_z(avg_x_e(u_s(im_of(ijk2),m),u_s(ijk2,m),i_of(ijk2)),&
502  avg_x_e(u_s(im_of(ijk1),m),u_s(ijk1,m),i_of(ijk1)),k_of(ijk2))
503 !
504  u_s_e = avg_z(avg_y(u_s(ijk2,m),u_s(ijpk2,m),j_of(ijk2)),avg_y(&
505  u_s(ijk1,m),u_s(ijpkp2,m),j_of(ijk1)),k_of(ijk2))
506 !
507  u_s_w = avg_z(avg_y(u_s(im_of(ijk2),m),u_s(im_of(ijpk2),m),j_of(&
508  im_of(ijk2))),avg_y(u_s(im_of(ijk1),m),u_s(im_of(ijpkp2),m),&
509  j_of(im_of(ijk1))),k_of(ijk2))
510 !
511  u_s_t = avg_y(avg_x_e(u_s(im_of(ijk1),m),u_s(ijk1,m),i_of(ijk1)),&
512  avg_x_e(u_s(im_of(ijpkp2),m),u_s(ijpkp2,m),i_of(ijpkp2)),j_of(&
513  ijk1))
514 !
515  u_s_b = avg_y(avg_x_e(u_s(im_of(ijk2),m),u_s(ijk2,m),i_of(ijk2)),&
516  avg_x_e(u_s(im_of(ijpk2),m),u_s(ijpk2,m),i_of(ijpk2)),j_of(&
517  ijk2))
518 !
519  v_s_n = avg_z(avg_y_n(v_s(ijk2,m),v_s(ijpk2,m)),avg_y_n(v_s(ijk1&
520  ,m),v_s(ijpkp2,m)),k_of(ijpk2))
521 !
522  v_s_s = avg_z(avg_y_n(v_s(jm_of(ijk2),m),v_s(ijk2,m)),avg_y_n(&
523  v_s(jm_of(ijk1),m),v_s(ijk1,m)),k_of(ijk2))
524 !
525  v_s_e = avg_z(avg_x(v_s(ijk2,m),v_s(ip_of(ijk2),m),i_of(ijk2)),&
526  avg_x(v_s(ijk1,m),v_s(ip_of(ijk1),m),i_of(ijk1)),k_of(ijk2))
527 !
528  v_s_w = avg_z(avg_x(v_s(im_of(ijk2),m),v_s(ijk2,m),i_of(im_of(&
529  ijk2))),avg_x(v_s(im_of(ijk1),m),v_s(ijk1,m),i_of(im_of(ijk1))&
530  ),k_of(ijk2))
531 !
532  v_s_t = v_s(ijk1,m)
533 !
534  v_s_b = v_s(ijk2,m)
535 !
536  w_s_n = zero
537 !
538  w_s_s = zero
539 !
540  w_s_e = zero
541 !
542  w_s_w = zero
543 !
544  w_s_t = avg_y(avg_z_t(w_s(ijk2,m),w_s(ijk1,m)),avg_z_t(w_s(ijpk2&
545  ,m),w_s(ijpkp2,m)),j_of(ijk1))
546 !
547  w_s_b = avg_y(avg_z_t(w_s(km_of(ijk2),m),w_s(ijk2,m)),avg_z_t(&
548  w_s(km_of(ijpk2),m),w_s(ijpk2,m)),j_of(ijk2))
549 !
550  IF (cylindrical) THEN
551  u_s_c = avg_z(u_s_b,u_s_t,k_of(ijk2))
552  w_s_c = avg_y(w_s(ijk2,m),w_s(ijpk2,m),j_of(ijk2))
553  ELSE
554  u_s_c = zero
555  w_s_c = zero
556  ENDIF
557 !
558  CALL sddots (ijk2, fcell, 'N', u_s_n, u_s_s, u_s_e, u_s_w, u_s_t, &
559  u_s_b, u_s_c, v_s_n, v_s_s, v_s_e, v_s_w, v_s_t, v_s_b, w_s_n, &
560  w_s_s, w_s_e, w_s_w, w_s_t, w_s_b, w_s_c, del_dot_u, s_ddot_s, &
561  s_dd)
562 !
563 !
564  CASE ('E')
565  ijpk2 = jp_of(ijk2)
566  imjpk2 = im_of(ijpk2)
567 !
568  u_s_n = zero
569 !
570  u_s_s = zero
571 !
572  u_s_e = avg_y(avg_x_e(u_s(ijk1,m),u_s(ijk2,m),i_of(ijk2)),avg_x_e&
573  (u_s(imjpk2,m),u_s(ijpk2,m),i_of(ijpk2)),j_of(ijk2))
574 !
575  u_s_w = avg_y(avg_x_e(u_s(im_of(ijk1),m),u_s(ijk1,m),i_of(ijk1)),&
576  avg_x_e(u_s(im_of(imjpk2),m),u_s(imjpk2,m),i_of(imjpk2)),j_of(&
577  ijk1))
578 !
579  u_s_t = zero
580 !
581  u_s_b = zero
582 !
583  v_s_n = avg_x(avg_y_n(v_s(ijk1,m),v_s(imjpk2,m)),avg_y_n(v_s(&
584  ijk2,m),v_s(ijpk2,m)),i_of(imjpk2))
585 !
586  v_s_s = avg_x(avg_y_n(v_s(jm_of(ijk1),m),v_s(ijk1,m)),avg_y_n(&
587  v_s(jm_of(ijk2),m),v_s(ijk2,m)),i_of(ijk1))
588 !
589  v_s_e = v_s(ijk2,m)
590 !
591  v_s_w = v_s(ijk1,m)
592 !
593  v_s_t = avg_x(avg_z(v_s(ijk1,m),v_s(kp_of(ijk1),m),k_of(ijk1)),&
594  avg_z(v_s(ijk2,m),v_s(kp_of(ijk2),m),k_of(ijk2)),i_of(ijk1))
595 !
596  v_s_b = avg_x(avg_z(v_s(km_of(ijk1),m),v_s(ijk1,m),k_of(km_of(&
597  ijk1))),avg_z(v_s(km_of(ijk2),m),v_s(ijk2,m),k_of(km_of(ijk2))&
598  ),i_of(ijk1))
599 !
600  w_s_n = avg_x(avg_z_t(w_s(km_of(imjpk2),m),w_s(imjpk2,m)),avg_z_t&
601  (w_s(km_of(ijpk2),m),w_s(ijpk2,m)),i_of(imjpk2))
602 !
603  w_s_s = avg_x(avg_z_t(w_s(km_of(ijk1),m),w_s(ijk1,m)),avg_z_t(&
604  w_s(km_of(ijk2),m),w_s(ijk2,m)),i_of(ijk1))
605 !
606  w_s_e = avg_y(avg_z_t(w_s(km_of(ijk2),m),w_s(ijk2,m)),avg_z_t(&
607  w_s(km_of(ijpk2),m),w_s(ijpk2,m)),j_of(ijk2))
608 !
609  w_s_w = avg_y(avg_z_t(w_s(km_of(ijk1),m),w_s(ijk1,m)),avg_z_t(&
610  w_s(km_of(imjpk2),m),w_s(imjpk2,m)),j_of(ijk1))
611 !
612  w_s_t = avg_x(avg_y(w_s(ijk1,m),w_s(imjpk2,m),j_of(ijk1)),avg_y(&
613  w_s(ijk2,m),w_s(ijpk2,m),j_of(ijk2)),i_of(ijk1))
614 !
615  w_s_b = avg_x(avg_y(w_s(km_of(ijk1),m),w_s(km_of(imjpk2),m),j_of(&
616  km_of(ijk1))),avg_y(w_s(km_of(ijk2),m),w_s(km_of(ijpk2),m),&
617  j_of(km_of(ijk2))),i_of(ijk1))
618 !
619  IF (cylindrical) THEN
620  u_s_c = avg_y(u_s(ijk1,m),u_s(imjpk2,m),j_of(ijk1))
621  w_s_c = avg_x(w_s_w,w_s_e,i_of(ijk1))
622  ELSE
623  u_s_c = zero
624  w_s_c = zero
625  ENDIF
626 !
627  CALL sddots (ijk1, fcell, 'Y', u_s_n, u_s_s, u_s_e, u_s_w, u_s_t, &
628  u_s_b, u_s_c, v_s_n, v_s_s, v_s_e, v_s_w, v_s_t, v_s_b, w_s_n, &
629  w_s_s, w_s_e, w_s_w, w_s_t, w_s_b, w_s_c, del_dot_u, s_ddot_s, &
630  s_dd)
631 !
632 !
633  CASE ('W')
634  ijpk2 = jp_of(ijk2)
635  ipjpk2 = ip_of(ijpk2)
636 !
637  u_s_n = zero
638 !
639  u_s_s = zero
640 !
641  u_s_e = avg_y(avg_x_e(u_s(ijk2,m),u_s(ijk1,m),i_of(ijk1)),avg_x_e&
642  (u_s(ijpk2,m),u_s(ipjpk2,m),i_of(ipjpk2)),j_of(ijk1))
643 !
644  u_s_w = avg_y(avg_x_e(u_s(im_of(ijk2),m),u_s(ijk2,m),i_of(ijk2)),&
645  avg_x_e(u_s(im_of(ijpk2),m),u_s(ijpk2,m),i_of(ijpk2)),j_of(&
646  ijk2))
647 !
648  u_s_t = zero
649 !
650  u_s_b = zero
651 !
652  v_s_n = avg_x(avg_y_n(v_s(ijk2,m),v_s(ijpk2,m)),avg_y_n(v_s(ijk1&
653  ,m),v_s(ipjpk2,m)),i_of(ijpk2))
654 !
655  v_s_s = avg_x(avg_y_n(v_s(jm_of(ijk2),m),v_s(ijk2,m)),avg_y_n(&
656  v_s(jm_of(ijk1),m),v_s(ijk1,m)),i_of(ijk2))
657 !
658  v_s_e = v_s(ijk1,m)
659 !
660  v_s_w = v_s(ijk2,m)
661 !
662  v_s_t = avg_x(avg_z(v_s(ijk2,m),v_s(kp_of(ijk2),m),k_of(ijk2)),&
663  avg_z(v_s(ijk1,m),v_s(kp_of(ijk1),m),k_of(ijk1)),i_of(ijk2))
664 !
665  v_s_b = avg_x(avg_z(v_s(km_of(ijk2),m),v_s(ijk2,m),k_of(km_of(&
666  ijk2))),avg_z(v_s(km_of(ijk1),m),v_s(ijk1,m),k_of(km_of(ijk1))&
667  ),i_of(ijk2))
668 !
669  w_s_n = avg_x(avg_z_t(w_s(km_of(ijpk2),m),w_s(ijpk2,m)),avg_z_t(&
670  w_s(km_of(ipjpk2),m),w_s(ipjpk2,m)),i_of(ijpk2))
671 !
672  w_s_s = avg_x(avg_z_t(w_s(km_of(ijk2),m),w_s(ijk2,m)),avg_z_t(&
673  w_s(km_of(ijk1),m),w_s(ijk1,m)),i_of(ijk2))
674 !
675  w_s_e = avg_y(avg_z_t(w_s(km_of(ijk1),m),w_s(ijk1,m)),avg_z_t(&
676  w_s(km_of(ipjpk2),m),w_s(ipjpk2,m)),j_of(ijk1))
677 !
678  w_s_w = avg_y(avg_z_t(w_s(km_of(ijk2),m),w_s(ijk2,m)),avg_z_t(&
679  w_s(km_of(ijpk2),m),w_s(ijpk2,m)),j_of(ijk2))
680 !
681  w_s_t = avg_x(avg_y(w_s(ijk2,m),w_s(ijpk2,m),j_of(ijk2)),avg_y(&
682  w_s(ijk1,m),w_s(ipjpk2,m),j_of(ijk1)),i_of(ijk2))
683 !
684  w_s_b = avg_x(avg_y(w_s(km_of(ijk2),m),w_s(km_of(ijpk2),m),j_of(&
685  km_of(ijk2))),avg_y(w_s(km_of(ijk1),m),w_s(km_of(ipjpk2),m),&
686  j_of(km_of(ijk1))),i_of(ijk2))
687 !
688  IF (cylindrical) THEN
689  u_s_c = avg_y(u_s(ijk2,m),u_s(ijpk2,m),j_of(ijk2))
690  w_s_c = avg_x(w_s_w,w_s_e,i_of(ijk2))
691  ELSE
692  u_s_c = zero
693  w_s_c = zero
694  ENDIF
695 !
696  CALL sddots (ijk2, fcell, 'Y', u_s_n, u_s_s, u_s_e, u_s_w, u_s_t, &
697  u_s_b, u_s_c, v_s_n, v_s_s, v_s_e, v_s_w, v_s_t, v_s_b, w_s_n, &
698  w_s_s, w_s_e, w_s_w, w_s_t, w_s_b, w_s_c, del_dot_u, s_ddot_s, &
699  s_dd)
700 !
701  END SELECT
702 !
703  CASE ('W')
704  SELECT CASE (trim(fcell))
705  CASE ('N')
706  ijkp2 = kp_of(ijk2)
707  ijmkp2 = jm_of(ijkp2)
708 !
709  u_s_n = avg_z(avg_x_e(u_s(im_of(ijk2),m),u_s(ijk2,m),i_of(ijk2)),&
710  avg_x_e(u_s(im_of(ijkp2),m),u_s(ijkp2,m),i_of(ijkp2)),k_of(&
711  ijk2))
712 !
713  u_s_s = avg_z(avg_x_e(u_s(im_of(ijk1),m),u_s(ijk1,m),i_of(ijk1)),&
714  avg_x_e(u_s(im_of(ijmkp2),m),u_s(ijmkp2,m),i_of(ijmkp2)),k_of(&
715  ijk1))
716 !
717  u_s_e = avg_z(avg_y(u_s(ijk1,m),u_s(ijk2,m),j_of(ijk1)),avg_y(&
718  u_s(ijmkp2,m),u_s(ijkp2,m),j_of(ijmkp2)),k_of(ijk1))
719 !
720  u_s_w = avg_z(avg_y(u_s(im_of(ijk1),m),u_s(im_of(ijk2),m),j_of(&
721  im_of(ijk1))),avg_y(u_s(im_of(ijmkp2),m),u_s(im_of(ijkp2),m),&
722  j_of(im_of(ijmkp2))),k_of(ijk1))
723 !
724  u_s_t = avg_y(avg_x_e(u_s(im_of(ijmkp2),m),u_s(ijmkp2,m),i_of(&
725  ijmkp2)),avg_x_e(u_s(im_of(ijkp2),m),u_s(ijkp2,m),i_of(ijkp2))&
726  ,j_of(ijmkp2))
727 !
728  u_s_b = avg_y(avg_x_e(u_s(im_of(ijk1),m),u_s(ijk1,m),i_of(ijk1)),&
729  avg_x_e(u_s(im_of(ijk2),m),u_s(ijk2,m),i_of(ijk2)),j_of(ijk1))
730 !
731  v_s_n = avg_z(avg_y_n(v_s(ijk1,m),v_s(ijk2,m)),avg_y_n(v_s(&
732  ijmkp2,m),v_s(ijkp2,m)),k_of(ijk2))
733 !
734  v_s_s = avg_z(avg_y_n(v_s(jm_of(ijk1),m),v_s(ijk1,m)),avg_y_n(&
735  v_s(jm_of(ijmkp2),m),v_s(ijmkp2,m)),k_of(ijk1))
736 !
737  v_s_e = zero
738 !
739  v_s_w = zero
740 !
741  v_s_t = zero
742 !
743  v_s_b = zero
744 !
745  w_s_n = w_s(ijk2,m)
746 !
747  w_s_s = w_s(ijk1,m)
748 !
749  w_s_e = avg_y(avg_x(w_s(ijk1,m),w_s(ip_of(ijk1),m),i_of(ijk1)),&
750  avg_x(w_s(ijk2,m),w_s(ip_of(ijk2),m),i_of(ijk2)),j_of(ijk1))
751 !
752  w_s_w = avg_y(avg_x(w_s(im_of(ijk1),m),w_s(ijk1,m),i_of(im_of(&
753  ijk1))),avg_x(w_s(im_of(ijk2),m),w_s(ijk2,m),i_of(im_of(ijk2))&
754  ),j_of(ijk1))
755 !
756  w_s_t = avg_y(avg_z_t(w_s(ijk1,m),w_s(ijmkp2,m)),avg_z_t(w_s(&
757  ijk2,m),w_s(ijkp2,m)),j_of(ijmkp2))
758 !
759  w_s_b = avg_y(avg_z_t(w_s(km_of(ijk1),m),w_s(ijk1,m)),avg_z_t(&
760  w_s(km_of(ijk2),m),w_s(ijk2,m)),j_of(ijk1))
761 !
762  IF (cylindrical) THEN
763  u_s_c = avg_z(u_s_b,u_s_t,k_of(ijk1))
764  w_s_c = avg_y(w_s(ijk1,m),w_s(ijk2,m),j_of(ijk1))
765  ELSE
766  u_s_c = zero
767  w_s_c = zero
768  ENDIF
769 !
770  CALL sddots (ijk1, fcell, 'N', u_s_n, u_s_s, u_s_e, u_s_w, u_s_t, &
771  u_s_b, u_s_c, v_s_n, v_s_s, v_s_e, v_s_w, v_s_t, v_s_b, w_s_n, &
772  w_s_s, w_s_e, w_s_w, w_s_t, w_s_b, w_s_c, del_dot_u, s_ddot_s, &
773  s_dd)
774 !
775 !
776  CASE ('S')
777  ijkp2 = kp_of(ijk2)
778  ijpkp2 = jp_of(ijkp2)
779 !
780  u_s_n = avg_z(avg_x_e(u_s(im_of(ijk1),m),u_s(ijk1,m),i_of(ijk1)),&
781  avg_x_e(u_s(im_of(ijpkp2),m),u_s(ijpkp2,m),i_of(ijpkp2)),k_of(&
782  ijk1))
783 !
784  u_s_s = avg_z(avg_x_e(u_s(im_of(ijk2),m),u_s(ijk2,m),i_of(ijk2)),&
785  avg_x_e(u_s(im_of(ijkp2),m),u_s(ijkp2,m),i_of(ijkp2)),k_of(&
786  ijk2))
787 !
788  u_s_e = avg_z(avg_y(u_s(ijk2,m),u_s(ijk1,m),j_of(ijk2)),avg_y(&
789  u_s(ijkp2,m),u_s(ijpkp2,m),j_of(ijkp2)),k_of(ijk2))
790 !
791  u_s_w = avg_z(avg_y(u_s(im_of(ijk2),m),u_s(im_of(ijk1),m),j_of(&
792  im_of(ijk2))),avg_y(u_s(im_of(ijkp2),m),u_s(im_of(ijpkp2),m),&
793  j_of(im_of(ijkp2))),k_of(ijk2))
794 !
795  u_s_t = avg_y(avg_x_e(u_s(im_of(ijkp2),m),u_s(ijkp2,m),i_of(ijkp2&
796  )),avg_x_e(u_s(im_of(ijpkp2),m),u_s(ijpkp2,m),i_of(ijpkp2)),&
797  j_of(ijkp2))
798 !
799  u_s_b = avg_y(avg_x_e(u_s(im_of(ijk2),m),u_s(ijk2,m),i_of(ijk2)),&
800  avg_x_e(u_s(im_of(ijk1),m),u_s(ijk1,m),i_of(ijk1)),j_of(ijk2))
801 !
802  v_s_n = avg_z(avg_y_n(v_s(ijk2,m),v_s(ijk1,m)),avg_y_n(v_s(ijkp2&
803  ,m),v_s(ijpkp2,m)),k_of(ijk1))
804 !
805  v_s_s = avg_z(avg_y_n(v_s(jm_of(ijk2),m),v_s(ijk2,m)),avg_y_n(&
806  v_s(jm_of(ijkp2),m),v_s(ijkp2,m)),k_of(ijk2))
807 !
808  v_s_e = zero
809 !
810  v_s_w = zero
811 !
812  v_s_t = zero
813 !
814  v_s_b = zero
815 !
816  w_s_n = w_s(ijk1,m)
817 !
818  w_s_s = w_s(ijk2,m)
819 !
820  w_s_e = avg_y(avg_x(w_s(ijk2,m),w_s(ip_of(ijk2),m),i_of(ijk2)),&
821  avg_x(w_s(ijk1,m),w_s(ip_of(ijk1),m),i_of(ijk1)),j_of(ijk2))
822 !
823  w_s_w = avg_y(avg_x(w_s(im_of(ijk2),m),w_s(ijk2,m),i_of(im_of(&
824  ijk2))),avg_x(w_s(im_of(ijk1),m),w_s(ijk1,m),i_of(im_of(ijk1))&
825  ),j_of(ijk2))
826 !
827  w_s_t = avg_y(avg_z_t(w_s(ijk2,m),w_s(ijkp2,m)),avg_z_t(w_s(ijk1&
828  ,m),w_s(ijpkp2,m)),j_of(ijkp2))
829 !
830  w_s_b = avg_y(avg_z_t(w_s(km_of(ijk2),m),w_s(ijk2,m)),avg_z_t(&
831  w_s(km_of(ijk1),m),w_s(ijk1,m)),j_of(ijk2))
832 !
833  IF (cylindrical) THEN
834  u_s_c = avg_z(u_s_b,u_s_t,k_of(ijk2))
835  w_s_c = avg_y(w_s(ijk2,m),w_s(ijk1,m),j_of(ijk2))
836  ELSE
837  u_s_c = zero
838  w_s_c = zero
839  ENDIF
840 !
841  CALL sddots (ijk2, fcell, 'N', u_s_n, u_s_s, u_s_e, u_s_w, u_s_t, &
842  u_s_b, u_s_c, v_s_n, v_s_s, v_s_e, v_s_w, v_s_t, v_s_b, w_s_n, &
843  w_s_s, w_s_e, w_s_w, w_s_t, w_s_b, w_s_c, del_dot_u, s_ddot_s, &
844  s_dd)
845 !
846 !
847  CASE ('E')
848  ijkp2 = kp_of(ijk2)
849  imjkp2 = im_of(ijkp2)
850 !
851  u_s_n = zero
852 !
853  u_s_s = zero
854 !
855  u_s_e = avg_z(avg_x_e(u_s(ijk1,m),u_s(ijk2,m),i_of(ijk2)),avg_x_e&
856  (u_s(imjkp2,m),u_s(ijkp2,m),i_of(ijkp2)),k_of(ijk2))
857 !
858  u_s_w = avg_z(avg_x_e(u_s(im_of(ijk1),m),u_s(ijk1,m),i_of(ijk1)),&
859  avg_x_e(u_s(im_of(imjkp2),m),u_s(imjkp2,m),i_of(imjkp2)),k_of(&
860  ijk1))
861 !
862  u_s_t = zero
863 !
864  u_s_b = zero
865 !
866  v_s_n = avg_x(avg_z(v_s(ijk1,m),v_s(imjkp2,m),k_of(ijk1)),avg_z(&
867  v_s(ijk2,m),v_s(ijkp2,m),k_of(ijk2)),i_of(ijk1))
868 !
869  v_s_s = avg_x(avg_z(v_s(jm_of(ijk1),m),v_s(jm_of(imjkp2),m),k_of(&
870  jm_of(ijk1))),avg_z(v_s(jm_of(ijk2),m),v_s(jm_of(ijkp2),m),&
871  k_of(jm_of(ijk2))),i_of(ijk1))
872 !
873  v_s_e = avg_z(avg_y_n(v_s(jm_of(ijk2),m),v_s(ijk2,m)),avg_y_n(&
874  v_s(jm_of(ijkp2),m),v_s(ijkp2,m)),k_of(ijk2))
875 !
876  v_s_w = avg_z(avg_y_n(v_s(jm_of(ijk1),m),v_s(ijk1,m)),avg_y_n(&
877  v_s(jm_of(imjkp2),m),v_s(imjkp2,m)),k_of(ijk1))
878 !
879  v_s_t = avg_x(avg_y_n(v_s(jm_of(imjkp2),m),v_s(imjkp2,m)),avg_y_n&
880  (v_s(jm_of(ijkp2),m),v_s(ijkp2,m)),i_of(imjkp2))
881 !
882  v_s_b = avg_x(avg_y_n(v_s(jm_of(ijk1),m),v_s(ijk1,m)),avg_y_n(&
883  v_s(jm_of(ijk2),m),v_s(ijk2,m)),i_of(ijk1))
884 !
885  w_s_n = avg_x(avg_y(w_s(ijk1,m),w_s(jp_of(ijk1),m),j_of(ijk1)),&
886  avg_y(w_s(ijk2,m),w_s(jp_of(ijk2),m),j_of(ijk2)),i_of(ijk1))
887 !
888  w_s_s = avg_x(avg_y(w_s(jm_of(ijk1),m),w_s(ijk1,m),j_of(jm_of(&
889  ijk1))),avg_y(w_s(jm_of(ijk2),m),w_s(ijk2,m),j_of(jm_of(ijk2))&
890  ),i_of(ijk1))
891 !
892  w_s_e = w_s(ijk2,m)
893 !
894  w_s_w = w_s(ijk1,m)
895 !
896  w_s_t = avg_x(avg_z_t(w_s(ijk1,m),w_s(imjkp2,m)),avg_z_t(w_s(&
897  ijk2,m),w_s(ijkp2,m)),i_of(imjkp2))
898 !
899  w_s_b = avg_x(avg_z_t(w_s(km_of(ijk1),m),w_s(ijk1,m)),avg_z_t(&
900  w_s(km_of(ijk2),m),w_s(ijk2,m)),i_of(ijk1))
901 !
902  IF (cylindrical) THEN
903  u_s_c = avg_z(u_s(ijk1,m),u_s(imjkp2,m),k_of(ijk1))
904  w_s_c = avg_x(w_s(ijk1,m),w_s(ijk2,m),i_of(ijk1))
905  ELSE
906  u_s_c = zero
907  w_s_c = zero
908  ENDIF
909 !
910  CALL sddots (ijk1, fcell, 'Y', u_s_n, u_s_s, u_s_e, u_s_w, u_s_t, &
911  u_s_b, u_s_c, v_s_n, v_s_s, v_s_e, v_s_w, v_s_t, v_s_b, w_s_n, &
912  w_s_s, w_s_e, w_s_w, w_s_t, w_s_b, w_s_c, del_dot_u, s_ddot_s, &
913  s_dd)
914 !
915 !
916  CASE ('W')
917  ijkp2 = kp_of(ijk2)
918  ipjkp2 = ip_of(ijkp2)
919 !
920  u_s_n = zero
921 !
922  u_s_s = zero
923 !
924  u_s_e = avg_z(avg_x_e(u_s(ijk2,m),u_s(ijk1,m),i_of(ijk1)),avg_x_e&
925  (u_s(ijkp2,m),u_s(ipjkp2,m),i_of(ipjkp2)),k_of(ijk1))
926 !
927  u_s_w = avg_z(avg_x_e(u_s(im_of(ijk2),m),u_s(ijk2,m),i_of(ijk2)),&
928  avg_x_e(u_s(im_of(ijkp2),m),u_s(ijkp2,m),i_of(ijkp2)),k_of(&
929  ijk2))
930 !
931  u_s_t = zero
932 !
933  u_s_b = zero
934 !
935  v_s_n = avg_x(avg_z(v_s(ijk2,m),v_s(ijkp2,m),k_of(ijk2)),avg_z(&
936  v_s(ijk1,m),v_s(ipjkp2,m),k_of(ijk1)),i_of(ijk2))
937 !
938  v_s_s = avg_x(avg_z(v_s(jm_of(ijk2),m),v_s(jm_of(ijkp2),m),k_of(&
939  jm_of(ijk2))),avg_z(v_s(jm_of(ijk1),m),v_s(jm_of(ipjkp2),m),&
940  k_of(jm_of(ijk1))),i_of(ijk2))
941 !
942  v_s_e = avg_z(avg_y_n(v_s(jm_of(ijk1),m),v_s(ijk1,m)),avg_y_n(&
943  v_s(jm_of(ipjkp2),m),v_s(ipjkp2,m)),k_of(ijk1))
944 !
945  v_s_w = avg_z(avg_y_n(v_s(jm_of(ijk2),m),v_s(ijk2,m)),avg_y_n(&
946  v_s(jm_of(ijkp2),m),v_s(ijkp2,m)),k_of(ijk2))
947 !
948  v_s_t = avg_x(avg_y_n(v_s(jm_of(ijkp2),m),v_s(ijkp2,m)),avg_y_n(&
949  v_s(jm_of(ipjkp2),m),v_s(ipjkp2,m)),i_of(ijkp2))
950 !
951  v_s_b = avg_x(avg_y_n(v_s(jm_of(ijk2),m),v_s(ijk2,m)),avg_y_n(&
952  v_s(jm_of(ijk1),m),v_s(ijk1,m)),i_of(ijk2))
953 !
954  w_s_n = avg_x(avg_y(w_s(ijk2,m),w_s(jp_of(ijk2),m),j_of(ijk2)),&
955  avg_y(w_s(ijk1,m),w_s(jp_of(ijk1),m),j_of(ijk1)),i_of(ijk2))
956 !
957  w_s_s = avg_x(avg_y(w_s(jm_of(ijk2),m),w_s(ijk2,m),j_of(jm_of(&
958  ijk2))),avg_y(w_s(jm_of(ijk1),m),w_s(ijk1,m),j_of(jm_of(ijk1))&
959  ),i_of(ijk2))
960 !
961  w_s_e = w_s(ijk1,m)
962 !
963  w_s_w = w_s(ijk2,m)
964 !
965  w_s_t = avg_x(avg_z_t(w_s(ijk2,m),w_s(ijkp2,m)),avg_z_t(w_s(ijk1&
966  ,m),w_s(ipjkp2,m)),i_of(ijkp2))
967 !
968  w_s_b = avg_x(avg_z_t(w_s(km_of(ijk2),m),w_s(ijk2,m)),avg_z_t(&
969  w_s(km_of(ijk1),m),w_s(ijk1,m)),i_of(ijk2))
970 !
971  IF (cylindrical) THEN
972  u_s_c = avg_z(u_s(ijk2,m),u_s(ijkp2,m),k_of(ijk2))
973  w_s_c = avg_x(w_s(ijk2,m),w_s(ijk1,m),i_of(ijk2))
974  ELSE
975  u_s_c = zero
976  w_s_c = zero
977  ENDIF
978 !
979  CALL sddots (ijk2, fcell, 'Y', u_s_n, u_s_s, u_s_e, u_s_w, u_s_t, &
980  u_s_b, u_s_c, v_s_n, v_s_s, v_s_e, v_s_w, v_s_t, v_s_b, w_s_n, &
981  w_s_s, w_s_e, w_s_w, w_s_t, w_s_b, w_s_c, del_dot_u, s_ddot_s, &
982  s_dd)
983 !
984  END SELECT
985  END SELECT
986  RETURN
987  END SUBROUTINE calc_s_ddot_s
988 !
989 !
990 !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
991 ! C
992 ! Module name: SUBROUTINE SDDOTS(IJK, NORMAL, DZVALUE, C
993 ! U_s_N, U_s_S, U_s_E, U_s_W, U_s_T, U_s_B, U_s_C, C
994 ! V_s_N, V_s_S, V_s_E, V_s_W, V_s_T, V_s_B, C
995 ! W_s_N, W_s_S, W_s_E, W_s_W, W_s_T, W_s_B, W_s_C, C
996 ! DEL_DOT_U, S_DDOT_S) C
997 ! C
998 ! Purpose: Calculate del.U (trace of D_s), S:S and S_xx, S_yy or S_zz C
999 ! at the boundary C
1000 ! C
1001 ! C
1002 ! Author: Anuj Srivastava, Princeton University Date: 4-APR-98 C
1003 ! Reviewer: Date: C
1004 ! C
1005 ! C
1006 ! Literature/Document References: C
1007 ! C
1008 ! Variables referenced: C
1009 ! Variables modified: C
1010 ! C
1011 ! Local variables: C
1012 ! C
1013 !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
1014 !
1015  SUBROUTINE sddots(IJK, NORMAL, DZVALUE, U_S_N, U_S_S, U_S_E, U_S_W, U_S_T&
1016  , u_s_b, u_s_c, v_s_n, v_s_s, v_s_e, v_s_w, v_s_t, v_s_b, w_s_n, w_s_s&
1017  , w_s_e, w_s_w, w_s_t, w_s_b, w_s_c, del_dot_u, s_ddot_s, s_dd)
1018 !...Translated by Pacific-Sierra Research VAST-90 2.06G5 12:17:31 12/09/98
1019 !...Switches: -xf
1020 !
1021 !
1022 !-----------------------------------------------
1023 ! M o d u l e s
1024 !-----------------------------------------------
1025  USE param
1026  USE param1
1027  USE constant
1028  USE fldvar
1029  USE geometry
1030  USE indices
1031  USE compar
1032  USE fun_avg
1033  USE functions
1034  IMPLICIT NONE
1035 !-----------------------------------------------
1036 ! G l o b a l P a r a m e t e r s
1037 !-----------------------------------------------
1038 !-----------------------------------------------
1039 ! D u m m y A r g u m e n t s
1040 !-----------------------------------------------
1041 ! Strain rate tensor components for mth solids phase
1042  DOUBLE PRECISION D_s(3,3)
1043 
1044 ! direction of normal to wall
1045  CHARACTER NORMAL
1046 
1047 ! tells what value should odelta_Z should take
1048  CHARACTER DZVALUE
1049 
1050 ! The location where D is calculated in located at the center of 4
1051 ! cells - 2 fluid cells and 2 wall cells. Coordinates of this are i,j,k
1052 
1053 ! U_s at the north (i, j+1/2, k)
1054  DOUBLE PRECISION U_s_N
1055 !
1056 ! U_s at the south (i, j-1/2, k)
1057  DOUBLE PRECISION U_s_S
1058 
1059 ! U_s at the east (i+1/2, j, k)
1060  DOUBLE PRECISION U_s_E
1061 !
1062 ! U_s at the west (i-1/2, j, k)
1063  DOUBLE PRECISION U_s_W
1064 !
1065 ! U_s at the top (i, j, k+1/2)
1066  DOUBLE PRECISION U_s_T
1067 !
1068 ! U_s at the bottom (i, j, k-1/2)
1069  DOUBLE PRECISION U_s_B
1070 !
1071 ! U_s at the center (i, j, k)
1072 ! Calculated for Cylindrical coordinates only.
1073  DOUBLE PRECISION U_s_C
1074 !
1075 ! V_s at the north (i, j+1/2, k)
1076  DOUBLE PRECISION V_s_N
1077 !
1078 ! V_s at the south (i, j-1/2, k)
1079  DOUBLE PRECISION V_s_S
1080 !
1081 ! V_s at the east (i+1/2, j, k)
1082  DOUBLE PRECISION V_s_E
1083 !
1084 ! V_s at the west (i-1/2, j, k)
1085  DOUBLE PRECISION V_s_W
1086 
1087 ! V_s at the top (i, j, k+1/2)
1088  DOUBLE PRECISION V_s_T
1089 !
1090 ! V_s at the bottom (i, j, k-1/2)
1091  DOUBLE PRECISION V_s_B
1092 
1093 ! W_s at the north (i, j+1/2, k)
1094  DOUBLE PRECISION W_s_N
1095 !
1096 ! W_s at the south (i, j-1/2, k)
1097  DOUBLE PRECISION W_s_S
1098 !
1099 ! W_s at the east (i+1/2, j, k)
1100  DOUBLE PRECISION W_s_E
1101 !
1102 ! W_s at the west (1-1/2, j, k)
1103  DOUBLE PRECISION W_s_W
1104 !
1105 ! W_s at the top (i, j, k+1/2)
1106  DOUBLE PRECISION W_s_T
1107 !
1108 ! W_s at the bottom (i, j, k-1/2)
1109  DOUBLE PRECISION W_s_B
1110 !
1111 ! W_s at the center (i, j, k).
1112 ! Calculated for Cylindrical coordinates only.
1113  DOUBLE PRECISION W_s_C
1114 
1115 ! del.u
1116  DOUBLE PRECISION DEL_DOT_U
1117 
1118 ! S:S
1119  DOUBLE PRECISION S_DDOT_S
1120 
1121 ! S_dd (d is x,y, or z)
1122  DOUBLE PRECISION S_dd
1123 
1124 ! trace of D
1125  DOUBLE PRECISION TRACE_D
1126 
1127 ! trace of the square of D
1128  DOUBLE PRECISION TRACE_sD
1129 
1130 ! Local indices
1131  INTEGER IJK, I, J, K, I1, I2
1132 
1133  DOUBLE PRECISION odelta_Z
1134 
1135 !-----------------------------------------------
1136 !
1137 ! Define I, J, K
1138 ! IJK is the cell whose coordinates are i-1/2, j-1/2 and k-1/2
1139 ! (this cell actually has two of the above coordinates. The other
1140 ! coordinate is i, j or k)
1141 ! This is necessary because we are using oDX_E, oDY_N, and oDZ_T to
1142 ! calculate D_s
1143 !
1144 !
1145  i = i_of(ijk)
1146  j = j_of(ijk)
1147  k = k_of(ijk)
1148 !
1149  IF (dzvalue == 'Y') THEN
1150  odelta_z = (odz(k)+odz(k_of(ip_of(ijk))))/2d0
1151  ELSE
1152  odelta_z = odz_t(k)
1153  ENDIF
1154 !
1155 !
1156 !
1157 ! Find components of Mth solids phase continuum strain rate
1158 ! tensor, D_s, at i,j,k
1159 !
1160  d_s(1,1) = (u_s_e - u_s_w)*odx_e(i)
1161  d_s(1,2) = half*((u_s_n - u_s_s)*ody_n(j)+(v_s_e-v_s_w)*odx_e(i))
1162  d_s(1,3) = half*((w_s_e - w_s_w)*odx_e(i)+(u_s_t-u_s_b)*(ox_e(i)*odelta_z&
1163  )-w_s_c*ox_e(i))
1164  d_s(2,1) = d_s(1,2)
1165  d_s(2,2) = (v_s_n - v_s_s)*ody_n(j)
1166  d_s(2,3)=half*((v_s_t-v_s_b)*(ox_e(i)*odelta_z)+(w_s_n-w_s_s)*ody_n(j))
1167  d_s(3,1) = d_s(1,3)
1168  d_s(3,2) = d_s(2,3)
1169  d_s(3,3) = (w_s_t - w_s_b)*(ox_e(i)*odelta_z) + u_s_c*ox_e(i)
1170 !
1171 ! Calculate the trace of D_s
1172  trace_d = d_s(1,1) + d_s(2,2) + d_s(3,3)
1173 !
1174 ! Calculate trace of the square of D_s
1175  trace_sd = 0.d0 !Initialize the totalizer
1176  DO i1 = 1, 3
1177  trace_sd = trace_sd + sum(d_s(i1,:)*d_s(i1,:))
1178  i2 = 4
1179  END DO
1180  del_dot_u = trace_d
1181 !
1182  s_ddot_s = trace_sd - (trace_d*trace_d)/3.d0
1183 !
1184  s_ddot_s = dmax1(1d-10,s_ddot_s)
1185 !
1186  IF (normal=='E' .OR. normal=='W') THEN
1187  s_dd = d_s(1,1) - trace_d/3d0
1188 !
1189  ELSE IF (normal=='N' .OR. normal=='S') THEN
1190  s_dd = d_s(2,2) - trace_d/3d0
1191 !
1192  ELSE
1193  s_dd = d_s(3,3) - trace_d/3d0
1194  ENDIF
1195 !
1196  RETURN
1197  END SUBROUTINE sddots
fldvar::v_s
double precision, dimension(:,:), allocatable v_s
Definition: fldvar_mod.f:105
indices::i_of
integer, dimension(:), allocatable i_of
Definition: indices_mod.f:45
geometry::ox_e
double precision, dimension(:), allocatable ox_e
Definition: geometry_mod.f:143
sddots
subroutine sddots(IJK, NORMAL, DZVALUE, U_S_N, U_S_S, U_S_E, U_S_W
Definition: calc_s_ddot_s.f:1016
fldvar::w_s
double precision, dimension(:,:), allocatable w_s
Definition: fldvar_mod.f:117
fldvar::u_s
double precision, dimension(:,:), allocatable u_s
Definition: fldvar_mod.f:93
indices::k_of
integer, dimension(:), allocatable k_of
Definition: indices_mod.f:47
geometry::ody_n
double precision, dimension(:), allocatable ody_n
Definition: geometry_mod.f:123
fun_avg
Definition: fun_avg_mod.f:1
calc_s_ddot_s
subroutine calc_s_ddot_s(IJK1, IJK2, FCELL, COM, M, DEL_DOT_U, S_DDOT_S,
Definition: calc_s_ddot_s.f:24
indices::j_of
integer, dimension(:), allocatable j_of
Definition: indices_mod.f:46
geometry::odx_e
double precision, dimension(:), allocatable odx_e
Definition: geometry_mod.f:121
param1::half
double precision, parameter half
Definition: param1_mod.f:28
param
Definition: param_mod.f:2
geometry::odz
double precision, dimension(:), allocatable odz
Definition: geometry_mod.f:118
geometry::cylindrical
logical cylindrical
Definition: geometry_mod.f:168
geometry::odz_t
double precision, dimension(:), allocatable odz_t
Definition: geometry_mod.f:125
param1::zero
double precision, parameter zero
Definition: param1_mod.f:27