File: /nfs/home/0/users/jenkins/mfix.git/model/cartesian_grid/CG_set_outflow.f

1     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
2     !                                                                      C
3     !  Module name: SET_OUTFLOW(BCV, I1, I2, J1, J2, K1, K2)               C
4     !  Purpose: Set specified pressure outflow bc for a specified range of C
5     !           cells                                                      C
6     !                                                                      C
7     !  Author: M. Syamlal                                 Date: 21-JAN-92  C
8     !  Reviewer:M. Syamlal, S. Venkatesan, P. Nicoletti,  Date: 29-JAN-92  C
9     !           W. Rogers                                                  C
10     !                                                                      C
11     !  Revision Number: 1                                                  C
12     !  Purpose: To incorporate Cartesian grid modifications                C
13     !  Author: Jeff Dietiker                              Date: 01-Jul-09  C
14     !                                                                      C
15     !  Literature/Document References:                                     C
16     !                                                                      C
17     !  Variables referenced: MMAX                                          C
18     !                                                                      C
19     !  Variables modified: I, J, K, RO_g, ROP_g,                           C
20     !                      EP_g, C
21     !                      T_g, T_s,  M, ROP_s, U_g, U_s, V_g, V_s,  C
22     !                      W_g, W_s,
23     !                                                                      C
24     !  Local variables: IJK, LFLUID                                      C
25     !                                                                      C
26     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
27     !
28           SUBROUTINE CG_SET_OUTFLOW
29     !...Translated by Pacific-Sierra Research VAST-90 2.06G5  12:17:31  12/09/98
30     !...Switches: -xf
31     !
32     !-----------------------------------------------
33     !   M o d u l e s
34     !-----------------------------------------------
35           USE bc
36           USE compar        !//d
37           USE cutcell
38           USE eos, ONLY: EOSG
39           USE fldvar
40           USE functions
41           USE geometry
42           USE indices
43           USE mflux
44           USE param
45           USE param1
46           USE physprop
47           USE quadric
48           USE run
49           USE scalars
50           IMPLICIT NONE
51     !-----------------------------------------------
52     !   G l o b a l   P a r a m e t e r s
53     !-----------------------------------------------
54     !-----------------------------------------------
55     !   D u m m y   A r g u m e n t s
56     !-----------------------------------------------
57     !
58     !                      indices
59           INTEGER          I, J, K, M, N
60     !
61     !                      Local index for boundary cell
62           INTEGER          IJK
63     !
64     !                      Boundary condition number
65           INTEGER          BCV
66     !
67     !                      Locall index for a fluid cell near the boundary cell
68           INTEGER          LFLUID
69     
70           INTEGER          IJKW,IJKWW,IJKS,IJKSS,IJKB
71     
72           CHARACTER(LEN=16) :: BCT1,BCT2,BCT3,BCT4
73     
74           LOGICAL :: TEST1,TEST2
75     !-----------------------------------------------
76     !
77     !      print*,'top of cg_set_outflow'
78           DO IJK = IJKSTART3, IJKEND3
79           IF(INTERIOR_CELL_AT(IJK)) THEN
80     
81     
82              I = I_OF(IJK)
83              J = J_OF(IJK)
84              K = K_OF(IJK)
85     
86     !//SP Check if current i,j,k resides on this PE
87              IF (.NOT.IS_ON_myPE_plus2layers(I,J,K)) CYCLE
88     !        IJK = FUNIJK(I,J,K)
89     !
90     ! Fluid cell at West
91     !
92     !      print*,'west'
93               BCV = BC_U_ID(IJK)
94     
95     
96               IF(BCV>0)THEN
97     
98                 IF (BC_TYPE(BCV) == 'CG_PO') THEN
99     
100                    IF (FLUID_AT(IM_OF(IJK))) THEN
101                       LFLUID = IM_OF(IJK)
102     !            print*,'west treatment:IJK,LFLUID=',IJK,LFLUID
103     !            read(*,*)
104     !                  IF (U_G(LFLUID)>=ZERO .OR. EP_G(IJK)==UNDEFINED) THEN
105     !                     IF (BC_TYPE(BCV) /= 'P_OUTFLOW') P_G(IJK) = P_G(LFLUID)
106                          T_G(IJK) = T_G(LFLUID)
107                          N = 1
108                          IF (NMAX(0) > 0) THEN
109                             X_G(IJK,:NMAX(0)) = X_G(LFLUID,:NMAX(0))
110                             N = NMAX(0) + 1
111                          ENDIF
112                          MW_MIX_G(IJK) = MW_MIX_G(LFLUID)
113                          IF (RO_G0 == UNDEFINED) RO_G(IJK) = EOSG(MW_MIX_G(IJK),P_G&
114                             (IJK),T_G(IJK))
115     !                  ENDIF
116                       P_STAR(IJK) = P_STAR(LFLUID)
117                       IF (BC_EP_G(BCV) == UNDEFINED) EP_G(IJK) = ONE
118     
119                       DO N = 1, NScalar
120                          M = Phase4Scalar(N)
121                          IF(M == 0)Then
122     !                      IF (U_G(LFLUID)>=ZERO) THEN
123                               Scalar(IJK, N) = Scalar(LFLUID, N)
124     !                      ENDIF
125                          Else
126     !                      IF (U_s(LFLUID, M)>=ZERO) THEN
127                               Scalar(IJK, N) = Scalar(LFLUID, N)
128     !                     ENDIF
129                          Endif
130                       END DO
131     
132                       IF(K_Epsilon) THEN
133     !                    IF (U_G(LFLUID) >= ZERO) THEN
134                          K_Turb_G(IJK) = K_Turb_G(LFLUID)
135                          E_Turb_G(IJK) = E_Turb_G(LFLUID)
136     !                   ENDIF
137                       ENDIF
138     
139     
140                       DO M = 1, MMAX
141                          P_S(IJK,M) = P_S(LFLUID,M)
142     !                     IF (U_S(LFLUID,M) >= ZERO) THEN
143                             ROP_S(IJK,M) = ROP_S(LFLUID,M)
144                             T_S(IJK,M) = T_S(LFLUID,M)
145     !                     ELSE
146     !                        ROP_S(IJK,M) = ZERO
147     !                     ENDIF
148     !
149                          IF(BC_ROP_S(BCV,M)/=UNDEFINED)ROP_S(IJK,M)=BC_ROP_S(BCV,M)
150     !
151                          IF(BC_EP_G(BCV)==UNDEFINED)EP_G(IJK)=EP_G(IJK)-EP_S(IJK,M)
152     !
153                          N = 1
154                          IF (NMAX(M) > 0) THEN
155                             X_S(IJK,M,:NMAX(M)) = X_S(LFLUID,M,:NMAX(M))
156                             N = NMAX(M) + 1
157                          ENDIF
158                       END DO
159                       ROP_G(IJK) = RO_G(IJK)*EP_G(IJK)
160     !                  IF (ROP_G(IJK) > ZERO) THEN
161     !                     U_G(IJK) = ROP_G(LFLUID)*U_G(LFLUID)/ROP_G(IJK)
162     !                  ELSE
163     !                     U_G(IJK) = ZERO
164     !                  ENDIF
165     !                  V_G(IJK) = V_G(LFLUID)
166     !                  W_G(IJK) = W_G(LFLUID)
167     !                  Flux_gE(IJK) = Flux_gE(LFLUID)
168     !                 Flux_gN(IJK) = Flux_gN(LFLUID)
169     !                 Flux_gT(IJK) = Flux_gT(LFLUID)
170     
171     !                  IF (MMAX > 0) THEN
172     !                     WHERE (ROP_S(IJK,:MMAX) > ZERO)
173     !                        U_S(IJK,:MMAX) = ROP_S(LFLUID,:MMAX)*U_S(LFLUID,:MMAX)/&
174     !                           ROP_S(IJK,:MMAX)
175     !                     ELSEWHERE
176     !                        U_S(IJK,:MMAX) = ZERO
177     !                    END WHERE
178     !                     V_S(IJK,:MMAX) = V_S(LFLUID,:MMAX)
179     !                     W_S(IJK,:MMAX) = W_S(LFLUID,:MMAX)
180     !                     Flux_sE(IJK,:MMAX) = Flux_sE(LFLUID,:MMAX)
181     !                    Flux_sN(IJK,:MMAX) = Flux_sN(LFLUID,:MMAX)
182     !                    Flux_sT(IJK,:MMAX) = Flux_sT(LFLUID,:MMAX)
183     !                  ENDIF
184                    ENDIF
185     
186     
187                 ENDIF
188              ENDIF
189     
190              IJKW = WEST_OF(IJK)
191              IJKWW = WEST_OF(IJKW)
192     
193              BCT1=''
194              BCT2=''
195              BCT3=''
196              BCV = BC_ID(IJK)
197              IF(BCV>0) BCT1 = BC_TYPE(BCV)
198              BCV = BC_U_ID(IJK)
199              IF(BCV>0) BCT2 = BC_TYPE(BCV)
200              BCV = BC_U_ID(IJKW)
201              IF(BCV>0) BCT3 = BC_TYPE(BCV)
202              BCV = BC_ID(IJKW)
203              IF(BCV>0) BCT4 = BC_TYPE(BCV)
204     
205              TEST1= ((BCT1 == 'CG_PO').AND.(BCT2 /= 'CG_PO').AND.(BCT3 == 'CG_PO'))
206     
207              TEST2 = (BCT2 == 'CG_PO').AND.(BCT4 /= 'CG_PO')
208     
209     !         TEST2 = (FLAG(IJK) == 11)
210     
211     
212              IF(TEST1) THEN
213                 BCV = BC_ID(IJK)
214     
215     !         IF(TEST1.OR.TEST2) THEN
216     
217     !            IF(FLUID_AT(IJKW)) THEN
218     
219     !            ELSEIF(FLUID_AT(IJKWW)) THEN
220     !               IJKW = IJKWW
221     !            ENDIF
222     
223     
224     !         IF((BCT1 == 'CG_PO').AND.(BCT2 /= 'CG_PO').AND.(BCT3 == 'CG_PO')) THEN
225     
226     
227     
228     !         IF(BCT2 == 'CG_PO')
229     
230     
231     !            print*,'IJK,IJKW=',IJK,IJKW
232     !            read(*,*)
233     
234                 T_G(IJK) = T_G(IJKW)
235                 N = 1
236                 IF (NMAX(0) > 0) THEN
237                    X_G(IJK,:NMAX(0)) = X_G(IJKW,:NMAX(0))
238                    N = NMAX(0) + 1
239                 ENDIF
240                 MW_MIX_G(IJK) = MW_MIX_G(IJKW)
241                 RO_G(IJK) = RO_G(IJKW)
242                 P_STAR(IJK) = P_STAR(IJKW)
243                 EP_G(IJK) = EP_G(IJKW)
244                 ROP_G(IJK) = ROP_G(IJKW)
245                 DO N = 1, NScalar
246                    M = Phase4Scalar(N)
247                    IF(M == 0)Then
248                       Scalar(IJK, N) = Scalar(IJKW, N)
249                    ELSE
250                       Scalar(IJK, N) = Scalar(IJKW, N)
251                    ENDIF
252                 ENDDO
253                 IF(K_Epsilon) THEN
254                    K_Turb_G(IJK) = K_Turb_G(IJKW)
255                    E_Turb_G(IJK) = E_Turb_G(IJKW)
256                 ENDIF
257     !            DO M = 1, MMAX
258     !               P_S(IJK,M) = P_S(IJKW,M)
259     !               ROP_S(IJK,M) = ROP_S(IJKW,M)
260     !               T_S(IJK,M) = T_S(IJKW,M)
261     !               N = 1
262     !               IF (NMAX(M) > 0) THEN
263     !                 X_S(IJK,M,:NMAX(M)) = X_S(IJKW,M,:NMAX(M))
264     !                 N = NMAX(M) + 1
265     !               ENDIF
266     !            ENDDO
267     
268                 DO M = 1, MMAX
269                    P_S(IJK,M) = P_S(IJKW,M)
270     !               IF (U_S(IJKW,M) >= ZERO) THEN
271                       ROP_S(IJK,M) = ROP_S(IJKW,M)
272                       T_S(IJK,M) = T_S(IJKW,M)
273     !               ELSE
274     !                  ROP_S(IJK,M) = ZERO
275     !               ENDIF
276     !
277                    IF(BC_ROP_S(BCV,M)/=UNDEFINED)ROP_S(IJK,M)=BC_ROP_S(BCV,M)
278     !
279                    IF(BC_EP_G(BCV)==UNDEFINED)EP_G(IJK)=EP_G(IJK)-EP_S(IJK,M)
280     !
281                    N = 1
282                    IF (NMAX(M) > 0) THEN
283                       X_S(IJK,M,:NMAX(M)) = X_S(IJKW,M,:NMAX(M))
284                       N = NMAX(M) + 1
285                    ENDIF
286                 END DO
287                 ROP_G(IJK) = RO_G(IJK)*EP_G(IJK)
288     !            IF (ROP_G(IJK) > ZERO) THEN
289     !               U_G(IJK) = ROP_G(IJKW)*U_G(IJKW)/ROP_G(IJK)
290     !            ELSE
291     !               U_G(IJK) = ZERO
292     !            ENDIF
293     !            V_G(IJK) = V_G(IJKW)
294     !            W_G(IJK) = W_G(IJKW)
295     !            Flux_gE(IJK) = Flux_gE(IJKW)
296     !            Flux_gN(IJK) = Flux_gN(IJKW)
297     !            Flux_gT(IJK) = Flux_gT(IJKW)
298     
299     !            IF (MMAX > 0) THEN
300     !               WHERE (ROP_S(IJK,:MMAX) > ZERO)
301     !                  U_S(IJK,:MMAX) = ROP_S(IJKW,:MMAX)*U_S(IJKW,:MMAX)/&
302     !                     ROP_S(IJK,:MMAX)
303     !               ELSEWHERE
304     !                  U_S(IJK,:MMAX) = ZERO
305     !               END WHERE
306     !               V_S(IJK,:MMAX) = V_S(IJKW,:MMAX)
307     !               W_S(IJK,:MMAX) = W_S(IJKW,:MMAX)
308     !               Flux_sE(IJK,:MMAX) = Flux_sE(IJKW,:MMAX)
309     !              Flux_sN(IJK,:MMAX) = Flux_sN(IJKW,:MMAX)
310     !              Flux_sT(IJK,:MMAX) = Flux_sT(IJKW,:MMAX)
311     !            ENDIF
312     
313     
314     
315     !            CYCLE
316     
317              ENDIF
318     
319     
320     
321     
322     !
323     ! Fluid cell at East
324     !      print*,'east'
325     
326               BCV = BC_U_ID(IJK)
327     
328               IF(BCV>0)THEN
329     
330                 IF (BC_TYPE(BCV) == 'CG_PO') THEN
331     
332     
333                    IF (FLUID_AT(IP_OF(IJK))) THEN
334                       LFLUID = IP_OF(IJK)
335     !            print*,'east treatment:IJK,LFLUID=',IJK,LFLUID
336     !            read(*,*)
337     !                  IF (U_G(IJK)<=ZERO .OR. EP_G(IJK)==UNDEFINED) THEN
338     !                     IF (BC_TYPE(BCV) /= 'P_OUTFLOW') P_G(IJK) = P_G(LFLUID)
339                          T_G(IJK) = T_G(LFLUID)
340                          N = 1
341                          IF (NMAX(0) > 0) THEN
342                             X_G(IJK,:NMAX(0)) = X_G(LFLUID,:NMAX(0))
343                             N = NMAX(0) + 1
344                          ENDIF
345                          MW_MIX_G(IJK) = MW_MIX_G(LFLUID)
346                          IF (RO_G0 == UNDEFINED) RO_G(IJK) = EOSG(MW_MIX_G(IJK),P_G&
347                             (IJK),T_G(IJK))
348     !                  ENDIF
349                       P_STAR(IJK) = P_STAR(LFLUID)
350                       IF (BC_EP_G(BCV) == UNDEFINED) EP_G(IJK) = ONE
351     
352                       DO N = 1, NScalar
353                          M = Phase4Scalar(N)
354                          IF(M == 0)Then
355     !                      IF (U_G(LFLUID) <= ZERO) THEN
356                               Scalar(IJK, N) = Scalar(LFLUID, N)
357     !                      ENDIF
358                          Else
359     !                      IF (U_s(LFLUID, M) <= ZERO) THEN
360                               Scalar(IJK, N) = Scalar(LFLUID, N)
361     !                      ENDIF
362                          Endif
363                       END DO
364     
365                       IF(K_Epsilon) THEN
366     !                    IF (U_G(LFLUID) <= ZERO) THEN
367                             K_Turb_G(IJK) = K_Turb_G(LFLUID)
368                             E_Turb_G(IJK) = E_Turb_G(LFLUID)
369     !                    ENDIF
370                       ENDIF
371     
372                       DO M = 1, MMAX
373                          P_S(IJK,M) = P_S(LFLUID,M)
374     !                     IF (U_S(IJK,M) <= ZERO) THEN
375                             ROP_S(IJK,M) = ROP_S(LFLUID,M)
376                             T_S(IJK,M) = T_S(LFLUID,M)
377     !                     ELSE
378     !                        ROP_S(IJK,M) = ZERO
379     !                     ENDIF
380     !
381                          IF(BC_ROP_S(BCV,M)/=UNDEFINED)ROP_S(IJK,M)=BC_ROP_S(BCV,M)
382     !
383                          IF(BC_EP_G(BCV)==UNDEFINED)EP_G(IJK)=EP_G(IJK)-EP_S(IJK,M)
384                          N = 1
385                          IF (NMAX(M) > 0) THEN
386                             X_S(IJK,M,:NMAX(M)) = X_S(LFLUID,M,:NMAX(M))
387                             N = NMAX(M) + 1
388                          ENDIF
389                       END DO
390                       ROP_G(IJK) = RO_G(IJK)*EP_G(IJK)
391     !                  IF (U_G(IJK) == UNDEFINED) THEN
392     !                     IF (ROP_G(IJK) > ZERO) THEN
393     !                        U_G(IJK) = ROP_G(LFLUID)*U_G(LFLUID)/ROP_G(IJK)
394     !                     ELSE
395     !                        U_G(IJK) = ZERO
396     !                     ENDIF
397     !                  ENDIF
398     !                  V_G(IJK) = V_G(LFLUID)
399     !                  W_G(IJK) = W_G(LFLUID)
400     !                  Flux_gE(IJK) = Flux_gE(LFLUID)
401     !                  Flux_gN(IJK) = Flux_gN(LFLUID)
402     !                  Flux_gT(IJK) = Flux_gT(LFLUID)
403     
404     !                  DO M = 1, MMAX
405     !                     IF (U_S(IJK,M) == UNDEFINED) THEN
406     !                        IF (ROP_S(IJK,M) > ZERO) THEN
407     !                           U_S(IJK,M) = ROP_S(LFLUID,M)*U_S(LFLUID,M)/ROP_S(IJK&
408     !                              ,M)
409     !                        ELSE
410     !                           U_S(IJK,M) = ZERO
411     !                        ENDIF
412     !                     ENDIF
413     !                     V_S(IJK,M) = V_S(LFLUID,M)
414     !                     W_S(IJK,M) = W_S(LFLUID,M)
415     !                     Flux_sE(IJK,M) = Flux_sE(LFLUID,M)
416     !                    Flux_sN(IJK,M) = Flux_sN(LFLUID,M)
417     !                    Flux_sT(IJK,M) = Flux_sT(LFLUID,M)
418     !                  END DO
419                    ENDIF
420     
421                 ENDIF
422              ENDIF
423     
424     !
425     ! Fluid cell at South
426     !
427     !      print*,'south'
428               BCV = BC_V_ID(IJK)
429     !      print*,ijk,I,J,K,bcv,BC_TYPE(BCV)
430               IF(BCV>0)THEN
431     
432                 IF (BC_TYPE(BCV) == 'CG_PO') THEN
433     
434     
435                    IF (FLUID_AT(JM_OF(IJK))) THEN
436                       LFLUID = JM_OF(IJK)
437     !            print*,'south treatment:IJK,LFLUID=',IJK,LFLUID
438     !            read(*,*)
439     
440     
441     !                  IF (V_G(LFLUID)>=ZERO .OR. EP_G(IJK)==UNDEFINED) THEN
442     !                     IF (BC_TYPE(BCV) /= 'P_OUTFLOW') P_G(IJK) = P_G(LFLUID)
443                          T_G(IJK) = T_G(LFLUID)
444                          N = 1
445                          IF (NMAX(0) > 0) THEN
446                             X_G(IJK,:NMAX(0)) = X_G(LFLUID,:NMAX(0))
447                             N = NMAX(0) + 1
448                          ENDIF
449                          MW_MIX_G(IJK) = MW_MIX_G(LFLUID)
450                          IF (RO_G0 == UNDEFINED) RO_G(IJK) = EOSG(MW_MIX_G(IJK),P_G&
451                             (IJK),T_G(IJK))
452     !                  ENDIF
453                       P_STAR(IJK) = P_STAR(LFLUID)
454                       IF (BC_EP_G(BCV) == UNDEFINED) EP_G(IJK) = ONE
455     
456                       DO N = 1, NScalar
457                          M = Phase4Scalar(N)
458                          IF(M == 0)Then
459     !                      IF (V_G(LFLUID) >= ZERO) THEN
460                               Scalar(IJK, N) = Scalar(LFLUID, N)
461     !                      ENDIF
462                          Else
463     !                      IF (V_s(LFLUID, M) >= ZERO) THEN
464                            Scalar(IJK, N) = Scalar(LFLUID, N)
465     !                     ENDIF
466                          Endif
467                       END DO
468     
469                       IF(K_Epsilon) THEN
470     !                    IF (V_G(LFLUID) >= ZERO) THEN
471                             K_Turb_G(IJK) = K_Turb_G(LFLUID)
472                             E_Turb_G(IJK) = E_Turb_G(LFLUID)
473     !                   ENDIF
474                       ENDIF
475     
476                       DO M = 1, MMAX
477                          P_S(IJK,M) = P_S(LFLUID,M)
478     !                     IF (V_S(LFLUID,M) >= 0.) THEN
479                             ROP_S(IJK,M) = ROP_S(LFLUID,M)
480                             T_S(IJK,M) = T_S(LFLUID,M)
481     !                     ELSE
482     !                        ROP_S(IJK,M) = ZERO
483     !                     ENDIF
484     !
485                          IF(BC_ROP_S(BCV,M)/=UNDEFINED)ROP_S(IJK,M)=BC_ROP_S(BCV,M)
486     !
487                          IF(BC_EP_G(BCV)==UNDEFINED)EP_G(IJK)=EP_G(IJK)-EP_S(IJK,M)
488                          N = 1
489                          IF (NMAX(M) > 0) THEN
490                             X_S(IJK,M,:NMAX(M)) = X_S(LFLUID,M,:NMAX(M))
491                             N = NMAX(M) + 1
492                          ENDIF
493                       END DO
494                       ROP_G(IJK) = RO_G(IJK)*EP_G(IJK)
495     !                  U_G(IJK) = U_G(LFLUID)
496     !                  IF (ROP_G(IJK) > ZERO) THEN
497     !                     V_G(IJK) = ROP_G(LFLUID)*V_G(LFLUID)/ROP_G(IJK)
498     !                  ELSE
499     !                     V_G(IJK) = ZERO
500     !                  ENDIF
501     !                  W_G(IJK) = W_G(LFLUID)
502     !                  Flux_gE(IJK) = Flux_gE(LFLUID)
503     !                  Flux_gN(IJK) = Flux_gN(LFLUID)
504     !                  Flux_gT(IJK) = Flux_gT(LFLUID)
505     
506     !                  IF (MMAX > 0) THEN
507     !                     U_S(IJK,:MMAX) = U_S(LFLUID,:MMAX)
508     !                     V_S(IJK,:MMAX) = V_S(LFLUID,:MMAX)
509     !                     W_S(IJK,:MMAX) = W_S(LFLUID,:MMAX)
510     !                     Flux_sE(IJK,:MMAX) = Flux_sE(LFLUID,:MMAX)
511     !                    Flux_sN(IJK,:MMAX) = Flux_sN(LFLUID,:MMAX)
512     !                    Flux_sT(IJK,:MMAX) = Flux_sT(LFLUID,:MMAX)
513     !                  ENDIF
514                    ENDIF
515     
516                 ENDIF
517              ENDIF
518     
519              IJKS = SOUTH_OF(IJK)
520              IJKSS = SOUTH_OF(IJKS)
521     !         print*,'ijk,ijks=',ijk,ijks,JM_OF(IJK)
522              BCT1=''
523              BCT2=''
524              BCT3=''
525              BCV = BC_ID(IJK)
526     !         print*,'bcv=',bcv
527              IF(BCV>0) BCT1 = BC_TYPE(BCV)
528              BCV = BC_V_ID(IJK)
529     !         print*,'bcv=',bcv
530              IF(BCV>0) BCT2 = BC_TYPE(BCV)
531              BCV = BC_V_ID(IJKS)
532     !         print*,'bcv=',bcv
533              IF(BCV>0) BCT3 = BC_TYPE(BCV)
534     
535     !         IF((BCT1 == 'CG_PO').AND.(BCT2 /= 'CG_PO').AND.(BCT3 == 'CG_PO')) THEN
536     
537     
538     
539              TEST1 = (BCT1 == 'CG_PO').AND.(BCT2 /= 'CG_PO').AND.(BCT3 == 'CG_PO')
540              TEST2 = (FLAG(IJK) == 11).AND.(BCT2 == 'CG_PO')
541     
542              TEST2 = (FLAG(IJK) == 11)
543     
544              IF(TEST1) THEN
545     
546                 BCV = BC_ID(IJK)
547     
548     !         IF(TEST1.OR.TEST2) THEN
549     
550     !            IF(FLUID_AT(IJKS)) THEN
551     
552     !            ELSEIF(FLUID_AT(IJKSS)) THEN
553     !               IJKS = IJKSS
554     !            ENDIF
555     
556     
557     !            print*,'IJK,IJKS=',IJK,IJKS
558     !            read(*,*)
559     
560                 T_G(IJK) = T_G(IJKS)
561                 N = 1
562                 IF (NMAX(0) > 0) THEN
563                    X_G(IJK,:NMAX(0)) = X_G(IJKS,:NMAX(0))
564                    N = NMAX(0) + 1
565                 ENDIF
566                 MW_MIX_G(IJK) = MW_MIX_G(IJKS)
567                 RO_G(IJK) = RO_G(IJKS)
568                 P_STAR(IJK) = P_STAR(IJKS)
569                 EP_G(IJK) = EP_G(IJKS)
570                 ROP_G(IJK) = ROP_G(IJKS)
571                 DO N = 1, NScalar
572                    M = Phase4Scalar(N)
573                    IF(M == 0)Then
574                       Scalar(IJK, N) = Scalar(IJKS, N)
575                    ELSE
576                       Scalar(IJK, N) = Scalar(IJKS, N)
577                    ENDIF
578                 ENDDO
579                 IF(K_Epsilon) THEN
580                    K_Turb_G(IJK) = K_Turb_G(IJKS)
581                    E_Turb_G(IJK) = E_Turb_G(IJKS)
582                 ENDIF
583     !            DO M = 1, MMAX
584     !               P_S(IJK,M) = P_S(IJKS,M)
585     !               ROP_S(IJK,M) = ROP_S(IJKS,M)
586     !               T_S(IJK,M) = T_S(IJKS,M)
587     !               N = 1
588     !               IF (NMAX(M) > 0) THEN
589     !                 X_S(IJK,M,:NMAX(M)) = X_S(IJKS,M,:NMAX(M))
590     !                 N = NMAX(M) + 1
591     !               ENDIF
592     !            ENDDO
593     
594                 DO M = 1, MMAX
595                    P_S(IJK,M) = P_S(IJKS,M)
596     !               IF (V_S(IJKS,M) >= 0.) THEN
597                       ROP_S(IJK,M) = ROP_S(IJKS,M)
598                       T_S(IJK,M) = T_S(IJKS,M)
599     !               ELSE
600     !                  ROP_S(IJK,M) = ZERO
601     !               ENDIF
602     !
603                    IF(BC_ROP_S(BCV,M)/=UNDEFINED)ROP_S(IJK,M)=BC_ROP_S(BCV,M)
604     !
605                    IF(BC_EP_G(BCV)==UNDEFINED)EP_G(IJK)=EP_G(IJK)-EP_S(IJK,M)
606                    N = 1
607                    IF (NMAX(M) > 0) THEN
608                       X_S(IJK,M,:NMAX(M)) = X_S(IJKS,M,:NMAX(M))
609                       N = NMAX(M) + 1
610                    ENDIF
611                 END DO
612                 ROP_G(IJK) = RO_G(IJK)*EP_G(IJK)
613     !            U_G(IJK) = U_G(IJKS)
614     !            IF (ROP_G(IJK) > ZERO) THEN
615     !               V_G(IJK) = ROP_G(IJKS)*V_G(IJKS)/ROP_G(IJK)
616     !            ELSE
617     !               V_G(IJK) = ZERO
618     !            ENDIF
619     !            W_G(IJK) = W_G(IJKS)
620     !            Flux_gE(IJK) = Flux_gE(IJKS)
621     !           Flux_gN(IJK) = Flux_gN(IJKS)
622     !           Flux_gT(IJK) = Flux_gT(IJKS)
623     
624     !            IF (MMAX > 0) THEN
625     !               U_S(IJK,:MMAX) = U_S(IJKS,:MMAX)
626     !               V_S(IJK,:MMAX) = V_S(IJKS,:MMAX)
627     !               W_S(IJK,:MMAX) = W_S(IJKS,:MMAX)
628     !               Flux_sE(IJK,:MMAX) = Flux_sE(IJKS,:MMAX)
629     !              Flux_sN(IJK,:MMAX) = Flux_sN(IJKS,:MMAX)
630     !              Flux_sT(IJK,:MMAX) = Flux_sT(IJKS,:MMAX)
631     !            ENDIF
632     
633     
634     
635     
636     !            CYCLE
637     
638              ENDIF
639     
640     
641     
642     !
643     ! Fluid cell at North
644     !
645     
646     !      print*,'north'
647               BCV = BC_V_ID(IJK)
648     
649               IF(BCV>0)THEN
650     
651                 IF (BC_TYPE(BCV) == 'CG_PO') THEN
652     
653     
654                    IF (FLUID_AT(JP_OF(IJK))) THEN
655                       LFLUID = JP_OF(IJK)
656     
657     !            print*,'north treatment:IJK,LFLUID=',IJK,LFLUID
658     !            read(*,*)
659     
660     !                  IF (V_G(IJK)<=ZERO .OR. EP_G(IJK)==UNDEFINED) THEN
661     !                     IF (BC_TYPE(BCV) /= 'P_OUTFLOW') P_G(IJK) = P_G(LFLUID)
662                          T_G(IJK) = T_G(LFLUID)
663                          N = 1
664                          IF (NMAX(0) > 0) THEN
665                             X_G(IJK,:NMAX(0)) = X_G(LFLUID,:NMAX(0))
666                             N = NMAX(0) + 1
667                          ENDIF
668                          MW_MIX_G(IJK) = MW_MIX_G(LFLUID)
669                          IF (RO_G0 == UNDEFINED) RO_G(IJK) = EOSG(MW_MIX_G(IJK),P_G&
670                             (IJK),T_G(IJK))
671     !                  ENDIF
672                       P_STAR(IJK) = P_STAR(LFLUID)
673                       IF (BC_EP_G(BCV) == UNDEFINED) EP_G(IJK) = ONE
674     
675                       DO N = 1, NScalar
676                          M = Phase4Scalar(N)
677                          IF(M == 0)Then
678     !                      IF (V_G(LFLUID) <= ZERO) THEN
679                              Scalar(IJK, N) = Scalar(LFLUID, N)
680     !                     ENDIF
681                          Else
682     !                      IF (V_s(LFLUID, M) <= ZERO) THEN
683                              Scalar(IJK, N) = Scalar(LFLUID, N)
684     !                     ENDIF
685                          Endif
686                       END DO
687     
688                       IF(K_Epsilon) THEN
689     !                    IF (V_G(LFLUID) <= ZERO) THEN
690                            K_Turb_G(IJK) = K_Turb_G(LFLUID)
691                            E_Turb_G(IJK) = E_Turb_G(LFLUID)
692     !                   ENDIF
693                       ENDIF
694     
695                       DO M = 1, MMAX
696                          P_S(IJK,M) = P_S(LFLUID,M)
697     !                     IF (V_S(IJK,M) <= ZERO) THEN
698                             ROP_S(IJK,M) = ROP_S(LFLUID,M)
699                             T_S(IJK,M) = T_S(LFLUID,M)
700     !                     ELSE
701     !                        ROP_S(IJK,M) = ZERO
702     !                     ENDIF
703     !
704                          IF(BC_ROP_S(BCV,M)/=UNDEFINED)ROP_S(IJK,M)=BC_ROP_S(BCV,M)
705     !
706                          IF(BC_EP_G(BCV)==UNDEFINED)EP_G(IJK)=EP_G(IJK)-EP_S(IJK,M)
707                          N = 1
708                          IF (NMAX(M) > 0) THEN
709                             X_S(IJK,M,:NMAX(M)) = X_S(LFLUID,M,:NMAX(M))
710                             N = NMAX(M) + 1
711                          ENDIF
712                       END DO
713                       ROP_G(IJK) = RO_G(IJK)*EP_G(IJK)
714     !                  U_G(IJK) = U_G(LFLUID)
715     !                  IF (V_G(IJK) == UNDEFINED) THEN
716     !                     IF (ROP_G(IJK) > ZERO) THEN
717     !                        V_G(IJK) = ROP_G(LFLUID)*V_G(LFLUID)/ROP_G(IJK)
718     !                     ELSE
719     !                        V_G(IJK) = ZERO
720     !                     ENDIF
721     !                  ENDIF
722     !                  W_G(IJK) = W_G(LFLUID)
723     !                  Flux_gE(IJK) = Flux_gE(LFLUID)
724     !                 Flux_gN(IJK) = Flux_gN(LFLUID)
725     !                 Flux_gT(IJK) = Flux_gT(LFLUID)
726     
727     !                  IF (MMAX > 0) THEN
728     !                     U_S(IJK,:MMAX) = U_S(LFLUID,:MMAX)
729     !                     WHERE (V_S(IJK,:MMAX) == UNDEFINED) V_S(IJK,:MMAX) = V_S(&
730     !                        LFLUID,:MMAX)
731     !                     W_S(IJK,:MMAX) = W_S(LFLUID,:MMAX)
732     !                     Flux_sE(IJK,:MMAX) = Flux_sE(LFLUID,:MMAX)
733     !                    Flux_sN(IJK,:MMAX) = Flux_sN(LFLUID,:MMAX)
734     !                    Flux_sT(IJK,:MMAX) = Flux_sT(LFLUID,:MMAX)
735     !                  ENDIF
736                    ENDIF
737     
738                 ENDIF
739              ENDIF
740     
741     !
742     ! Fluid cell at Bottom
743     !
744     !      print*,'bottom'
745     
746               BCV = BC_W_ID(IJK)
747     
748               IF(BCV>0)THEN
749     
750                 IF (BC_TYPE(BCV) == 'CG_PO') THEN
751     
752                    IF (FLUID_AT(KM_OF(IJK))) THEN
753                       LFLUID = KM_OF(IJK)
754     !                  IF (W_G(LFLUID)>=ZERO .OR. EP_G(IJK)==UNDEFINED) THEN
755     !                     IF (BC_TYPE(BCV) /= 'P_OUTFLOW') P_G(IJK) = P_G(LFLUID)
756                          T_G(IJK) = T_G(LFLUID)
757                          N = 1
758                          IF (NMAX(0) > 0) THEN
759                             X_G(IJK,:NMAX(0)) = X_G(LFLUID,:NMAX(0))
760                             N = NMAX(0) + 1
761                          ENDIF
762                          MW_MIX_G(IJK) = MW_MIX_G(LFLUID)
763                          IF (RO_G0 == UNDEFINED) RO_G(IJK) = EOSG(MW_MIX_G(IJK),P_G&
764                             (IJK),T_G(IJK))
765     !                  ENDIF
766                       P_STAR(IJK) = P_STAR(LFLUID)
767                       IF (BC_EP_G(BCV) == UNDEFINED) EP_G(IJK) = ONE
768     
769                       DO N = 1, NScalar
770                          M = Phase4Scalar(N)
771                          IF(M == 0)Then
772     !                      IF (W_G(LFLUID) >= ZERO) THEN
773                              Scalar(IJK, N) = Scalar(LFLUID, N)
774     !                      ENDIF
775                          Else
776     !                      IF (W_s(LFLUID, M) >= ZERO) THEN
777                              Scalar(IJK, N) = Scalar(LFLUID, N)
778     !                     ENDIF
779                          Endif
780                       END DO
781     
782                       IF(K_Epsilon) THEN
783     !                    IF (W_G(LFLUID) >= ZERO) THEN
784                            K_Turb_G(IJK) = K_Turb_G(LFLUID)
785                            E_Turb_G(IJK) = E_Turb_G(LFLUID)
786     !                   ENDIF
787                       ENDIF
788     
789                       DO M = 1, MMAX
790                          P_S(IJK,M) = P_S(LFLUID,M)
791     !                     IF (W_S(LFLUID,M) >= 0.) THEN
792                             ROP_S(IJK,M) = ROP_S(LFLUID,M)
793                             T_S(IJK,M) = T_S(LFLUID,M)
794     !                     ELSE
795                             ROP_S(IJK,M) = ZERO
796     !                     ENDIF
797     !
798                          IF(BC_ROP_S(BCV,M)/=UNDEFINED)ROP_S(IJK,M)=BC_ROP_S(BCV,M)
799     !
800                          IF(BC_EP_G(BCV)==UNDEFINED)EP_G(IJK)=EP_G(IJK)-EP_S(IJK,M)
801                          N = 1
802                          IF (NMAX(M) > 0) THEN
803                             X_S(IJK,M,:NMAX(M)) = X_S(LFLUID,M,:NMAX(M))
804                             N = NMAX(M) + 1
805                          ENDIF
806                       END DO
807                       ROP_G(IJK) = RO_G(IJK)*EP_G(IJK)
808     !                  U_G(IJK) = U_G(LFLUID)
809     !                  V_G(IJK) = V_G(LFLUID)
810     !                  IF (ROP_G(IJK) > ZERO) THEN
811     !                     W_G(IJK) = ROP_G(LFLUID)*W_G(LFLUID)/ROP_G(IJK)
812     !                  ELSE
813     !                     W_G(IJK) = ZERO
814     !                  ENDIF
815     !                  Flux_gE(IJK) = Flux_gE(LFLUID)
816     !                 Flux_gN(IJK) = Flux_gN(LFLUID)
817     !                 Flux_gT(IJK) = Flux_gT(LFLUID)
818     
819     !                  IF (MMAX > 0) THEN
820     !                     U_S(IJK,:MMAX) = U_S(LFLUID,:MMAX)
821     !                     V_S(IJK,:MMAX) = V_S(LFLUID,:MMAX)
822     !                     W_S(IJK,:MMAX) = W_S(LFLUID,:MMAX)
823     !                     Flux_sE(IJK,:MMAX) = Flux_sE(LFLUID,:MMAX)
824     !                    Flux_sN(IJK,:MMAX) = Flux_sN(LFLUID,:MMAX)
825     !                    Flux_sT(IJK,:MMAX) = Flux_sT(LFLUID,:MMAX)
826     !                  ENDIF
827                    ENDIF
828     
829     
830                 ENDIF
831              ENDIF
832     
833              IJKB = BOTTOM_OF(IJK)
834              BCT1=''
835              BCT2=''
836              BCT3=''
837              BCV = BC_ID(IJK)
838              IF(BCV>0) BCT1 = BC_TYPE(BCV)
839              BCV = BC_W_ID(IJK)
840              IF(BCV>0) BCT2 = BC_TYPE(BCV)
841              BCV = BC_W_ID(IJKB)
842              IF(BCV>0) BCT3 = BC_TYPE(BCV)
843     
844              IF((BCT1 == 'CG_PO').AND.(BCT2 /= 'CG_PO').AND.(BCT3 == 'CG_PO')) THEN
845     !            print*,'IJK,IJKB=',IJK,IJKB
846     !            read(*,*)
847     
848                 BCV = BC_ID(IJK)
849     
850                 T_G(IJK) = T_G(IJKB)
851                 N = 1
852                 IF (NMAX(0) > 0) THEN
853                    X_G(IJK,:NMAX(0)) = X_G(IJKB,:NMAX(0))
854                    N = NMAX(0) + 1
855                 ENDIF
856                 MW_MIX_G(IJK) = MW_MIX_G(IJKB)
857                 RO_G(IJK) = RO_G(IJKB)
858                 P_STAR(IJK) = P_STAR(IJKB)
859                 EP_G(IJK) = EP_G(IJKB)
860                 ROP_G(IJK) = ROP_G(IJKB)
861                 DO N = 1, NScalar
862                    M = Phase4Scalar(N)
863                    IF(M == 0)Then
864                       Scalar(IJK, N) = Scalar(IJKB, N)
865                    ELSE
866                       Scalar(IJK, N) = Scalar(IJKB, N)
867                    ENDIF
868                 ENDDO
869                 IF(K_Epsilon) THEN
870                    K_Turb_G(IJK) = K_Turb_G(IJKB)
871                    E_Turb_G(IJK) = E_Turb_G(IJKB)
872                 ENDIF
873     !            DO M = 1, MMAX
874     !               P_S(IJK,M) = P_S(IJKB,M)
875     !               ROP_S(IJK,M) = ROP_S(IJKB,M)
876     !               T_S(IJK,M) = T_S(IJKB,M)
877     !               N = 1
878     !               IF (NMAX(M) > 0) THEN
879     !                 X_S(IJK,M,:NMAX(M)) = X_S(IJKB,M,:NMAX(M))
880     !                 N = NMAX(M) + 1
881     !               ENDIF
882     !            ENDDO
883     
884                 DO M = 1, MMAX
885                    P_S(IJK,M) = P_S(IJKB,M)
886     !               IF (W_S(IJKB,M) >= 0.) THEN
887                       ROP_S(IJK,M) = ROP_S(IJKB,M)
888                       T_S(IJK,M) = T_S(IJKB,M)
889     !               ELSE
890     !                  ROP_S(IJK,M) = ZERO
891     !               ENDIF
892     !
893                    IF(BC_ROP_S(BCV,M)/=UNDEFINED)ROP_S(IJK,M)=BC_ROP_S(BCV,M)
894     !
895                    IF(BC_EP_G(BCV)==UNDEFINED)EP_G(IJK)=EP_G(IJK)-EP_S(IJK,M)
896                    N = 1
897                    IF (NMAX(M) > 0) THEN
898                       X_S(IJK,M,:NMAX(M)) = X_S(IJKB,M,:NMAX(M))
899                       N = NMAX(M) + 1
900                    ENDIF
901                 END DO
902                 ROP_G(IJK) = RO_G(IJK)*EP_G(IJK)
903     !            U_G(IJK) = U_G(IJKB)
904     !            V_G(IJK) = V_G(IJKB)
905     !            IF (ROP_G(IJK) > ZERO) THEN
906     !               W_G(IJK) = ROP_G(IJKB)*W_G(IJKB)/ROP_G(IJK)
907     !            ELSE
908     !               W_G(IJK) = ZERO
909     !            ENDIF
910     !            Flux_gE(IJK) = Flux_gE(IJKB)
911     !           Flux_gN(IJK) = Flux_gN(IJKB)
912     !           Flux_gT(IJK) = Flux_gT(IJKB)
913     
914     !            IF (MMAX > 0) THEN
915     !               U_S(IJK,:MMAX) = U_S(IJKB,:MMAX)
916     !               V_S(IJK,:MMAX) = V_S(IJKB,:MMAX)
917     !               W_S(IJK,:MMAX) = W_S(IJKB,:MMAX)
918     !               Flux_sE(IJK,:MMAX) = Flux_sE(IJKB,:MMAX)
919     !              Flux_sN(IJK,:MMAX) = Flux_sN(IJKB,:MMAX)
920     !              Flux_sT(IJK,:MMAX) = Flux_sT(IJKB,:MMAX)
921     !            ENDIF
922     
923     !            CYCLE
924     
925              ENDIF
926     
927     
928     !
929     ! Fluid cell at Top
930     !
931     !      print*,'top'
932               BCV = BC_W_ID(IJK)
933     
934               IF(BCV>0)THEN
935     
936                 IF (BC_TYPE(BCV) == 'CG_PO') THEN
937     
938     
939                    IF (FLUID_AT(KP_OF(IJK))) THEN
940                       LFLUID = KP_OF(IJK)
941     !                  IF (W_G(IJK)<=ZERO .OR. EP_G(IJK)==UNDEFINED) THEN
942     !                     IF (BC_TYPE(BCV) /= 'P_OUTFLOW') P_G(IJK) = P_G(LFLUID)
943                          T_G(IJK) = T_G(LFLUID)
944                          N = 1
945                          IF (NMAX(0) > 0) THEN
946                             X_G(IJK,:NMAX(0)) = X_G(LFLUID,:NMAX(0))
947                             N = NMAX(0) + 1
948                          ENDIF
949                          MW_MIX_G(IJK) = MW_MIX_G(LFLUID)
950                          IF (RO_G0 == UNDEFINED) RO_G(IJK) = EOSG(MW_MIX_G(IJK),P_G&
951                             (IJK),T_G(IJK))
952     !                  ENDIF
953                       P_STAR(IJK) = P_STAR(LFLUID)
954                       IF (BC_EP_G(BCV) == UNDEFINED) EP_G(IJK) = ONE
955     
956                       DO N = 1, NScalar
957                          M = Phase4Scalar(N)
958                          IF(M == 0)Then
959     !                      IF (W_G(LFLUID) <= ZERO) THEN
960                              Scalar(IJK, N) = Scalar(LFLUID, N)
961     !                     ENDIF
962                          Else
963     !                      IF (W_s(LFLUID, M) <= ZERO) THEN
964                              Scalar(IJK, N) = Scalar(LFLUID, N)
965     !                      ENDIF
966                          Endif
967                       END DO
968     
969                       IF(K_Epsilon) THEN
970     !                    IF (W_G(LFLUID) <= ZERO) THEN
971                            K_Turb_G(IJK) = K_Turb_G(LFLUID)
972                            E_Turb_G(IJK) = E_Turb_G(LFLUID)
973     !                   ENDIF
974                       ENDIF
975     
976                       DO M = 1, MMAX
977                          P_S(IJK,M) = P_S(LFLUID,M)
978     !                     IF (W_S(IJK,M) <= ZERO) THEN
979                             ROP_S(IJK,M) = ROP_S(LFLUID,M)
980                             T_S(IJK,M) = T_S(LFLUID,M)
981     !                     ELSE
982     !                        ROP_S(IJK,M) = ZERO
983     !                     ENDIF
984     !
985                          IF(BC_ROP_S(BCV,M)/=UNDEFINED)ROP_S(IJK,M)=BC_ROP_S(BCV,M)
986     !
987                          IF(BC_EP_G(BCV)==UNDEFINED)EP_G(IJK)=EP_G(IJK)-EP_S(IJK,M)
988                          N = 1
989                          IF (NMAX(M) > 0) THEN
990                             X_S(IJK,M,:NMAX(M)) = X_S(LFLUID,M,:NMAX(M))
991                             N = NMAX(M) + 1
992                          ENDIF
993                       END DO
994                       ROP_G(IJK) = RO_G(IJK)*EP_G(IJK)
995     !                  U_G(IJK) = U_G(LFLUID)
996     !                  V_G(IJK) = V_G(LFLUID)
997     !                  IF (W_G(IJK) == UNDEFINED) THEN
998     !                     IF (ROP_G(IJK) > ZERO) THEN
999     !                        W_G(IJK) = ROP_G(LFLUID)*W_G(LFLUID)/ROP_G(IJK)
1000     !                     ELSE
1001     !                        W_G(IJK) = ZERO
1002     !                     ENDIF
1003     !                  ENDIF
1004     !                  Flux_gE(IJK) = Flux_gE(LFLUID)
1005     !                 Flux_gN(IJK) = Flux_gN(LFLUID)
1006     !                 Flux_gT(IJK) = Flux_gT(LFLUID)
1007     
1008     !                  IF (MMAX > 0) THEN
1009     !                     U_S(IJK,:MMAX) = U_S(LFLUID,:MMAX)
1010     !                     V_S(IJK,:MMAX) = V_S(LFLUID,:MMAX)
1011     !                     WHERE (W_S(IJK,:MMAX) == UNDEFINED) W_S(IJK,:MMAX) = W_S(&
1012     !                        LFLUID,:MMAX)
1013     !                     Flux_sE(IJK,:MMAX) = Flux_sE(LFLUID,:MMAX)
1014     !                    Flux_sN(IJK,:MMAX) = Flux_sN(LFLUID,:MMAX)
1015     !                    Flux_sT(IJK,:MMAX) = Flux_sT(LFLUID,:MMAX)
1016     !                  ENDIF
1017                    ENDIF
1018     
1019                 ENDIF
1020              ENDIF
1021     
1022           ENDIF
1023           END DO
1024     
1025     !      print*,'bottom of cg_set_outflow'
1026           RETURN
1027           END SUBROUTINE CG_SET_OUTFLOW
1028