Télécharger racsti.eso

Retour à la liste

racsti
C RACSTI    SOURCE    BECC      11/05/26    21:16:16     6981
      SUBROUTINE RACSTI(
     &     PC_L, GAM_L,
     &     PC_R, GAM_R,
     &     PMIN,
     &     RHO_L, P_L, U_L,
     &     RHO_R, P_R, U_R,
     &     D_L, D_B, D_A, D_R,
     &     RHO_B, P_B, U_B,
     &     RHO_A, P_A, U_A,
     &     LOGDEB, LOGAN, LOGNC, LOGVAC)
*
*************************************************************************
*
* project           :
*
* name              :  racsti
*
* description       :  euler equations for a mixture of stiffened gases
*                      flux in the non-reactive case.
*
*                      field-by-field decomposition (entropy-respecting).
*
* language          :  fortran 77
*
* author            :  a. beccantini den/dm2s/sfme/ltmf
*
*************************************************************************
*
* called by         :
*
*
*************************************************************************
*
***** input
*
*     pc_l, gam_l     =  properties of the gas in the left
*
*     pc_r, gam_r     =  properties of the gas in the right
*
*     rho_lr, p_lr, w_lr
*                     = density, pressure, velocity
*
*     pmin            = if p < pmin => vacuum
*
*     logdeb          = debugging ?
*
***** input / output
*
*     logvac          = vacuum ?
*
***** output
*
*     rho_ba, p_ba, w_ba
*                     = density, pressure, velocity
*
*     d_l, d_b, d_a, d_r
*                     = wave speeds d_l <= d_b <= d_a <= d_r
*
*     logan           = anomaly
*
*     lognc               =  if true, the Newton for the Riemann problem
*                            solution did not converge
*
*************************************************************************
*
*     26/11/2009      created
*     25/05/2011      evolution in CAST3M
*
*************************************************************************
*
* n.b.: all variables are declared
*
C      IMPLICIT NONE
      IMPLICIT INTEGER(I-N)
      INTEGER ITER, NITER
      PARAMETER (NITER=50)
      REAL*8
     &     PMIN,
     &     PC_L, GAM_L,
     &     PC_R, GAM_R,
     &     D_L, D_B, D_A, D_R,
     &     RHO_L, P_L, U_L,
     &     RHO_B, P_B, U_B,
     &     RHO_A, P_A, U_A,
     &     RHO_R, P_R, U_R
     &     , EPSERR
      PARAMETER (EPSERR = 1.0D-6)
*
      REAL*8 P_VAC, A_L, A_R
     &     , P_ES, ROES_B, UES_B, DES_B, DER_B
     &     , ROES_A, UES_A, DES_A, DER_A
     &     , P_INT
*
*     debugging ?
*
      LOGICAL LOGDEB, LOGAN, LOGNC, LOGVAC
*
      LOGNC = .TRUE.
*
***** initialisation of the speed of the wave on the left
*     and on the right
*
      CALL ASTIFF(RHO_L, P_L, GAM_L, PC_L, A_L, LOGDEB, LOGAN)
      CALL ASTIFF(RHO_R, P_R, GAM_R, PC_R, A_R, LOGDEB, LOGAN)
      IF ( LOGAN ) THEN
         WRITE(*,*) 'SUBROUTINE RACSTI.F'
         WRITE(*,*) 'INITIALISATION OF THE WAVE SPEEDS'
         WRITE(*,*) 'ANOMALY DETECTED'
         GOTO 9999
      ENDIF
      D_L = U_L - A_L
      D_R = U_R + A_R
