get_eq.f

Go to the documentation of this file.

!vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
!                                                                      C
!  Module name: GET_EQ(A, B, VEL, DTxFA, IJK1)                         C
!  Purpose: Complete the filling of coefficient matrix A               C
!                                                                      C
!  Author: M. Syamlal                                 Date: 23-JUL-92  C
!  Reviewer: W. Rogers                                Date: 11-DEC-92  C
!                                                                      C
!  Revision Number:                                                    C
!  Purpose:                                                            C
!  Author:                                            Date: dd-mmm-yy  C
!  Reviewer:                                          Date: dd-mmm-yy  C
!                                                                      C
!  Literature/Document References:                                     C
!                                                                      C
!  Variables referenced:                                               C
!  Variables modified:                                                 C
!                                                                      C
!  Local variables:                                                    C
!                                                                      C
!^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
!
      SUBROUTINE get_eq(A, BB, VEL, DTXFA, IJK1)
!...Translated by Pacific-Sierra Research VAST-90 2.06G5  12:17:31  12/09/98
!...Switches: -xf
!
!  Include param.inc file to specify parameter values
!
!-----------------------------------------------
!   M o d u l e s
!-----------------------------------------------
      USE param
      USE param1
      USE physprop
      USE geometry
      USE fldvar
      USE indices
      IMPLICIT NONE
!-----------------------------------------------
!   G l o b a l   P a r a m e t e r s
!-----------------------------------------------
!-----------------------------------------------
!   D u m m y   A r g u m e n t s
!-----------------------------------------------
!
!                      Local index of the solids phases
      INTEGER          L, M
!
!                      IMJK etc.
      INTEGER          IJK1
!
!                      Average DTxF_gs and DTxF_ss
      DOUBLE PRECISION DTxFA(0:dimension_m, 0:dimension_m)
!
!                      Coefficients in the discretized momentum equations
!                        Av = B, where v(0) = V_g, v(1) = V_s1 etc.
      DOUBLE PRECISION A(0:dimension_m, 0:dimension_m), BB(0:dimension_m)
!
!                      Solids velocity array passed
      DOUBLE PRECISION VEL (dimension_3, dimension_m)
!
!-----------------------------------------------
!
!
!     Compute off-diagonal coefficients
!
      DO m = 1, mmax
         IF (a(m,m) /= zero) THEN
            a(0,m) = -dtxfa(0,m)
            a(m,0) = -dtxfa(0,m)
            DO l = m + 1, mmax
               IF (a(l,l) /= zero) THEN
                  a(l,m) = -dtxfa(l,m)
                  a(m,l) = -dtxfa(l,m)
               ELSE
                  a(l,m) = zero
                  a(m,l) = zero
               ENDIF
            END DO
         ELSE
            a(0,m) = zero
            a(m,0) = zero
            l = m + 1
            IF (mmax - m > 0) THEN
               a(m+1:mmax,m) = zero
               a(m,m+1:mmax) = zero
               l = mmax + 1
            ENDIF
            IF (bb(m) == undefined) THEN
               a(0,0) = a(0,0) + dtxfa(0,m)
               bb(0) = bb(0) + dtxfa(0,m)*vel(ijk1,m)
               DO l = 1, mmax
                  IF (l/=m .AND. a(l,l)/=zero) THEN
                     a(l,l) = a(l,l) + dtxfa(l,m)
                     bb(l) = bb(l) + dtxfa(l,m)*vel(ijk1,m)
                  ENDIF
               END DO
            ENDIF
         ENDIF
      END DO
      DO m = 0, mmax
         DO l = 0, mmax
            IF (l /= m) a(m,m) = a(m,m) - a(m,l)
         END DO
      END DO
      RETURN
      END SUBROUTINE get_eq