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