check_convergence.f

Go to the documentation of this file.

!vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
!                                                                      C
!  Subroutine: CHECK_CONVERGENCE                                       C
!  Purpose: Monitor convergence                                        C
!                                                                      C
!                                                                      C
!  Author: M. Syamlal                                 Date: 8-JUL-96   C
!  Reviewer:                                          Date:            C
!                                                                      C
!                                                                      C
!  Literature/Document References:                                     C
!                                                                      C
!  Variables referenced:                                               C
!  Variables modified:                                                 C
!  Local variables:                                                    C
!                                                                      C
!^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C

      SUBROUTINE check_convergence(NIT, errorpercent, MUSTIT)

!-----------------------------------------------
! Modules
!-----------------------------------------------
      USE geometry
      USE indices
      USE mpi_utility
      USE param
      USE param1
      USE physprop
      USE residual
      USE run
      USE scalars, only :nscalar
      USE toleranc
      USE utilities, ONLY: check_vel_bound
      IMPLICIT NONE
!-----------------------------------------------
! Local parameters
!-----------------------------------------------
! Maximum % error allowed in fluid continuity
!      DOUBLE PRECISION, PARAMETER :: MaxErrorPercent = 1.0E-6
!-----------------------------------------------
! Dummy arguments
!-----------------------------------------------
! Iteration number
      INTEGER, INTENT(IN) :: NIT
! %error in fluid mass balance
      DOUBLE PRECISION, INTENT(IN) :: errorpercent
! value tells whether to iterate (1) or not (0).
      INTEGER, INTENT(INOUT) :: MUSTIT
!-----------------------------------------------
! Local variables
!-----------------------------------------------
! sum of residuals
      DOUBLE PRECISION :: SUM, SUM_T, SUM_X, SUM_Th
! max of residuals
      DOUBLE PRECISION :: maxres
! index
      INTEGER :: L, M, maxL, maxM, NN, maxN
! to indicate undefined residual in species eq at the
! beginning of iterations
      LOGICAL :: NO_RESID
!-----------------------------------------------

! sum the residuals from correction equation (pressure and/or
! solids), continuity equations (gas and/or solids) and momentum
! equations (gas and solids)

! add pressure correction residual
!      if(abs(errorpercent) > MaxErrorPercent)then
        sum = resid(resid_p,0)
!      else
!        SUM = zero
!      endif


! add solids correction residual
      IF(mmax > 0) sum = sum + resid(resid_p,1)

! add continuity equation residuals
      DO m = 0, mmax
         sum = sum + resid(resid_ro,m)
      ENDDO
! add momentum equation residuals
      DO m = 0, mmax
         sum = sum + resid(resid_u,m)
      ENDDO
      DO m = 0, mmax
         sum = sum + resid(resid_v,m)
      ENDDO
      IF (do_k) THEN
         DO m = 0, mmax
            sum = sum + resid(resid_w,m)
         ENDDO
      ENDIF
!      call global_all_sum(SUM)


! sum the granular energy equation residuals
      sum_th = zero
      IF (granular_energy) THEN
         DO m = 1, mmax
            sum_th = sum_th + resid(resid_th,m)
         ENDDO
      ENDIF


! sum the energy equation residuals
      sum_t = zero
      IF (energy_eq) THEN
         DO m = 0, mmax
            sum_t = sum_t + resid(resid_t,m)
         END DO
      ENDIF
!      call global_all_sum(SUM_T)

! sum the species equation residuals
      sum_x = zero
      no_resid = .false.
      DO m = 0, mmax
         IF (species_eq(m)) THEN
            DO nn = 1, nmax(m)
               IF (resid(resid_x+(nn-1),m) == undefined) no_resid = .true.
               sum_x = sum_x + resid(resid_x+(nn-1),m)
            ENDDO
         ENDIF
      ENDDO
!      call global_all_sum(SUM_X)
      IF (no_resid) sum_x = tol_resid_x + one


! find the variable with maximum residual
      IF (resid_index(max_resid_index,1) == undefined_i) THEN
         maxres = zero
         DO l = 1, nresid
            DO m = 0, mmax
               IF (resid(l,m) >= maxres) THEN
                  maxres = resid(l,m)
                  maxl = l
                  maxm = m
                  IF (l >= resid_x) THEN
                     maxn = l - resid_x + 1
                  ELSE
                     maxn = undefined_i
                  ENDIF
               ENDIF
            END DO
         END DO
         IF (maxn == undefined_i) THEN
            WRITE (resid_string(max_resid_index), '(A1,I1)') resid_prefix(maxl)&
               , maxm
         ELSE
            WRITE (resid_string(max_resid_index), '(A1,I1,I2.0)') 'X', maxm, &
               maxn
         ENDIF
      ENDIF

      IF (group_resid) THEN
         resid_grp(hydro_grp) = sum
         IF(granular_energy) resid_grp(theta_grp) = sum_th
         IF(energy_eq) resid_grp(energy_grp) = sum_t
         IF(any_species_eq) resid_grp(species_grp) = sum_x
         IF(nscalar > 0) resid_grp(scalar_grp) = resid(resid_sc,0)
         IF(k_epsilon) resid_grp(ke_grp) = resid(resid_ke,0)
      ENDIF

! Every 5 iterations detect whether the run is stalled by checking
! that the total residual has decreased.
      IF(detect_stall .AND. mod(nit,5) == 0) THEN
         IF(nit > 10) THEN
            IF(sum5_resid <= sum) THEN
! The run is stalled. Reduce the time step.
               IF(.NOT.persistent_mode) THEN
                  mustit = 2
                  RETURN
! Forces the max number of iterations for DT=DT_MIN
               ELSEIF(dt > dt_min) THEN
                  mustit = 1
                  RETURN
               ENDIF
            ENDIF
         ENDIF
         sum5_resid = sum
      ENDIF

! Require at least two iterations.
      IF(nit == 1) THEN
         mustit = 1
         RETURN
      ENDIF

! total residual
      IF(sum<=tol_resid .AND. sum_t<=tol_resid_t .AND. &
         resid(resid_sc,0)<=tol_resid_scalar .AND. sum_x<=tol_resid_x &
        .AND. resid(resid_ke,0)<=tol_resid_k_epsilon &
        .AND. sum_th <=tol_resid_th)THEN
         mustit = 0                              !converged
      ELSEIF (sum>=tol_diverge .OR. sum_t>=tol_diverge .OR.&
              resid(resid_sc,0)>= tol_diverge .OR. sum_x>=tol_diverge&
              .OR. resid(resid_ke,0)>= tol_diverge &
              .OR. sum_th >= tol_diverge ) THEN
         IF (nit /= 1) THEN
            mustit = 2                           !diverged
         ELSE
            mustit = 1                           !not converged
         ENDIF
      ELSE
         mustit = 1                              !not converged
      ENDIF


! to check upper bound (speed of sound) limit for gas and
! solids velocity components
! only check velocity if any of the momentum equations are solved
      DO m = 0,mmax
         IF (momentum_x_eq(m) .OR. momentum_y_eq(m) .OR. &
             momentum_z_eq(m)) THEN
            IF(check_vel_bound()) mustit = 2     !divergence
            EXIT  ! only need to call check_vel_bound once
         ENDIF
      ENDDO


      RETURN
      END SUBROUTINE check_convergence