File: RELATIVE:/../../../mfix.git/model/check_convergence.f

1     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
2     !                                                                      C
3     !  Subroutine: CHECK_CONVERGENCE                                       C
4     !  Purpose: Monitor convergence                                        C
5     !                                                                      C
6     !                                                                      C
7     !  Author: M. Syamlal                                 Date: 8-JUL-96   C
8     !  Reviewer:                                          Date:            C
9     !                                                                      C
10     !                                                                      C
11     !  Literature/Document References:                                     C
12     !                                                                      C
13     !  Variables referenced:                                               C
14     !  Variables modified:                                                 C
15     !  Local variables:                                                    C
16     !                                                                      C
17     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
18     
19           SUBROUTINE CHECK_CONVERGENCE(NIT, errorpercent, MUSTIT)
20     
21     !-----------------------------------------------
22     ! Modules
23     !-----------------------------------------------
24           USE geometry
25           USE indices
26           USE mpi_utility
27           USE param
28           USE param1
29           USE physprop
30           USE residual
31           USE run
32           USE scalars, only :NScalar
33           USE toleranc
34           USE utilities, ONLY: check_vel_bound
35           IMPLICIT NONE
36     !-----------------------------------------------
37     ! Local parameters
38     !-----------------------------------------------
39     ! Maximum % error allowed in fluid continuity
40     !      DOUBLE PRECISION, PARAMETER :: MaxErrorPercent = 1.0E-6
41     !-----------------------------------------------
42     ! Dummy arguments
43     !-----------------------------------------------
44     ! Iteration number
45           INTEGER, INTENT(IN) :: NIT
46     ! %error in fluid mass balance
47           DOUBLE PRECISION, INTENT(IN) :: errorpercent
48     ! value tells whether to iterate (1) or not (0).
49           INTEGER, INTENT(INOUT) :: MUSTIT
50     !-----------------------------------------------
51     ! Local variables
52     !-----------------------------------------------
53     ! sum of residuals
54           DOUBLE PRECISION :: SUM, SUM_T, SUM_X, SUM_Th
55     ! max of residuals
56           DOUBLE PRECISION :: maxres
57     ! index
58           INTEGER :: L, M, maxL, maxM, N, maxN
59     ! to indicate undefined residual in species eq at the
60     ! beginning of iterations
61           LOGICAL :: NO_RESID
62     !-----------------------------------------------
63     
64     ! sum the residuals from correction equation (pressure and/or
65     ! solids), continuity equations (gas and/or solids) and momentum
66     ! equations (gas and solids)
67     
68     ! add pressure correction residual
69     !      if(abs(errorpercent) > MaxErrorPercent)then
70             SUM = RESID(RESID_P,0)
71     !      else
72     !        SUM = zero
73     !      endif
74     
75     
76     ! add solids correction residual
77           IF(MMAX > 0) SUM = SUM + RESID(RESID_P,1)
78     
79     ! add continuity equation residuals
80           DO M = 0, MMAX
81              SUM = SUM + RESID(RESID_RO,M)
82           ENDDO
83     ! add momentum equation residuals
84           DO M = 0, MMAX
85              SUM = SUM + RESID(RESID_U,M)
86           ENDDO
87           DO M = 0, MMAX
88              SUM = SUM + RESID(RESID_V,M)
89           ENDDO
90           IF (DO_K) THEN
91              DO M = 0, MMAX
92                 SUM = SUM + RESID(RESID_W,M)
93              ENDDO
94           ENDIF
95     !      call global_all_sum(SUM)
96     
97     
98     ! sum the granular energy equation residuals
99           SUM_Th = zero
100           IF (GRANULAR_ENERGY) THEN
101              DO M = 1, MMAX
102                 SUM_Th = SUM_Th + RESID(RESID_TH,M)
103              ENDDO
104           ENDIF
105     
106     
107     ! sum the energy equation residuals
108           SUM_T = ZERO
109           IF (ENERGY_EQ) THEN
110              DO M = 0, MMAX
111                 SUM_T = SUM_T + RESID(RESID_T,M)
112              END DO
113           ENDIF
114     !      call global_all_sum(SUM_T)
115     
116     ! sum the species equation residuals
117           SUM_X = ZERO
118           NO_RESID = .FALSE.
119           DO M = 0, MMAX
120              IF (SPECIES_EQ(M)) THEN
121                 DO N = 1, NMAX(M)
122                    IF (RESID(RESID_X+(N-1),M) == UNDEFINED) NO_RESID = .TRUE.
123                    SUM_X = SUM_X + RESID(RESID_X+(N-1),M)
124                 ENDDO
125              ENDIF
126           ENDDO
127     !      call global_all_sum(SUM_X)
128           IF (NO_RESID) SUM_X = TOL_RESID_X + ONE
129     
130     
131     ! find the variable with maximum residual
132           IF (RESID_INDEX(MAX_RESID_INDEX,1) == UNDEFINED_I) THEN
133              MAXRES = ZERO
134              DO L = 1, NRESID
135                 DO M = 0, MMAX
136                    IF (RESID(L,M) >= MAXRES) THEN
137                       MAXRES = RESID(L,M)
138                       MAXL = L
139                       MAXM = M
140                       IF (L >= RESID_X) THEN
141                          MAXN = L - RESID_X + 1
142                       ELSE
143                          MAXN = UNDEFINED_I
144                       ENDIF
145                    ENDIF
146                 END DO
147              END DO
148              IF (MAXN == UNDEFINED_I) THEN
149                 WRITE (RESID_STRING(MAX_RESID_INDEX), '(A1,I1)') RESID_PREFIX(MAXL)&
150                    , MAXM
151              ELSE
152                 WRITE (RESID_STRING(MAX_RESID_INDEX), '(A1,I1,I2.0)') 'X', MAXM, &
153                    MAXN
154              ENDIF
155           ENDIF
156     
157           IF (GROUP_RESID) THEN
158              RESID_GRP(HYDRO_GRP) = SUM
159              IF(GRANULAR_ENERGY) RESID_GRP(THETA_GRP) = SUM_TH
160              IF(ENERGY_EQ) RESID_GRP(ENERGY_GRP) = SUM_T
161              IF(ANY_SPECIES_EQ) RESID_GRP(SPECIES_GRP) = SUM_X
162              IF(NScalar > 0) RESID_GRP(SCALAR_GRP) = RESID(RESID_sc,0)
163              IF(K_EPSILON) RESID_GRP(KE_GRP) = RESID(RESID_ke,0)
164           ENDIF
165     
166     ! Flag set in interactive mode to bypass time-step advancement.
167           IF(INTERACTIVE_NITS /= UNDEFINED_I) THEN
168              MUSTIT = 1
169              RETURN
170           ENDIF
171     
172     ! Every 5 iterations detect whether the run is stalled by checking
173     ! that the total residual has decreased.
174           IF(DETECT_STALL .AND. MOD(NIT,5) == 0) THEN
175              IF(NIT > 10) THEN
176                 IF(SUM5_RESID <= SUM) THEN
177     ! The run is stalled. Reduce the time step.
178                    IF(.NOT.PERSISTENT_MODE) THEN
179                       MUSTIT = 2
180                       RETURN
181     ! Forces the max number of iterations for DT=DT_MIN
182                    ELSEIF(DT > DT_MIN) THEN
183                       MUSTIT = 1
184                       RETURN
185                    ENDIF
186                 ENDIF
187              ENDIF
188              SUM5_RESID = SUM
189           ENDIF
190     
191     ! Require at least two iterations.
192           IF(NIT == 1) THEN
193              MUSTIT = 1
194              RETURN
195           ENDIF
196     
197     ! total residual
198           IF(SUM<=TOL_RESID .AND. SUM_T<=TOL_RESID_T .AND. &
199              RESID(RESID_sc,0)<=TOL_RESID_Scalar .AND. SUM_X<=TOL_RESID_X &
200             .AND. RESID(RESID_ke,0)<=TOL_RESID_K_Epsilon &
201             .AND. SUM_Th <=TOL_RESID_Th)THEN
202              MUSTIT = 0                              !converged
203           ELSEIF (SUM>=TOL_DIVERGE .OR. SUM_T>=TOL_DIVERGE .OR.&
204                   RESID(RESID_sc,0)>= TOL_DIVERGE .OR. SUM_X>=TOL_DIVERGE&
205                   .OR. RESID(RESID_ke,0)>= TOL_DIVERGE &
206                   .OR. SUM_Th >= TOL_DIVERGE ) THEN
207              IF (NIT /= 1) THEN
208                 MUSTIT = 2                           !diverged
209              ELSE
210                 MUSTIT = 1                           !not converged
211              ENDIF
212           ELSE
213              MUSTIT = 1                              !not converged
214           ENDIF
215     
216     
217     ! to check upper bound (speed of sound) limit for gas and
218     ! solids velocity components
219     ! only check velocity if any of the momentum equations are solved
220           DO M = 0,MMAX
221              IF (MOMENTUM_X_EQ(M) .OR. MOMENTUM_Y_EQ(M) .OR. &
222                  MOMENTUM_Z_EQ(M)) THEN
223                 IF(CHECK_VEL_BOUND()) MUSTIT = 2     !divergence
224                 EXIT  ! only need to call check_vel_bound once
225              ENDIF
226           ENDDO
227     
228     
229           RETURN
230           END SUBROUTINE CHECK_CONVERGENCE
231     
232     
233