*
******** states
*
*         l
*        ***
*          |*                              r
*          | *   b                       *********
*          |  ********                  *|
*          |  |      |  *   a          * |
*          |  |      |    *************  |
*          |  |      |    |           |  |
*          |  |      |    |           |  |
*        d_l d_b   un_b   un_a       d_a d_r
*
*
*        if no vacuum, un_b = un_a
*
      P_VAC = MIN( P_L, P_R )
      P_VAC = MAX( P_VAC, PMIN)
      P_VAC = 1.0D-8 * P_VAC
      CALL FSTRAL( PC_L, GAM_L, P_VAC, RHO_L, P_L, U_L,
     &     U_B, RHO_B, D_B, LOGDEB, LOGAN)
      CALL FSTRAR(PC_R, GAM_R, P_VAC, RHO_R, P_R, U_R,
     &     U_A, RHO_A, D_A, LOGDEB, LOGAN)
      IF ( LOGAN ) THEN
         WRITE(*,*) 'SUBROUTINE RACSTI.F'
         WRITE(*,*) 'COMPUTATION OF THE VACUUM STATE'
         WRITE(*,*) 'ANOMALY DETECTED'
         GOTO 9999
      ENDIF
*
      IF (LOGVAC) THEN
*
*        From elsewhere, we know that we are in the vacuum...
*
         P_B = P_VAC
         P_A = P_VAC
         LOGNC = .FALSE.
 
C
      ELSEIF ( U_B .LT. U_A) THEN
*
*        vacuum
*
         P_B = P_VAC
         P_A = P_VAC
         LOGVAC = .TRUE.
         LOGNC = .FALSE.
      ELSE
         LOGVAC = .FALSE.
*
*************************************************************************
***** intersection (states l, b, a, r) **********************************
*************************************************************************
*
*
***** initialization of p_int and of the states b, a
*
*
         P_INT = 0.5D0 * (P_R + P_L)
*
         CALL FSTERL( PC_L, GAM_L, P_INT, RHO_L, P_L, U_L,
     &        U_B, RHO_B, D_B, LOGDEB, LOGAN)
         CALL FSTERR( PC_R, GAM_R, P_INT, RHO_R, P_R, U_R,
     &        U_A, RHO_A, D_A, LOGDEB, LOGAN)
*
         IF ( LOGAN ) THEN
            WRITE(*,*) 'SUBROUTINE RACSTI.F'
            WRITE(*,*) 'INITIALIZATION OF THE STATES'
            WRITE(*,*) 'ANOMALY DETECTED'
            GOTO 9999
         ENDIF
*
         DO ITER = 1, NITER, 1
*
*           pmin
*
            P_VAC = 1D-8 * MIN (P_L, P_R, P_INT)
            P_VAC = MAX( P_VAC, PMIN)
*
            IF (.TRUE.) THEN
*
*           exact evaluation of the derivative
*
               CALL DERL(PC_L, GAM_L, P_INT, RHO_L, P_L, U_L, DER_B,
     $              LOGDEB, LOGAN)
               IF ( LOGAN ) THEN
                  WRITE(*,*) 'SUBROUTINE RACSTI.F'
                  WRITE(*,*) 'EVALUATION OF THE DERIVATES'
                  WRITE(*,*) 'ANOMALY DETECTED'
                  GOTO 9999
               ENDIF
               CALL DERR(PC_R, GAM_R, P_INT, RHO_R, P_R, U_R, DER_A,
     $              LOGDEB, LOGAN)
*
C               WRITE(*,*) 'DER_B, DER_A ex = ', DER_B, DER_A
C
               IF ( LOGAN ) THEN
                  WRITE(*,*) 'SUBROUTINE RACSTI.F'
                  WRITE(*,*) 'EVALUATION OF THE DERIVATES'
                  WRITE(*,*) 'ANOMALY DETECTED'
                  GOTO 9999
               ENDIF
            ELSE
