File: RELATIVE:/../../../mfix.git/model/cartesian_grid/CG_source_w_g.f

1     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
2     !                                                                      C
3     !  Module name: CG_SOURCE_W_g(A_m, B_m, IER)                           C
4     !  Purpose: Determine contribution of cut-cell to source terms         C
5     !  for W_g momentum eq.                                                C
6     !                                                                      C
7     !                                                                      C
8     !  Author: Jeff Dietiker                              Date: 01-MAY-09  C
9     !  Reviewer:                                          Date:            C
10     !                                                                      C
11     !                                                                      C
12     !  Literature/Document References:                                     C
13     !                                                                      C
14     !  Variables referenced:                                               C
15     !  Variables modified:                                                 C
16     !                                                                      C
17     !  Local variables:                                                    C
18     !                                                                      C
19     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
20     !
21           SUBROUTINE CG_SOURCE_W_G(A_M, B_M)
22     !...Translated by Pacific-Sierra Research VAST-90 2.06G5  12:17:31  12/09/98
23     !...Switches: -xf
24     !
25     !  Include param.inc file to specify parameter values
26     !
27     !-----------------------------------------------
28     !   M o d u l e s
29     !-----------------------------------------------
30           USE param
31           USE param1
32           USE parallel
33           USE matrix
34           USE scales
35           USE constant
36           USE physprop
37           USE fldvar
38           USE visc_g
39           USE rxns
40           USE run
41           USE toleranc
42           USE geometry
43           USE indices
44           USE is
45           USE tau_g
46           USE bc
47           USE compar
48           USE sendrecv
49           USE ghdtheory
50           USE drag
51           USE fun_avg
52           USE functions
53     !=======================================================================
54     ! JFD: START MODIFICATION FOR CARTESIAN GRID IMPLEMENTATION
55     !=======================================================================
56           USE cutcell
57           USE quadric
58     !=======================================================================
59     ! JFD: END MODIFICATION FOR CARTESIAN GRID IMPLEMENTATION
60     !=======================================================================
61     
62           IMPLICIT NONE
63     !-----------------------------------------------
64     !   G l o b a l   P a r a m e t e r s
65     !-----------------------------------------------
66     !-----------------------------------------------
67     !   D u m m y   A r g u m e n t s
68     !-----------------------------------------------
69     !
70     !                      Indices
71           INTEGER          I, J, K, IJK, IJKT, IMJK, IJKP, IMJKP,&
72                            IJKE, IJKW, IJKTE, IM, IPJK
73     !
74     !                      Phase index
75           INTEGER          M
76     !
77     !                      Average volume fraction
78           DOUBLE PRECISION EPGA
79     !
80     !                      Septadiagonal matrix A_m
81           DOUBLE PRECISION A_m(DIMENSION_3, -3:3, 0:DIMENSION_M)
82     !
83     !                      Vector b_m
84           DOUBLE PRECISION B_m(DIMENSION_3, 0:DIMENSION_M)
85     !
86     !=======================================================================
87     ! JFD: START MODIFICATION FOR CARTESIAN GRID IMPLEMENTATION
88     !=======================================================================
89           INTEGER :: JM,IP,JP,IJMK,IJPK,IJKC,IJKN,IJKNE,IJKS,IJKSE,IPJMK,IJKM,KM,KP,IJMKP
90           INTEGER :: IJKTN,IJKWT,IJKST
91           DOUBLE PRECISION :: We,Ww,Wn,Ws,Wt,Wb
92           DOUBLE PRECISION :: B_NOC
93           DOUBLE PRECISION :: MU_GT_E,MU_GT_W,MU_GT_N,MU_GT_S,MU_GT_T,MU_GT_B,MU_GT_CUT
94           DOUBLE PRECISION :: WW_g
95           INTEGER :: BCV
96           CHARACTER(LEN=9) :: BCT
97     !                       virtual (added) mass
98           DOUBLE PRECISION F_vir, ROP_MA, U_se, Usw, Wse, Wsw, Wsn, Wss, Wst, Wsb, Usc,Vsc,Vsn,Vss
99     ! Wall function
100           DOUBLE PRECISION :: W_F_Slip
101     !=======================================================================
102     ! JFD: END MODIFICATION FOR CARTESIAN GRID IMPLEMENTATION
103     !=======================================================================
104     !-----------------------------------------------
105     
106           IF(CG_SAFE_MODE(5)==1) RETURN
107     
108     !
109           M = 0
110           IF (.NOT.MOMENTUM_Z_EQ(0)) RETURN
111     !
112     !
113     !!!$omp  parallel do private( I, J, K, IJK, IJKT, ISV, Sdp, V0, Vpm, Vmt, Vbf, &
114     !!!$omp&  PGT, ROGA, IMJK, IJKP, IMJKP, IJKW, IJKTE, IJKTW, IM, IPJK,  &
115     !!!$omp&  CTE, CTW, SXZB, EPMUOX, VXZA, VXZB, UGT, VCOA, VCOB, IJKE,&
116     !!!$omp&  MUGA, ROPGA, EPGA, LINE)
117           DO IJK = ijkstart3, ijkend3
118              I = I_OF(IJK)
119              J = J_OF(IJK)
120              K = K_OF(IJK)
121              IJKT = TOP_OF(IJK)
122              EPGA = AVG_Z(EP_G(IJK),EP_G(IJKT),K)
123              IF (IP_AT_T(IJK)) THEN
124     !
125     !        do nothing
126     !
127     !       dilute flow
128              ELSE IF (EPGA <= DIL_EP_S) THEN
129     !
130     !        do nothing
131     !
132     !         ELSE
133              ELSEIF(INTERIOR_CELL_AT(IJK)) THEN
134     !
135                 BCV = BC_W_ID(IJK)
136     
137                 IF(BCV > 0 ) THEN
138                    BCT = BC_TYPE(BCV)
139                 ELSE
140                    BCT = 'NONE'
141                 ENDIF
142     
143                 SELECT CASE (BCT)
144                    CASE ('CG_NSW')
145                       NOC_WG = .TRUE.
146                       WW_g = ZERO
147                       MU_GT_CUT =  (VOL(IJK)*MU_GT(IJK) + VOL(IJKT)*MU_GT(IJKT))/(VOL(IJK) + VOL(IJKT))
148     
149                       IF(.NOT.K_EPSILON) THEN
150                          A_M(IJK,0,M) = A_M(IJK,0,M) - MU_GT_CUT * Area_W_CUT(IJK)/DELH_W(IJK)
151                       ELSE
152                          CALL Wall_Function(IJK,IJK,ONE/DELH_W(IJK),W_F_Slip)
153                          A_M(IJK,0,M) = A_M(IJK,0,M)  - MU_GT_CUT * Area_W_CUT(IJK)*(ONE-W_F_Slip)/DELH_W(IJK)
154                       ENDIF
155     
156                    CASE ('CG_FSW')
157                       NOC_WG = .FALSE.
158                       WW_g = ZERO
159                    CASE('CG_PSW')
160                       IF(BC_HW_G(BCV)==UNDEFINED) THEN   ! same as NSW
161                          NOC_WG = .TRUE.
162                          WW_g = BC_WW_G(BCV)
163                          MU_GT_CUT = (VOL(IJK)*MU_GT(IJK) + VOL(IJKT)*MU_GT(IJKT))/(VOL(IJK) + VOL(IJKT))
164                          A_M(IJK,0,M) = A_M(IJK,0,M) - MU_GT_CUT * Area_W_CUT(IJK)/DELH_W(IJK)
165                          B_M(IJK,M) = B_M(IJK,M) - MU_GT_CUT * WW_g * Area_W_CUT(IJK)/DELH_W(IJK)
166                       ELSEIF(BC_HW_G(BCV)==ZERO) THEN   ! same as FSW
167                          NOC_WG = .FALSE.
168                          WW_g = ZERO
169                       ELSE                              ! partial slip
170                          NOC_WG = .FALSE.
171                          WW_g = BC_WW_G(BCV)
172                          MU_GT_CUT = (VOL(IJK)*MU_GT(IJK) + VOL(IJKT)*MU_GT(IJKT))/(VOL(IJK) + VOL(IJKT))
173                          A_M(IJK,0,M) = A_M(IJK,0,M) - MU_GT_CUT * Area_W_CUT(IJK)*(BC_HW_G(BCV))
174                          B_M(IJK,M) = B_M(IJK,M) - MU_GT_CUT * WW_g * Area_W_CUT(IJK)*(BC_HW_G(BCV))
175     
176                       ENDIF
177                    CASE ('NONE')
178                       NOC_WG = .FALSE.
179                 END SELECT
180     
181                 IF(NOC_WG) THEN
182     
183                    IMJK = IM_OF(IJK)
184                    IJMK = JM_OF(IJK)
185                    IJKM = KM_OF(IJK)
186                    IPJK = IP_OF(IJK)
187                    IJPK = JP_OF(IJK)
188                    IJKP = KP_OF(IJK)
189     
190                    We = Theta_We_bar(IJK)  * W_g(IJK)  + Theta_We(IJK)  * W_g(IPJK)
191                    Ww = Theta_We_bar(IMJK) * W_g(IMJK) + Theta_We(IMJK) * W_g(IJK)
192     
193                    Wn = Theta_Wn_bar(IJK)  * W_g(IJK)  + Theta_Wn(IJK)  * W_g(IJPK)
194                    Ws = Theta_Wn_bar(IJMK) * W_g(IJMK) + Theta_Wn(IJMK) * W_g(IJK)
195     
196                    Wt = Theta_Wt_bar(IJK)  * W_g(IJK)  + Theta_Wt(IJK)  * W_g(IJKP)
197                    Wb = Theta_Wt_bar(IJKM) * W_g(IJKM) + Theta_Wt(IJKM) * W_g(IJK)
198     
199                    IJKE = EAST_OF(IJK)
200     
201                    ijkt = top_of(ijk)
202     
203                    IF (WALL_AT(IJK)) THEN
204                       IJKC = IJKT
205                    ELSE
206                       IJKC = IJK
207                    ENDIF
208     
209                    IP = IP1(I)
210                    IM = IM1(I)
211                    IJKN = NORTH_OF(IJK)
212                    IJKNE = EAST_OF(IJKN)
213     
214                    JM = JM1(J)
215                    IPJMK = IP_OF(IJMK)
216                    IJKS = SOUTH_OF(IJK)
217                    IJKSE = EAST_OF(IJKS)
218     
219                    KP = KP1(K)
220                    IJKT = TOP_OF(IJK)
221                    IJKE = EAST_OF(IJK)
222                    IJKP = KP_OF(IJK)
223                    IJKTN = NORTH_OF(IJKT)
224                    IJKTE = EAST_OF(IJKT)
225                    IJKW = WEST_OF(IJK)
226                    IJKWT = TOP_OF(IJKW)
227                    IJKS = SOUTH_OF(IJK)
228                    IJKST = TOP_OF(IJKS)
229     
230                    MU_GT_E = AVG_Z_H(AVG_X_H(MU_GT(IJKC),MU_GT(IJKE),I),&
231                              AVG_X_H(MU_GT(IJKT),MU_GT(IJKTE),I),K)
232     
233                    MU_GT_W = AVG_Z_H(AVG_X_H(MU_GT(IJKW),MU_GT(IJKC),IM),&
234                              AVG_X_H(MU_GT(IJKWT),MU_GT(IJKT),IM),K)
235     
236                    MU_GT_N = AVG_Z_H(AVG_Y_H(MU_GT(IJKC),MU_GT(IJKN),J),&
237                              AVG_Y_H(MU_GT(IJKT),MU_GT(IJKTN),J),K)
238     
239                    MU_GT_S = AVG_Z_H(AVG_Y_H(MU_GT(IJKS),MU_GT(IJKC),JM),&
240                              AVG_Y_H(MU_GT(IJKST),MU_GT(IJKT),JM),K)
241     
242                    MU_GT_T = MU_GT(IJKT)
243                    MU_GT_B = MU_GT(IJKC)
244     
245                    B_NOC =     MU_GT_E * Ayz_W(IJK)  * (We-WW_g) * NOC_W_E(IJK)  &
246                            -   MU_GT_W * Ayz_W(IMJK) * (Ww-WW_g) * NOC_W_E(IMJK) &
247                            +   MU_GT_N * Axz_W(IJK)  * (Wn-WW_g) * NOC_W_N(IJK)  &
248                            -   MU_GT_S * Axz_W(IJMK) * (Ws-WW_g) * NOC_W_N(IJMK) &
249                            +   MU_GT_T * Axy_W(IJK)  * (Wt-WW_g) * NOC_W_T(IJK)  *2.0d0&
250                            -   MU_GT_B * Axy_W(IJKM) * (Wb-WW_g) * NOC_W_T(IJKM) *2.0D0
251     
252                    B_M(IJK,M) = B_M(IJK,M)   +  B_NOC
253                 ENDIF
254     
255                 IF(CUT_W_TREATMENT_AT(IJK)) THEN
256     !
257     !!! BEGIN VIRTUAL MASS SECTION (explicit terms)
258     ! adding transient term  dWs/dt to virtual mass term
259                    F_vir = ZERO
260                    IF(Added_Mass) THEN
261     
262                       F_vir = ( (W_s(IJK,M_AM) - W_sO(IJK,M_AM)) )*ODT*VOL_W(IJK)
263     
264                       I = I_OF(IJK)
265                       J = J_OF(IJK)
266                       K = K_OF(IJK)
267     
268                       IM = I - 1
269                       JM = J - 1
270                       KM = K - 1
271     
272                       IP = I + 1
273                       JP = J + 1
274                       KP = K + 1
275     
276                       IMJK = FUNIJK(IM,J,K)
277                       IJMK = FUNIJK(I,JM,K)
278                       IPJK = FUNIJK(IP,J,K)
279                       IJPK = FUNIJK(I,JP,K)
280                       IJKP = FUNIJK(I,J,KP)
281                       IJKM = FUNIJK(I,J,KM)
282     
283                       IMJKP = KP_OF(IMJK)
284                       IJMKP = KP_OF(IJMK)
285     
286                       IJKE = EAST_OF(IJK)
287     !
288     ! defining gas-particles velocity at momentum cell faces (or scalar cell center)
289     
290     
291                       Wse = Theta_We_bar(IJK) * W_s(IJK,M_AM) + Theta_We(IJK) * W_s(IPJK,M_AM)
292                       Wsw = Theta_We_bar(IMJK) * W_s(IMJK,M_AM) + Theta_We(IMJK) * W_s(IJK,M_AM)
293     
294                       U_se = Theta_W_te(IJK) * U_s(IJK,M_AM) + Theta_W_be(IJK) * U_s(IJKP,M_AM)
295                       Usw = Theta_W_te(IMJK) * U_s(IMJK,M_AM) + Theta_W_be(IMJK) * U_s(IMJKP,M_AM)
296     
297                       Usc = (DELX_we(IJK) * Usw + DELX_ww(IJK) * U_se) / (DELX_we(IJK) + DELX_ww(IJK))
298     
299     
300                       Wsn = Theta_Wn_bar(IJK) * W_s(IJK,M_AM) + Theta_Wn(IJK) * W_s(IJPK,M_AM)
301                       Wss = Theta_Wn_bar(IJMK) * W_s(IJMK,M_AM) + Theta_Wn(IJMK) * W_s(IJK,M_AM)
302     
303                       Vsn =  Theta_W_tn(IJK)  * V_s(IJK,M_AM)  + Theta_W_bn(IJK)  * V_s(IJKP,M_AM)
304                       Vss =  Theta_W_tn(IJMK) * V_s(IJMK,M_AM) + Theta_W_bn(IJMK) * V_s(IJMKP,M_AM)
305     
306                       Vsc = (DELY_wn(IJK) * Vss + DELY_ws(IJK) * Vsn) / (DELY_wn(IJK) + DELY_ws(IJK))
307     
308     
309                       Wst = Theta_Wt_bar(IJK)  * W_s(IJK,M_AM)  + Theta_Wt(IJK)  * W_s(IJKP,M_AM)
310                       Wsb = Theta_Wt_bar(IMJK) * W_s(IMJK,M_AM) + Theta_Wt(IMJK) * W_s(IMJKP,M_AM)
311     
312     !
313     ! adding convective terms (U dW/dx + V dW/dy + W dW/dz) to virtual mass
314     
315                       F_vir = F_vir +  Usc * (Wse - Wsw)*AYZ(IJK)    + &
316                                        Vsc * (Wsn - Wss)*AXZ(IJK)  + &
317                                        W_s(IJK,M_AM)*(Wst - Wsb)*AXY(IJK)
318     
319     
320                       ROP_MA = (VOL(IJK)*ROP_g(IJK)*EP_s(IJK,M_AM) + VOL(IJKT)*ROP_g(IJKT)*EP_s(IJKT,M_AM))/(VOL(IJK) + VOL(IJKT))
321     
322                       F_vir = F_vir * Cv * ROP_MA
323     
324                       B_M(IJK,M) = B_M(IJK,M) - F_vir ! explicit part of virtual mass force
325     
326                    ENDIF
327     !
328     !!! END VIRTUAL MASS SECTION
329     
330                 ENDIF
331     
332              ENDIF
333           END DO
334     
335           RETURN
336           END SUBROUTINE CG_SOURCE_W_G
337     
338     !
339     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
340     !                                                                      C
341     !  Module name: CG_SOURCE_W_g_BC(A_m, B_m, IER)                        C
342     !  Purpose: Determine contribution of cut-cell to source terms         C
343     !  for W_g momentum eq.                                                C
344     !                                                                      C
345     !                                                                      C
346     !  Author: Jeff Dietiker                              Date: 01-MAY-09  C
347     !  Reviewer:                                          Date:            C
348     !                                                                      C
349     !                                                                      C
350     !  Literature/Document References:                                     C
351     !                                                                      C
352     !  Variables referenced:                                               C
353     !  Variables modified:                                                 C
354     !                                                                      C
355     !  Local variables:                                                    C
356     !                                                                      C
357     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
358     !
359           SUBROUTINE CG_SOURCE_W_G_BC(A_M, B_M)
360     !...Translated by Pacific-Sierra Research VAST-90 2.06G5  12:17:31  12/09/98
361     !...Switches: -xf
362     !
363     !  Include param.inc file to specify parameter values
364     !
365     !-----------------------------------------------
366     !   M o d u l e s
367     !-----------------------------------------------
368           USE param
369           USE param1
370           USE parallel
371           USE matrix
372           USE scales
373           USE constant
374           USE physprop
375           USE fldvar
376           USE visc_g
377           USE rxns
378           USE run
379           USE toleranc
380           USE geometry
381           USE indices
382           USE is
383           USE tau_g
384           USE bc
385           USE output
386           USE compar
387           USE fun_avg
388           USE functions
389     
390     !=======================================================================
391     ! JFD: START MODIFICATION FOR CARTESIAN GRID IMPLEMENTATION
392     !=======================================================================
393           USE cutcell
394           USE quadric
395     !=======================================================================
396     ! JFD: END MODIFICATION FOR CARTESIAN GRID IMPLEMENTATION
397     !=======================================================================
398     
399     
400           IMPLICIT NONE
401     !-----------------------------------------------
402     !   G l o b a l   P a r a m e t e r s
403     !-----------------------------------------------
404     !-----------------------------------------------
405     !   D u m m y   A r g u m e n t s
406     !-----------------------------------------------
407     !
408     !                      Indices
409           INTEGER          IJK, IJKB
410     !
411     !                      Solids phase
412           INTEGER          M
413     !
414     !                      Septadiagonal matrix A_m
415           DOUBLE PRECISION A_m(DIMENSION_3, -3:3, 0:DIMENSION_M)
416     !
417     !                      Vector b_m
418           DOUBLE PRECISION B_m(DIMENSION_3, 0:DIMENSION_M)
419     !
420     
421     !=======================================================================
422     ! JFD: START MODIFICATION FOR CARTESIAN GRID IMPLEMENTATION
423     !=======================================================================
424           INTEGER :: BCV
425           CHARACTER(LEN=9) :: BCT
426     
427     !-----------------------------------------------
428     
429           IF(CG_SAFE_MODE(5)==1) RETURN
430     
431     !
432           M = 0
433     
434           DO IJK = ijkstart3, ijkend3
435     
436              BCV = BC_W_ID(IJK)
437     
438              IF(BCV > 0 ) THEN
439                 BCT = BC_TYPE(BCV)
440              ELSE
441                 BCT = 'NONE'
442              ENDIF
443     
444              SELECT CASE (BCT)
445     
446                 CASE ('CG_NSW')
447     
448                    IF(WALL_W_AT(IJK)) THEN
449     
450                       A_M(IJK,E,M) = ZERO
451                       A_M(IJK,W,M) = ZERO
452                       A_M(IJK,N,M) = ZERO
453                       A_M(IJK,S,M) = ZERO
454                       A_M(IJK,T,M) = ZERO
455                       A_M(IJK,B,M) = ZERO
456                       A_M(IJK,0,M) = -ONE
457     
458                       B_M(IJK,M) = ZERO
459     
460                    ENDIF
461     
462                 CASE ('CG_FSW')
463     
464                    IF(WALL_W_AT(IJK)) THEN
465     
466                       A_M(IJK,E,M) = ZERO
467                       A_M(IJK,W,M) = ZERO
468                       A_M(IJK,N,M) = ZERO
469                       A_M(IJK,S,M) = ZERO
470                       A_M(IJK,T,M) = ZERO
471                       A_M(IJK,B,M) = ZERO
472                       A_M(IJK,0,M) = -ONE
473     
474     !                  B_M(IJK,M) = - W_g(W_MASTER_OF(IJK))  ! Velocity of master node
475     
476                       B_M(IJK,M) = ZERO
477     
478                       IF(DABS(NORMAL_W(IJK,3))/=ONE) THEN
479     
480                          IF (W_MASTER_OF(IJK) == EAST_OF(IJK)) THEN
481                             A_M(IJK,E,M) = ONE
482                          ELSEIF (W_MASTER_OF(IJK) == WEST_OF(IJK)) THEN
483                             A_M(IJK,W,M) = ONE
484                          ELSEIF (W_MASTER_OF(IJK) == NORTH_OF(IJK)) THEN
485                             A_M(IJK,N,M) = ONE
486                          ELSEIF (W_MASTER_OF(IJK) == SOUTH_OF(IJK)) THEN
487                             A_M(IJK,S,M) = ONE
488                          ELSEIF (W_MASTER_OF(IJK) == TOP_OF(IJK)) THEN
489                             A_M(IJK,T,M) = ONE
490                          ELSEIF (W_MASTER_OF(IJK) == BOTTOM_OF(IJK)) THEN
491                             A_M(IJK,B,M) = ONE
492                          ENDIF
493     
494                       ENDIF
495     
496                    ENDIF
497     
498                 CASE ('CG_PSW')
499     
500                    IF(WALL_W_AT(IJK)) THEN
501     
502                       A_M(IJK,E,M) = ZERO
503                       A_M(IJK,W,M) = ZERO
504                       A_M(IJK,N,M) = ZERO
505                       A_M(IJK,S,M) = ZERO
506                       A_M(IJK,T,M) = ZERO
507                       A_M(IJK,B,M) = ZERO
508                       A_M(IJK,0,M) = -ONE
509     
510     
511                       IF(BC_HW_G(BCV)==UNDEFINED) THEN   ! same as NSW
512                          B_M(IJK,M) = -BC_WW_G(BCV)
513                       ELSEIF(BC_HW_G(BCV)==ZERO) THEN   ! same as FSW
514                          B_M(IJK,M) = ZERO
515     
516                          IF(DABS(NORMAL_W(IJK,3))/=ONE) THEN
517     
518                             IF (W_MASTER_OF(IJK) == EAST_OF(IJK)) THEN
519                                A_M(IJK,E,M) = ONE
520                             ELSEIF (W_MASTER_OF(IJK) == WEST_OF(IJK)) THEN
521                                A_M(IJK,W,M) = ONE
522                             ELSEIF (W_MASTER_OF(IJK) == NORTH_OF(IJK)) THEN
523                                A_M(IJK,N,M) = ONE
524                             ELSEIF (W_MASTER_OF(IJK) == SOUTH_OF(IJK)) THEN
525                                A_M(IJK,S,M) = ONE
526                             ELSEIF (W_MASTER_OF(IJK) == TOP_OF(IJK)) THEN
527                                A_M(IJK,T,M) = ONE
528                             ELSEIF (W_MASTER_OF(IJK) == BOTTOM_OF(IJK)) THEN
529                                A_M(IJK,B,M) = ONE
530                             ENDIF
531     
532                          ENDIF
533     
534                       ELSE                              ! partial slip
535     
536     
537     
538     
539     
540                       ENDIF
541     
542                    ENDIF
543     
544     
545                 CASE ('CG_MI')
546     
547                    A_M(IJK,E,M) = ZERO
548                    A_M(IJK,W,M) = ZERO
549                    A_M(IJK,N,M) = ZERO
550                    A_M(IJK,S,M) = ZERO
551                    A_M(IJK,T,M) = ZERO
552                    A_M(IJK,B,M) = ZERO
553                    A_M(IJK,0,M) = -ONE
554     
555                    IF(BC_W_g(BCV)/=UNDEFINED) THEN
556                       B_M(IJK,M) = - BC_W_g(BCV)
557                    ELSE
558                       B_M(IJK,M) = - BC_VELMAG_g(BCV)*NORMAL_W(IJK,3)
559                    ENDIF
560     
561     
562                    IJKB = BOTTOM_OF(IJK)
563                    IF(FLUID_AT(IJKB)) THEN
564     
565                       A_M(IJKB,E,M) = ZERO
566                       A_M(IJKB,W,M) = ZERO
567                       A_M(IJKB,N,M) = ZERO
568                       A_M(IJKB,S,M) = ZERO
569                       A_M(IJKB,T,M) = ZERO
570                       A_M(IJKB,B,M) = ZERO
571                       A_M(IJKB,0,M) = -ONE
572     
573                       IF(BC_W_g(BCV)/=UNDEFINED) THEN
574                          B_M(IJKB,M) = - BC_W_g(BCV)
575                       ELSE
576                          B_M(IJKB,M) = - BC_VELMAG_g(BCV)*NORMAL_W(IJK,3)
577                       ENDIF
578     
579     
580                    ENDIF
581     
582                 CASE ('CG_PO')
583     
584                    A_M(IJK,E,M) = ZERO
585                    A_M(IJK,W,M) = ZERO
586                    A_M(IJK,N,M) = ZERO
587                    A_M(IJK,S,M) = ZERO
588                    A_M(IJK,T,M) = ZERO
589                    A_M(IJK,B,M) = ZERO
590                    A_M(IJK,0,M) = -ONE
591                    B_M(IJK,M) = ZERO
592     
593                    IJKB = BOTTOM_OF(IJK)
594                    IF(FLUID_AT(IJKB)) THEN
595     
596                       A_M(IJK,B,M) = ONE
597                       A_M(IJK,0,M) = -ONE
598     
599                    ENDIF
600     
601              END SELECT
602     
603              BCV = BC_ID(IJK)
604     
605              IF(BCV > 0 ) THEN
606                 BCT = BC_TYPE(BCV)
607              ELSE
608                 BCT = 'NONE'
609              ENDIF
610     
611              SELECT CASE (BCT)
612     
613                 CASE ('CG_MI')
614     
615                    A_M(IJK,E,M) = ZERO
616                    A_M(IJK,W,M) = ZERO
617                    A_M(IJK,N,M) = ZERO
618                    A_M(IJK,S,M) = ZERO
619                    A_M(IJK,T,M) = ZERO
620                    A_M(IJK,B,M) = ZERO
621                    A_M(IJK,0,M) = -ONE
622     
623                    IF(BC_W_g(BCV)/=UNDEFINED) THEN
624                       B_M(IJK,M) = - BC_W_g(BCV)
625                    ELSE
626                       B_M(IJK,M) = - BC_VELMAG_g(BCV)*NORMAL_S(IJK,3)
627                    ENDIF
628     
629     
630                    IJKB = BOTTOM_OF(IJK)
631                    IF(FLUID_AT(IJKB)) THEN
632     
633                       A_M(IJKB,E,M) = ZERO
634                       A_M(IJKB,W,M) = ZERO
635                       A_M(IJKB,N,M) = ZERO
636                       A_M(IJKB,S,M) = ZERO
637                       A_M(IJKB,T,M) = ZERO
638                       A_M(IJKB,B,M) = ZERO
639                       A_M(IJKB,0,M) = -ONE
640     
641                       IF(BC_W_g(BCV)/=UNDEFINED) THEN
642                          B_M(IJKB,M) = - BC_W_g(BCV)
643                       ELSE
644                          B_M(IJKB,M) = - BC_VELMAG_g(BCV)*NORMAL_S(IJK,3)
645                       ENDIF
646     
647     
648                    ENDIF
649     
650                 CASE ('CG_PO')
651     
652                    A_M(IJK,E,M) = ZERO
653                    A_M(IJK,W,M) = ZERO
654                    A_M(IJK,N,M) = ZERO
655                    A_M(IJK,S,M) = ZERO
656                    A_M(IJK,T,M) = ZERO
657                    A_M(IJK,B,M) = ZERO
658                    A_M(IJK,0,M) = -ONE
659                    B_M(IJK,M) = ZERO
660     
661                    IJKB = BOTTOM_OF(IJK)
662                    IF(FLUID_AT(IJKB)) THEN
663     
664                       A_M(IJK,B,M) = ONE
665                       A_M(IJK,0,M) = -ONE
666     
667                    ENDIF
668     
669              END SELECT
670     
671     
672           ENDDO
673     
674           RETURN
675     
676     !=======================================================================
677     ! JFD: END MODIFICATION FOR CARTESIAN GRID IMPLEMENTATION
678     !=======================================================================
679     
680           END SUBROUTINE CG_SOURCE_W_G_BC
681     
682     !// Comments on the modifications for DMP version implementation
683     !// 001 Include header file and common declarations for parallelization
684     !// 350 Changed do loop limits: 1,kmax2->kmin3,kmax3
685     !// 360 Check if i,j,k resides on current processor
686