Esope source files

   1 : C FRUSBM    SOURCE    CHAT      05/01/13    00:09:18     5004           
   2 :       SUBROUTINE FRUSBM(NESP,
   3 :      &                  GAMG,ROG,PG,UNG,UTG,
   4 :      &                  GAMD,ROD,PD,UND,UTD,
   5 :      &                  YG,YD,VINF,FLU,
   6 :      &                  CELLT)
   7 : C
   8 : C PROJET            :  CASTEM 2000
   9 : C
  10 : C NOM               :  FRUSBM
  11 : C
  12 : C DESCRIPTION       :  Formulation Volumes Finis  pour les Equations 
  13 : C                      d'Euler Multi-Especes relatives à un melange 
  14 : C                      de gaz parfaits.
  15 : C
  16 : C                      Calcul du flux aux interfaces avec la methode
  17 : C                      de Rusanov preconditionée
  18 : C
  19 : C LANGUAGE          :  FORTRAN 77 
  20 : C      
  21 : C AUTEUR            :  A. BECCANTINI SFME/LTMF
  22 : C
  23 : C************************************************************************
  24 : C
  25 : C APPELES           :  FRBM1
  26 : C
  27 : C************************************************************************
  28 : C
  29 : C**** Entrées:
  30 : C
  31 : C     NESP            =  nombre d'especes considérées dans les Equations
  32 : C                        d'Euler
  33 : C
  34 : C     GAMG, GAMD      =  les "gamma" du gaz (gauche et droite)
  35 : C
  36 : C     ROG, ROD        =  les densités 
  37 : C
  38 : C     PG, PD          =  les pressions
  39 : C
  40 : C     UNG, UND        =  vitesses normales
  41 : C
  42 : C     UTG, UTD        =  vitesses tangentielles
  43 : C
  44 : C     VINF            = vitesse de cut-off
  45 : C
  46 : C**** Sorties:
  47 : C
  48 : C     FLU             =  table du flux a l'interface dans le repaire
  49 : C                        (n,t), i.e.
  50 : C                        (rho*un, rho*un*un + p, rho*un*ut, rho*un*ht,
  51 : C                         rho*un*y1, ...)
  52 : C
  53 : C     CELLT           =  condition de stabilité, i.e.
  54 : C
  55 : C                        dT/diamax < cellt
  56 : C      
  57 : C************************************************************************
  58 : C
  59 : C HISTORIQUE (Anomalies et modifications éventuelles)
  60 : C
  61 : C HISTORIQUE :  Créé le 08.02.00
  62 : C
  63 : C               21.02.00 transport de scalaires passifs
  64 : C
  65 : C************************************************************************
  66 : C
  67 : C N.B.: Toutes les variables sont DECLAREES
  68 : C                          
  69 : C
  70 :       IMPLICIT INTEGER(I-N)
  71 :       INTEGER NESP,I1
  72 :       REAL*8  GAMG,ROG,PG,UNG,UTG
  73 :      &     ,GAMD,ROD,PD,UND,UTD
  74 :      &     ,YG(*),YD(*),VINF
  75 :      &     ,FLU(*),CELLT, FLUM
  76 :      &     ,UNF,UTF,ROF,PF,GAMF
  77 :      &     ,AMAT(4,4)
  78 :      &     ,CSON,RETG,RETD,DRO,DRUN,DRUT,DRET,RLAMMA
  79 : C
  80 : C**** Etat moyen
  81 : C
  82 : C     SIGN(4.0D0,-1.0D0)=-4
  83 :       UNF=0.5D0*(UNG+UND)
  84 :       UTF=0.5D0*(UTG+UTD)
  85 :       ROF=0.5D0*(ROG+ROD)
  86 :       PF=0.5D0*(PG+PD)
  87 :       GAMF=0.5D0*(GAMG+GAMD)
  88 : C
  89 : C**** We include in the cut-off UNG,UND,UTG,UTD in order to
  90 : C     avoid low diffusivity in stagnation regions
  91 : C
  92 :       VINF=MAX(VINF,((UNG*UNG+UTG*UTG)**0.5D0))
  93 :       VINF=MAX(VINF,((UND*UND+UTD*UTD)**0.5D0))
  94 : C
  95 : C**** Calcule de la matrice de diffusivité pour le variables
  96 : C     conservatives.
  97 : C     NB: F(UL,Ur) = 0.5 (F(UL)+F(UR)+AMAT*RLAMMA*(UL-UR))
  98 : C
  98 :       CALL FRBM1(VINF, ROF,UNF,UTF,PF,GAMF,AMAT, RLAMMA)
  99 : C
 100 : C**** Vitesse de l'onde la plus rapide à l'interface (la vrai!)
 101 : C
 102 :       CSON= (GAMF*PF/ROF)**0.5D0
 103 :       CELLT=1.0D0/(CSON+ABS(UNF))
 104 : C
 105 : C**** DELTA
 106 : C
 107 :       RETG = ((1.0D0 / (GAMG - 1.0D0)) * PG) 
 108 :      &     + (0.5D0*ROG*((UNG*UNG)+(UTG*UTG)))
 109 :       RETD = ((1.0D0 / (GAMD - 1.0D0)) * PD) 
 110 :      &     + (0.5D0*ROD*((UND*UND)+(UTD*UTD)))
 111 :       DRO = ROG - ROD
 112 :       DRUN = (ROG * UNG) - (ROD * UND)
 113 :       DRUT = (ROG * UTG) - (ROD * UTD)
 114 :       DRET = RETG - RETD
 115 : C
 116 :       FLU(1) = 0.5D0 * ((ROG * UNG) + (ROD * UND) + 
 117 :      &     (((AMAT(1,1) * DRO) + (AMAT(1,2) * DRUN) + (AMAT(1,3) * DRUT) 
 118 :      &     + (AMAT(1,4) * DRET)) * RLAMMA))
 119 :       FLU(2) = 0.5D0 * ((ROG * UNG * UNG) + (ROD * UND * UND) + 
 120 :      &     (PG + PD) + 
 121 :      &     (((AMAT(2,1) * DRO) + (AMAT(2,2) * DRUN) + (AMAT(2,3) * DRUT) 
 122 :      &     +(AMAT(2,4) * DRET)) * RLAMMA))
 123 :       FLU(3) = 0.5D0 * ((ROG * UNG * UTG) + (ROD * UND * UTD) + 
 124 :      &     (((AMAT(3,1) * DRO) + (AMAT(3,2) * DRUN) + (AMAT(3,3) * DRUT) 
 125 :      &     +(AMAT(3,4) * DRET))*RLAMMA))
 126 :       FLU(4) = 0.5D0 * ((UNG * (RETG + PG)) + (UND * (RETD + PD)) +
 127 :      &     (((AMAT(4,1) * DRO) + (AMAT(4,2) * DRUN) + (AMAT(4,3) * DRUT) 
 128 :      &     +(AMAT(4,4) * DRET))*RLAMMA))
 129 : C
 130 :       FLUM = FLU(1)
 131 :       IF(FLUM .GT. 0)THEN
 132 :          DO I1 = 1, NESP, 1
 133 :             FLU(4+I1)=FLUM * YG(I1)
 134 :          ENDDO
 135 :       ELSE
 136 :          DO I1 = 1, NESP, 1
 137 :             FLU(4+I1)=FLUM * YD(I1)
 138 :          ENDDO
 139 :       ENDIF
 140 : C      
 141 :       RETURN
 142 :       END
 143 :  
 144 :