*
*           numerical evaluation of the derivatives
*
               P_ES = 1.001D0 * P_INT
               CALL FSTERL( PC_L, GAM_L, P_ES, RHO_L, P_L, U_L,
     &              UES_B, ROES_B, DES_B, LOGDEB, LOGAN)
               IF ( LOGAN ) THEN
                  WRITE(*,*) 'SUBROUTINE RACSTI.F'
                  WRITE(*,*) 'EVALUATION OF THE DERIVATES'
                  WRITE(*,*) 'ANOMALY DETECTED'
                  GOTO 9999
               ENDIF
               DER_B =  (UES_B - U_B) / (P_ES - P_INT)
               CALL FSTERR( PC_R, GAM_R, P_ES, RHO_R, P_R, U_R,
     &              UES_A, ROES_A, DES_A, LOGDEB, LOGAN)
               IF ( LOGAN ) THEN
                  WRITE(*,*) 'SUBROUTINE RACSTI.F'
                  WRITE(*,*) 'EVALUATION OF THE DERIVATES'
                  WRITE(*,*) 'ANOMALY DETECTED'
                  GOTO 9999
               ENDIF
               DER_A =  (UES_A - U_A) / (P_ES - P_INT)
C               WRITE(*,*) 'DER_B, DER_A num = ', DER_B, DER_A
C               STOP
            ENDIF
            IF ((DER_B .GT. 0.0D0)) THEN
               WRITE(*,*) 'SUBROUTINE RACSTI.F'
               WRITE (*,*) 'DER_B =',  DER_B
               LOGAN = .TRUE.
               GOTO 9999
            ENDIF
            IF ((DER_A .LT. 0.0D0)) THEN
               WRITE(*,*) 'SUBROUTINE RACSTI.F'
               WRITE (*,*) 'DER_A =',  DER_A
               LOGAN = .TRUE.
               GOTO 9999
            ENDIF
            IF ((DER_B - DER_A) .GE. 0.0D0) THEN
               WRITE(*,*) 'SUBROUTINE RACSTI.F'
               WRITE (*,*) 'DER_B - DER_A =',  (DER_B - DER_A)
               LOGAN = .TRUE.
               GOTO 9999
            ENDIF
C
            P_ES = P_INT - (U_B - U_A) / (DER_B - DER_A)
            P_ES = MAX (P_ES, P_VAC)
            CALL FSTERL( PC_L, GAM_L, P_ES, RHO_L, P_L, U_L,
     &           U_B, RHO_B, D_B, LOGDEB, LOGAN)
            CALL FSTERR( PC_R, GAM_R, P_ES, RHO_R, P_R, U_R,
     &           U_A, RHO_A, D_A, LOGDEB, LOGAN)
            IF ( LOGAN ) THEN
               WRITE(*,*) 'SUBROUTINE RACSTI.F'
               WRITE(*,*) 'COMPUTATION OF THE INTERSECTION'
               WRITE(*,*) 'ANOMALY DETECTED'
               GOTO 9999
            ENDIF
            IF ((ABS (P_ES - P_INT) / P_INT) .LT. EPSERR) THEN
               LOGNC = .FALSE.
               P_A = P_ES
               P_B = P_ES
               IF (P_B .GT. P_L) THEN
                  D_L = D_B
               ENDIF
               IF (P_A .GT. P_R) THEN
                  D_R = D_A
               ENDIF
               IF (LOGDEB) THEN
                  WRITE(*,*) 'SUBROUTINE RACSTI.F'
                  WRITE(*,*) 'CONVERGENCE ACHIEVED'
                  WRITE(*,*) 'ITER =', ITER
               ENDIF
               GOTO 9999
            ENDIF
            P_INT = P_ES
            IF (LOGDEB) THEN
               WRITE(*,*) 'SUBROUTINE RACSTI'
               WRITE(*,*) 'P_INT = ', P_INT
               WRITE(*,*) 'RHO_B, U_B', RHO_B, U_B
               WRITE(*,*) 'RHO_A, U_A', RHO_A, U_A
            ENDIF
         ENDDO
         P_A = P_ES
         P_B = P_ES
         IF (LOGNC) THEN
            WRITE(*,*) 'SUBROUTINE RACSTI.F'
            WRITE(*,*) 'CONVERGENCE NOT ACHIEVED'
            WRITE(*,*) 'ITER =', ITER
            WRITE(*,*) 'ERROR =', (ABS (P_ES - P_INT) / P_INT)
         ENDIF
      ENDIF
*
 9999 RETURN
      END
 
 

© Cast3M 2003 - Tous droits réservés.
Mentions légales