Télécharger la plaquette Cast3M


fr en

Télécharger cli272.eso

Retour à la liste

Numérotation des lignes : 

cli272

  1. C CLI272    SOURCE    CB215821  20/11/25    13:20:53     10792          
  2.       SUBROUTINE CLI272(NSP,MELEMF,MELEMC,MELECB,MELEFC,MELRES,INORM,
  3.      &     ICHPVO,ICHPSU,LRECP,LRECV,
  4.      &     IROC,IVITC,IPC,IYC,ICHLIM,ILIINC,ILIINP,IJAC,IJACO)
  5. C************************************************************************
  6. C
  7. C PROJET            :  CASTEM 2000
  8. C
  9. C NOM               :  CLI112
  10. C
  11. C DESCRIPTION       :  Subroutine appellée par CLIM22
  12. C                      OPTION: 'INJE' -- Jacobian
  13. C
  14. C LANGAGE           :  FORTRAN 77 + ESOPE 2000 (avec estensions CISI)
  15. C
  16. C AUTEUR            :  S. Kudriakov, DEN/DM2S/SFME/LTMF
  17. C
  18. C************************************************************************
  19. C
  20. C APPELES (Calcul)  :
  21. C
  22. C************************************************************************
  23. C
  24. C HISTORIQUE (Anomalies et modifications éventuelles)
  25. C
  26. C HISTORIQUE :
  27. C
  28. C************************************************************************
  29. C
  30. C----------------------------------------------------
  31. C**** Variables de COOPTIO
  32. C----------------------------------------------------
  33. C      INTEGER IPLLB, IERPER, IERMAX, IERR, INTERR
  34. C     &        ,IOTER,   IOLEC,  IOIMP, IOCAR,  IOACQ
  35. C     &        ,IOPER,   IOSGB,  IOGRA, IOSAU,  IORES
  36. C     &        ,IECHO,   IIMPI,  IOSPI
  37. C     &        ,IDIM,    IFICLE, IPREFI
  38. C     &        ,MCOORD
  39. C     &        ,IFOMOD, NIFOUR, IFOUR, NSDPGE, IONIVE
  40. C     &        ,NGMAXY, IZROSF, ISOTYP, IOSCR,LTEXLU
  41. C     &        ,NORINC,NORVAL,NORIND,NORVAD
  42. C     &        ,NUCROU, IPSAUV
  43. C
  44.       IMPLICIT INTEGER(I-N)
  45.       INTEGER MELEMF,MELEMC,MELECB,INORM,ICHPVO,ICHPSU, IROC,IVITC,IPC
  46.      &     ,IGAMC,ICHLIM,ICEL,NFAC,IFAC,MELRES,IJACO
  47.      &     ,NGF,NGC,NLF,NLC,NLCB
  48.      &     ,ILIINC,ILIINP,IJAC,II,JJ,K
  49.      &     ,MP, NBEL, NBME, NBSOUS, NKID, NKMT, NMATRI, NP, NRIGE
  50.      &     ,NSP,I, IYC,J, LRECP,LRECV,KV
  51.       REAL*8 VOLU,SURF,RC,PC,UXC,UYC,GAMC,CNX,CNY
  52.      &     ,PSRF,RHOUF,DECIDP,DUNDPC,DUXDPC,DUYDPC
  53.      &     ,UN,RHO,P,COEF
  54.      &     ,BR1,BOT,TOP
  55.  
  56.       CHARACTER*(8) TYPE
  57.  
  58. -INC PPARAM
  59. -INC CCOPTIO
  60. -INC SMLMOTS
  61. -INC SMELEME
  62.       POINTEUR MELEFC.MELEME
  63. -INC SMLENTI
  64.       POINTEUR MLEMC.MLENTI, MLEMCB.MLENTI,MLEMF.MLENTI
  65. -INC SMCHPOI
  66.       POINTEUR MPNORM.MPOVAL, MPVOL.MPOVAL, MPSURF.MPOVAL, MPRC.MPOVAL,
  67.      &     MPVC.MPOVAL, MPPC.MPOVAL, MPLIM.MPOVAL, MPYC.MPOVAL
  69.       POINTEUR CELL.IZAFM
  70. C-------------------------------------------------------
  71. -INC SMLREEL
  72.       POINTEUR MLRECP.MLREEL, MLRECV.MLREEL
  73. C-------------------------------------------------------
  74. C*********  Les Jacobians ******************************
  75. C-------------------------------------------------------
  76.       SEGMENT JACEL
  77.          REAL*8 JAC(3+NSP,3+NSP)
  78.       ENDSEGMENT
  79.       POINTEUR JTL.JACEL, WL.JACEL, JTT.JACEL
  80. C-------------------------------------------------------------
  81. C******* Les fractionines massiques **************************
  82. C-------------------------------------------------------------
  83.        SEGMENT FRAMAS
  84.          REAL*8 YET(NSP)
  85.       ENDSEGMENT
  86.       POINTEUR YC.FRAMAS, YF.FRAMAS
  87. C-------------------------------------------------------
  88. C**********  Les CP's and CV's   ***********************
  89. C-------------------------------------------------------
  90.       SEGMENT GCONST
  91.          REAL*8 GC(NSP)
  92.       ENDSEGMENT
  93.       POINTEUR CP.GCONST, CV.GCONST
  94. C-------------------------------------------------------------
  95. C**********  Segments for the vectors  ***********************
  96. C-------------------------------------------------------------
  97.        SEGMENT VECEL
  98.          REAL*8 VV(NSP)
  99.        ENDSEGMENT
  100.        POINTEUR DGDYC.VECEL
  101. C--------------------------------------------------
  102.       SEGINI JTL
  103.       SEGINI JTT
  104.       SEGINI WL
  105. C----------------------------------------------------
  106. C**** KRIPAD pour la correspondance global/local
  107. C----------------------------------------------------
  108.       CALL KRIPAD(MELEMC,MLEMC)
  109.       CALL KRIPAD(MELECB,MLEMCB)
  110.       CALL KRIPAD(MELEMF,MLEMF)
  111. C----------------------------------------------------
  112. C**** CHPOINTs de la table DOMAINE
  113. C----------------------------------------------------
  114.       CALL LICHT(INORM,MPNORM,TYPE,ICEL)
  115.       CALL LICHT(ICHPVO,MPVOL,TYPE,ICEL)
  116.       CALL LICHT(ICHPSU,MPSURF,TYPE,ICEL)
  117. C----------------------------------------------------
  118. C**** CHPOINTs des variables
  119. C----------------------------------------------------
  120.       CALL LICHT(IROC,MPRC,TYPE,ICEL)
  121.       CALL LICHT(IVITC,MPVC,TYPE,ICEL)
  122.       CALL LICHT(IPC,MPPC,TYPE,ICEL)
  123.       CALL LICHT(IYC,MPYC,TYPE,ICEL)
  124.       CALL LICHT(ICHLIM,MPLIM,TYPE,ICEL)
  125. C--------------------------------------------------------
  126. C**** Boucle sur le face pour le calcul des invariants de
  127. C     Riemann et du flux
  128. C--------------------------------------------------------
  129.       SEGACT MELEFC
  130.       NFAC=MELEFC.NUM(/2)
  131. C---------------------------------
  132. C**** Objet MATRIK
  133. C---------------------------------
  134.       NRIGE = 7
  135.       NMATRI = 1
  136.       NKID = 9
  137.       NKMT = 7
  138. C---------------------------------
  139.       SEGINI MATRIK
  140.       IJACO = MATRIK
  141.       MATRIK.IRIGEL(1,1) = MELRES
  142.       MATRIK.IRIGEL(2,1) = MELRES
  143. C---------------------------------
  144. C**** Matrice non symetrique
  145. C---------------------------------
  146.       MATRIK.IRIGEL(7,1) = 2
  147. C---------------------------------
  148.       NBME = (3+NSP)*(3+NSP)
  149.       NBSOUS = 1
  150.       SEGINI IMATRI
  151.       IF(IJAC.EQ.1)THEN
  152.          MLMOTS=ILIINC
  153.       ELSEIF(IJAC.EQ.2)THEN
  154.          MLMOTS=ILIINP
  155.       ENDIF
  156.       SEGACT MLMOTS
  157.       MATRIK.IRIGEL(4,1) = IMATRI
  158. C-------------------------------------------
  159.       DO 1 J=1,(NSP+3)
  160.         KV=(J-1)*(3+NSP)
  161.         IMATRI.LISPRI(KV+1) = MLMOTS.MOTS(1)
  162.         IMATRI.LISPRI(KV+2) = MLMOTS.MOTS(2)
  163.         IMATRI.LISPRI(KV+3) = MLMOTS.MOTS(3)
  164.         IMATRI.LISPRI(KV+4) = MLMOTS.MOTS(4)
  165.         DO 2 I=1,(NSP-1)
  166.           IMATRI.LISPRI(KV+4+I) = MLMOTS.MOTS(4+I)
  167.  2      CONTINUE
  168.  1    CONTINUE
  169. C-----------------------------------------------
  170.       SEGDES MLMOTS
  171.       MLMOTS=ILIINC
  172.       SEGACT MLMOTS
  173. C-----------------------------------------------
  174.       DO 3 J=1,(NSP+3)
  175.         KV=(J-1)*(3+NSP)
  176.         IMATRI.LISDUA(KV+1) = MLMOTS.MOTS(j)
  177.         IMATRI.LISDUA(KV+2) = MLMOTS.MOTS(j)
  178.         IMATRI.LISDUA(KV+3) = MLMOTS.MOTS(j)
  179.         IMATRI.LISDUA(KV+4) = MLMOTS.MOTS(j)
  180.         DO 4 I=1,(NSP-1)
  181.           IMATRI.LISDUA(KV+4+I) = MLMOTS.MOTS(j)
  182.  4     CONTINUE
  183.  3    CONTINUE
  184. C-----------------------------------------------
  185. C-----------------------------------------------
  186.       SEGDES MLMOTS
  187.       NBEL = NFAC
  188.       NBSOUS = 1
  189.       NP = 1
  190.       MP = 1
  191. C-----------------------------------------------------------
  192. C-----------------------------------------------------------
  193.       DO 5 I=1,NBME
  194.         SEGINI CELL
  195.         IMATRI.LIZAFM(1,I) = CELL
  196.  5    CONTINUE
  197. C---------------------------------
  198. C---------------------------------
  199. C**** Fin definition MATRIK
  200. C---------------------------------
  201.       DO IFAC=1,NFAC,1
  202.          NGF=MELEFC.NUM(1,IFAC)
  203.          NGC=MELEFC.NUM(2,IFAC)
  204.          NLF=MLEMF.LECT(NGF)
  205.          NLC=MLEMC.LECT(NGC)
  206.          NLCB=MLEMCB.LECT(NGF)
  207.          VOLU=MPVOL.VPOCHA(NLC,1)
  208.          SURF=MPSURF.VPOCHA(NLF,1)
  209. C        In CASTEM les normales sont sortantes
  210.          CNX=-1*MPNORM.VPOCHA(NLF,1)
  211.          CNY=-1*MPNORM.VPOCHA(NLF,2)
  212. C----------------------------------------------
  213.          SEGINI CP, CV
  214.          MLRECP = LRECP
  215.          MLRECV = LRECV
  216.          SEGACT MLRECP, MLRECV
  217.          DO 10 I=1,(NSP-1)
  218.            CP.GC(I)=MLRECP.PROG(I)
  219.            CV.GC(I)=MLRECV.PROG(I)
  220.  10      CONTINUE
  221.          CP.GC(NSP)=MLRECP.PROG(NSP)
  222.          CV.GC(NSP)=MLRECV.PROG(NSP)
  223. C---------------------------------
  224. C        Variables au centre
  225. C---------------------------------
  226.          RC=MPRC.VPOCHA(NLC,1)
  227.          PC=MPPC.VPOCHA(NLC,1)
  228.          UXC=MPVC.VPOCHA(NLC,1)
  229.          UYC=MPVC.VPOCHA(NLC,2)
  230.          SEGINI YC
  231.          SEGACT MPYC
  232.          DO 100 I=1,(NSP-1)
  233.            YC.YET(I)=MPYC.VPOCHA(NLC,I)
  234.  100     CONTINUE
  235. C---------------------------------
  236. C        Variables à la face
  237. C---------------------------------
  238.          RHOUF=MPLIM.VPOCHA(NLCB,1)
  239.          PSRF=MPLIM.VPOCHA(NLCB,2)
  240.          SEGINI YF
  241.          SEGACT MPLIM
  242.          DO 101 I=1,(NSP-1)
  243.            YF.YET(I)=MPLIM.VPOCHA(NLCB,2+I)
  244.  101     CONTINUE
  245. c-------------------------------------------------------------
  246. c  Computing GAMMA at the cell-center
  247. c-------------------------------------------------------------
  248.          top=0.0D0
  249.          bot=0.0D0
  250.          do 102  i=1,(nsp-1)
  251.           top=top+yc.yet(i)*(cp.gc(i)-cp.gc(nsp))
  252.           bot=bot+yc.yet(i)*(cv.gc(i)-cv.gc(nsp))
  253.  102     continue
  254.          top=cp.gc(nsp)+top
  255.          bot=cv.gc(nsp)+bot
  256.          GAMC=top/bot
  257. C-------------------------------------------------------------
  258.          SEGINI DGDYC
  259.          do 41 i=1,(nsp-1)
  260.            dgdyc.vv(i)=(cp.gc(i)-cp.gc(nsp)-
  261.      &       GAMC*(cv.gc(i)-cv.gc(nsp)))/bot
  262.  41    continue
  263. c-------------------------------------------------------------
  264. C------------------------------------------------
  265. C******* Densite, vitesse, pression sur le bord
  266. C------------------------------------------------
  267.          P=PC
  268.          RHO=P/PSRF
  269.          UN=RHOUF/RHO
  270. C------------------------------
  271. C******* Derivatives w.r.t. PC
  272. C------------------------------
  273.          DUNDPC=-1*RHOUF/(RHO*RHO*PSRF)
  274.          DUXDPC=CNX*DUNDPC
  275.          DUYDPC=CNY*DUNDPC
  276.          DECIDP=UN*DUNDPC
  277. C--------------------------------------------------------------
  278.          COEF=SURF/VOLU
  279. C-------------------------------------------------------------
  280.          JTL.JAC(1,1) = 0.0D0
  281.          JTL.JAC(1,2) = 0.0D0
  282.          JTL.JAC(1,3) = 0.0D0
  283.          JTL.JAC(1,4) = 0.0D0
  284.          DO 107 I=1,(NSP-1)
  285.            JTL.JAC(1,4+I) = 0.0D0
  286.  107     CONTINUE
  287. C----------------------------------------
  288.          JTL.JAC(2,1) = 0.0D0
  289.          JTL.JAC(2,2) = 0.0D0
  290.          JTL.JAC(2,3) = 0.0D0
  291.          JTL.JAC(2,4) = (RHOUF*DUXDPC+CNX)*COEF
  292.          DO 108 I=1,(NSP-1)
  293.            JTL.JAC(2,4+I) = 0.0D0
  294.  108     CONTINUE
  295. C----------------------------------------
  296.          JTL.JAC(3,1) = 0.0D0
  297.          JTL.JAC(3,2) = 0.0D0
  298.          JTL.JAC(3,3) = 0.0D0
  299.          JTL.JAC(3,4) = (RHOUF*DUYDPC+CNY)*COEF
  300.          DO 109 I=1,(NSP-1)
  301.            JTL.JAC(3,4+I) = 0.0D0
  302.  109     CONTINUE
  303. C----------------------------------------
  304.          JTL.JAC(4,1) = 0.0D0
  305.          JTL.JAC(4,2) = 0.0D0
  306.          JTL.JAC(4,3) = 0.0D0
  307.          JTL.JAC(4,4) = (RHOUF*DECIDP)*COEF
  308.          DO 110 I=1,(NSP-1)
  309.            JTL.JAC(4,4+I) = 0.0D0
  310.  110     CONTINUE
  311. C----------------------------------------
  312.          DO 111 I=1,(NSP-1)
  313.           JTL.JAC(4+I,1) = 0.0D0
  314.           JTL.JAC(4+I,2) = 0.0D0
  315.           JTL.JAC(4+I,3) = 0.0D0
  316.           JTL.JAC(4+I,4) = 0.0D0
  317.           DO 112 J=1,(NSP-1)
  318.            JTL.JAC(4+I,4+J) = 0.0D0
  319.  112      CONTINUE
  320.  111     CONTINUE
  321. C---------------------------------------------------------
  322. c  matrix wl(i,j) represents the derivative of the i-component
  323. c  of the vector of primitive variables of the left state with
  324. c  respect to the j-component of the vector of the conservative
  325. c  variables of the left state.
  326. c
  327. c  Here: (rho, ux, uy, p, Y_1,...,Y_(nsp-1)) -
  328. c                               vector of primitive variables;
  329. c  (rho, rho ux, rho uy, rho e, rho Y_1,..., rho Y_(nsp-1)) -
  330. c                            vector of conservative variables.
  331. c-------------------------------------------------------------
  332.        wl.jac(1,1)=1.0d0
  333.        wl.jac(1,2)=0.0d0
  334.        wl.jac(1,3)=0.0d0
  335.        wl.jac(1,4)=0.0d0
  336.        do 83 i=1,(nsp-1)
  337.          wl.jac(1,4+i)=0.0d0
  338.  83    continue
  339. c------------------------------
  340.        wl.jac(2,1)=-UXC/RC
  341.        wl.jac(2,2)=1.0d0/RC
  342.        wl.jac(2,3)=0.0d0
  343.        wl.jac(2,4)=0.0d0
  344.        do 84 i=1,(nsp-1)
  345.          wl.jac(2,4+i)=0.0d0
  346.  84    continue
  347. c------------------------------
  348.        wl.jac(3,1)=-UYC/RC
  349.        wl.jac(3,2)=0.0d0
  350.        wl.jac(3,3)=1.0d0/RC
  351.        wl.jac(3,4)=0.0d0
  352.        do 85 i=1,(nsp-1)
  353.          wl.jac(3,4+i)=0.0d0
  354.  85    continue
  355. c------------------------------
  356.        br1=0.0d0
  357.        do 86 i=1,(nsp-1)
  358.          br1=br1+dgdyc.vv(i)*yc.yet(i)
  359.  86    continue
  360.        br1=br1*PC/(RC*(GAMC-1.0D0))
  361.        wl.jac(4,1)=(GAMC-1.0D0)*(UXC*UXC+UYC*UYC)/2.0d0-br1
  362.        wl.jac(4,2)=-UXC*(GAMC-1.0D0)
  363.        wl.jac(4,3)=-UYC*(GAMC-1.0D0)
  364.        wl.jac(4,4)=(GAMC-1.0D0)
  365.        do 87 i=1,(nsp-1)
  366.          wl.jac(4,4+i)=dgdyc.vv(i)*PC/(RC*(GAMC-1.0D0))
  367.  87    continue
  368. c------------------------------
  369.        do 88 i=1,(nsp-1)
  370.         do 89 j=1,4
  371.           wl.jac(4+i,j)=0.0d0
  372.           if(j.eq.1) wl.jac(4+i,j)=-yc.yet(i)/RC
  373.  89     continue
  374. c------------
  375.         do 890 j=5,(4+nsp-1)
  376.           wl.jac(4+i,j)=0.0d0
  377.           if(4+i.eq.j)  then
  378.             wl.jac(4+i,j)=1.0d0/RC
  379.           endif
  380.  890    continue
  381.  88    continue
  382. c------------------------------------------------
  383. C------------------------------------------------
  384.        do 114 i=1,(3+nsp)
  385.         do 115 j=1,(3+nsp)
  386.          jtt.jac(i,j)=0.0d0
  387.          do 116 k=1,(3+nsp)
  388.           jtt.jac(i,j)=jtt.jac(i,j)+jtl.jac(i,k)*wl.jac(k,j)
  389.  116     continue
  390.  115    continue
  391.  114  continue
  392. C----------------------------------------------------------------
  393. C******* Jacobian with respect to conservative variables
  394. C----------------------------------------------------------------
  395.          IF(IJAC.EQ.1)THEN
  396.            DO 9 II = 1,(3+NSP)
  397.               DO 15 JJ = 1,(3+NSP)
  398.                 KV = (II-1)*(3+NSP)
  399. C----------------------------------
  400.                 CELL = IMATRI.LIZAFM(1,KV+JJ)
  401.                 CELL.AM(IFAC,1,1) = JTT.JAC(II,JJ)
  402.  15           CONTINUE
  403.  9         CONTINUE
  404.          ELSEIF(IJAC.EQ.2)THEN
  405.            DO 20 II = 1,(3+NSP)
  406.               DO 25 JJ = 1,(3+NSP)
  407.                 KV = (II-1)*(3+NSP)
  408. C----------------------------------
  409.                 CELL = IMATRI.LIZAFM(1,KV+JJ)
  410.                 CELL.AM(IFAC,1,1) = JTL.JAC(II,JJ)
  411.  25           CONTINUE
  412.  20         CONTINUE
  413.          ENDIF
  414. c--------------------------------------------------
  415.       ENDDO
  416. C
  417.       SEGDES MELEFC
  418. C
  419.       SEGSUP MLEMC
  420.       SEGSUP MLEMCB
  421.       SEGSUP MLEMF
  422. C
  423.       SEGDES MPNORM
  424.       SEGDES MPVOL
  425.       SEGDES MPSURF
  426.       SEGDES MPRC
  427.       SEGDES MPPC
  428.       SEGDES MPVC
  429.       SEGDES MPYC
  430.       SEGDES MPLIM
  431.       SEGDES YC
  432.       SEGDES YF
  433.       SEGDES CP
  434.       SEGDES CV
  435.       SEGDES JTL
  436.       SEGDES JTT
  437.       SEGDES WL
  438.       SEGDES  DGDYC
  439.       SEGDES MATRIK
  440.       DO 80 II=1,NBME
  441.        CELL = IMATRI.LIZAFM(1,II)
  442.        SEGDES CELL
  443.  80   CONTINUE
  444.       SEGDES IMATRI
  445. C---------------------------------------------
  446.  9999 CONTINUE
  447.       RETURN
  448.       END
  449.  
  450.  
  451.  
  452.  
  453.  
  454.  
  455.  
  456.  
  457.  
  458.  
  459.  
  460.  
  461.  
  462.  
  463.  

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