Télécharger cli262.eso

Retour à la liste

Numérotation des lignes :

  1. C CLI262 SOURCE PV 16/11/17 21:58:45 9180
  2. SUBROUTINE CLI262(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: 'INSU' -- 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,CTX,CTY
  52. & ,PSRF,COEF
  53. & ,BR1,BOT,TOP,DUNDUY,GAMF,UNC
  54. & ,COEF5,DPSRUN,DUXDUN,DRDUN,DUYDUN,DUNDUX,DPDUN
  55. & ,GM1,USGM1,RF,UNF,UTF,UXF,UYF,HTF,PF,SF,ECIN
  56. & ,TVECT(2),NVECT(2),WVEC_L(4),WVEC_R(4)
  57. CHARACTER*(8) TYPE
  58. -INC CCOPTIO
  59. -INC SMLMOTS
  60. -INC SMELEME
  61. POINTEUR MELEFC.MELEME
  62. -INC SMLENTI
  63. POINTEUR MLEMC.MLENTI, MLEMCB.MLENTI,MLEMF.MLENTI
  64. -INC SMCHPOI
  65. POINTEUR MPNORM.MPOVAL, MPVOL.MPOVAL, MPSURF.MPOVAL, MPRC.MPOVAL,
  66. & MPVC.MPOVAL, MPPC.MPOVAL, MPLIM.MPOVAL, MPYC.MPOVAL
  67. POINTEUR CELL.IZAFM
  68. C-------------------------------------------------------
  69. -INC SMLREEL
  70. POINTEUR MLRECP.MLREEL, MLRECV.MLREEL
  71. C-------------------------------------------------------
  72. C********* Les Jacobians ******************************
  73. C-------------------------------------------------------
  74. SEGMENT JACEL
  75. REAL*8 JAC(3+NSP,3+NSP)
  76. ENDSEGMENT
  77. POINTEUR JLL.JACEL, WL.JACEL, JRR.JACEL, JTT.JACEL,
  78. & JFL.JACEL, JFR.JACEL, JR.JACEL
  79. C-------------------------------------------------------------
  80. C******* Les fractionines massiques **************************
  81. C-------------------------------------------------------------
  82. SEGMENT FRAMAS
  83. REAL*8 YET(NSP)
  84. ENDSEGMENT
  85. POINTEUR YC.FRAMAS, YF.FRAMAS
  86. C-------------------------------------------------------
  87. C********** Les CP's and CV's ***********************
  88. C-------------------------------------------------------
  89. SEGMENT GCONST
  90. REAL*8 GC(NSP)
  91. ENDSEGMENT
  92. POINTEUR CP.GCONST, CV.GCONST
  93. C-------------------------------------------------------------
  94. C********** Segments for the vectors ***********************
  95. C-------------------------------------------------------------
  96. SEGMENT VECEL
  97. REAL*8 VV(NSP)
  98. ENDSEGMENT
  99. POINTEUR DGDYC.VECEL
  100. C--------------------------------------------------
  101. SEGINI JTT
  102. SEGINI JLL
  103. SEGINI JRR
  104. SEGINI JFL
  105. SEGINI JFR
  106. SEGINI WL
  107. C----------------------------------------------------
  108. C**** KRIPAD pour la correspondance global/local
  109. C----------------------------------------------------
  110. CALL KRIPAD(MELEMC,MLEMC)
  111. CALL KRIPAD(MELECB,MLEMCB)
  112. CALL KRIPAD(MELEMF,MLEMF)
  113. C----------------------------------------------------
  114. C**** CHPOINTs de la table DOMAINE
  115. C----------------------------------------------------
  116. CALL LICHT(INORM,MPNORM,TYPE,ICEL)
  117. CALL LICHT(ICHPVO,MPVOL,TYPE,ICEL)
  118. CALL LICHT(ICHPSU,MPSURF,TYPE,ICEL)
  119. C----------------------------------------------------
  120. C**** CHPOINTs des variables
  121. C----------------------------------------------------
  122. CALL LICHT(IROC,MPRC,TYPE,ICEL)
  123. CALL LICHT(IVITC,MPVC,TYPE,ICEL)
  124. CALL LICHT(IPC,MPPC,TYPE,ICEL)
  125. CALL LICHT(IYC,MPYC,TYPE,ICEL)
  126. CALL LICHT(ICHLIM,MPLIM,TYPE,ICEL)
  127. C--------------------------------------------------------
  128. C**** Boucle sur le face pour le calcul des invariants de
  129. C Riemann et du flux
  130. C--------------------------------------------------------
  131. SEGACT MELEFC
  132. NFAC=MELEFC.NUM(/2)
  133. C---------------------------------
  134. C**** Objet MATRIK
  135. C---------------------------------
  136. NRIGE = 7
  137. NMATRI = 1
  138. NKID = 9
  139. NKMT = 7
  140. C---------------------------------
  141. SEGINI MATRIK
  142. IJACO = MATRIK
  143. MATRIK.IRIGEL(1,1) = MELRES
  144. MATRIK.IRIGEL(2,1) = MELRES
  145. C---------------------------------
  146. C**** Matrice non symetrique
  147. C---------------------------------
  148. MATRIK.IRIGEL(7,1) = 2
  149. C---------------------------------
  150. NBME = (3+NSP)*(3+NSP)
  151. NBSOUS = 1
  152. SEGINI IMATRI
  153. IF(IJAC.EQ.1)THEN
  154. MLMOTS=ILIINC
  155. ELSEIF(IJAC.EQ.2)THEN
  156. MLMOTS=ILIINP
  157. ENDIF
  158. SEGACT MLMOTS
  159. MATRIK.IRIGEL(4,1) = IMATRI
  160. C-------------------------------------------
  161. DO 1 J=1,(NSP+3)
  162. KV=(J-1)*(3+NSP)
  163. IMATRI.LISPRI(KV+1) = MLMOTS.MOTS(1)
  164. IMATRI.LISPRI(KV+2) = MLMOTS.MOTS(2)
  165. IMATRI.LISPRI(KV+3) = MLMOTS.MOTS(3)
  166. IMATRI.LISPRI(KV+4) = MLMOTS.MOTS(4)
  167. DO 2 I=1,(NSP-1)
  168. IMATRI.LISPRI(KV+4+I) = MLMOTS.MOTS(4+I)
  169. 2 CONTINUE
  170. 1 CONTINUE
  171. C-----------------------------------------------
  172. SEGDES MLMOTS
  173. MLMOTS=ILIINC
  174. SEGACT MLMOTS
  175. C-----------------------------------------------
  176. DO 3 J=1,(NSP+3)
  177. KV=(J-1)*(3+NSP)
  178. IMATRI.LISDUA(KV+1) = MLMOTS.MOTS(j)
  179. IMATRI.LISDUA(KV+2) = MLMOTS.MOTS(j)
  180. IMATRI.LISDUA(KV+3) = MLMOTS.MOTS(j)
  181. IMATRI.LISDUA(KV+4) = MLMOTS.MOTS(j)
  182. DO 4 I=1,(NSP-1)
  183. IMATRI.LISDUA(KV+4+I) = MLMOTS.MOTS(j)
  184. 4 CONTINUE
  185. 3 CONTINUE
  186. C-----------------------------------------------
  187. C-----------------------------------------------
  188. SEGDES MLMOTS
  189. NBEL = NFAC
  190. NBSOUS = 1
  191. NP = 1
  192. MP = 1
  193. C-----------------------------------------------------------
  194. C-----------------------------------------------------------
  195. DO 5 I=1,NBME
  196. SEGINI CELL
  197. IMATRI.LIZAFM(1,I) = CELL
  198. 5 CONTINUE
  199. C---------------------------------
  200. C---------------------------------
  201. C**** Fin definition MATRIK
  202. C---------------------------------
  203. DO IFAC=1,NFAC,1
  204. NGF=MELEFC.NUM(1,IFAC)
  205. NGC=MELEFC.NUM(2,IFAC)
  206. NLF=MLEMF.LECT(NGF)
  207. NLC=MLEMC.LECT(NGC)
  208. NLCB=MLEMCB.LECT(NGF)
  209. VOLU=MPVOL.VPOCHA(NLC,1)
  210. SURF=MPSURF.VPOCHA(NLF,1)
  211. C In CASTEM les normales sont sortantes
  212. CNX=-1*MPNORM.VPOCHA(NLF,1)
  213. CNY=-1*MPNORM.VPOCHA(NLF,2)
  214. CTX=-1.0D0*CNY
  215. CTY=CNX
  216. C----------------------------------------------
  217. SEGINI CP, CV
  218. MLRECP = LRECP
  219. MLRECV = LRECV
  220. SEGACT MLRECP, MLRECV
  221. DO 10 I=1,(NSP-1)
  222. CP.GC(I)=MLRECP.PROG(I)
  223. CV.GC(I)=MLRECV.PROG(I)
  224. 10 CONTINUE
  225. CP.GC(NSP)=MLRECP.PROG(NSP)
  226. CV.GC(NSP)=MLRECV.PROG(NSP)
  227. C---------------------------------
  228. C Variables au centre
  229. C---------------------------------
  230. RC=MPRC.VPOCHA(NLC,1)
  231. PC=MPPC.VPOCHA(NLC,1)
  232. UXC=MPVC.VPOCHA(NLC,1)
  233. UYC=MPVC.VPOCHA(NLC,2)
  234. SEGINI YC
  235. SEGACT MPYC
  236. DO 100 I=1,(NSP-1)
  237. YC.YET(I)=MPYC.VPOCHA(NLC,I)
  238. 100 CONTINUE
  239. UNC=(UXC*CNX)+(UYC*CNY)
  240. C---------------------------------
  241. C Variables à la face
  242. C---------------------------------
  243. SEGACT MPLIM
  244. HTF=MPLIM.VPOCHA(NLCB,1)
  245. SF=MPLIM.VPOCHA(NLCB,2)
  246. SEGINI YF
  247. DO 101 I=1,(NSP-1)
  248. YF.YET(I)=MPLIM.VPOCHA(NLCB,2+I)
  249. 101 CONTINUE
  250. UNF=UNC
  251. UTF=0.0D0
  252. UXF=UNF*CNX+UTF*CTX
  253. UYF=UNF*CNY+UTF*CTY
  254. c-------------------------------------------------------------
  255. c Computing GAMMA at the face
  256. c-------------------------------------------------------------
  257. top=0.0D0
  258. bot=0.0D0
  259. do 102 i=1,(nsp-1)
  260. top=top+yf.yet(i)*(cp.gc(i)-cp.gc(nsp))
  261. bot=bot+yf.yet(i)*(cv.gc(i)-cv.gc(nsp))
  262. 102 continue
  263. top=cp.gc(nsp)+top
  264. bot=cv.gc(nsp)+bot
  265. GAMF=top/bot
  266. GM1=GAMF-1.0D0
  267. USGM1=1.0D0/GM1
  268. c-------------------------------------------------------------
  269. c Computing GAMMA at the cell-centre
  270. c-------------------------------------------------------------
  271. top=0.0D0
  272. bot=0.0D0
  273. do 103 i=1,(nsp-1)
  274. top=top+yc.yet(i)*(cp.gc(i)-cp.gc(nsp))
  275. bot=bot+yc.yet(i)*(cv.gc(i)-cv.gc(nsp))
  276. 103 continue
  277. top=cp.gc(nsp)+top
  278. bot=cv.gc(nsp)+bot
  279. GAMC=top/bot
  280. C-------------------------------------------------------------
  281. SEGINI DGDYC
  282. do 41 i=1,(nsp-1)
  283. dgdyc.vv(i)=(cp.gc(i)-cp.gc(nsp)-
  284. & GAMC*(cv.gc(i)-cv.gc(nsp)))/bot
  285. 41 continue
  286. c-------------------------------------------------------------
  287. C------------------------------------------------
  288. C******* Densite, vitesse, pression sur le bord
  289. C------------------------------------------------
  290. ECIN=0.5D0*((UXF*UXF)+(UYF*UYF))
  291. PSRF=(GM1/GAMF)*(HTF-ECIN)
  292. RF=PSRF/SF
  293. RF=RF**(1.0D0/GM1)
  294. PF=SF*(RF**GAMF)
  295. C------------------------------
  296. C******* Derivatives w.r.t. PC
  297. C------------------------------
  298. DPSRUN=-1*(GM1/GAMF)*UNC
  299. DRDUN=USGM1*RF/PSRF*DPSRUN
  300. DPDUN=GAMF*PSRF*DRDUN
  301. DUXDUN=CNX
  302. DUYDUN=CNY
  303. C-------------------------------
  304. DUNDUX=CNX
  305. DUNDUY=CNY
  306. C--------------------------------------------------------------
  307. C------------------------------
  308. C******* Derivatives
  309. C------------------------------
  310. wvec_l(1)=RF
  311. wvec_l(2)=UXF
  312. wvec_l(3)=UYF
  313. wvec_l(4)=PF
  314. C--------------------------
  315. wvec_r(1)=RC
  316. wvec_r(2)=UXC
  317. wvec_r(3)=UYC
  318. wvec_r(4)=PC
  319. C--------------------------
  320. nvect(1)=CNX
  321. nvect(2)=CNY
  322. tvect(1)=CTX
  323. tvect(2)=CTY
  324. call copms2(nsp,jfl,jfr,wvec_l,wvec_r,nvect,tvect,
  325. & yf,mpyc,lrecp,lrecv,nlc)
  326. C--------------------------------------------------------------
  327. COEF=SURF/VOLU
  328. JLL=JFL
  329. JRR=JFR
  330. SEGACT JLL
  331. SEGACT JRR
  332. SEGINI JR
  333. C-------------------------------------------------------------
  334. C Taking in to account the dependance of the variables
  335. C at the face on the variables at the centre (UNC)
  336. C-------------------------------------------------------------
  337. DO 105 I=1,(3+NSP)
  338. DO 1050 J=1,(3+NSP)
  339. IF(J .EQ. 2) THEN
  340. COEF5=(JLL.JAC(I,1)*DRDUN)+(JLL.JAC(I,2)*DUXDUN)+
  341. & (JLL.JAC(I,3)*DUYDUN)+(JLL.JAC(I,4)*DPDUN)
  342. JR.JAC(I,J)=(COEF5*DUNDUX+JRR.JAC(I,J))*COEF
  343. ELSEIF(J .EQ. 3) THEN
  344. COEF5=(JLL.JAC(I,1)*DRDUN)+(JLL.JAC(I,2)*DUXDUN)+
  345. & (JLL.JAC(I,3)*DUYDUN)+(JLL.JAC(I,4)*DPDUN)
  346. JR.JAC(I,J)=(COEF5*DUNDUY+JRR.JAC(I,J))*COEF
  347. ELSE
  348. JR.JAC(I,J)=JRR.JAC(I,J)*COEF
  349. ENDIF
  350. 1050 CONTINUE
  351. 105 CONTINUE
  352. C-------------------------------------------------------------
  353. c matrix wl(i,j) represents the derivative of the i-component
  354. c of the vector of primitive variables of the left state with
  355. c respect to the j-component of the vector of the conservative
  356. c variables of the left state.
  357. c
  358. c Here: (rho, ux, uy, p, Y_1,...,Y_(nsp-1)) -
  359. c vector of primitive variables;
  360. c (rho, rho ux, rho uy, rho e, rho Y_1,..., rho Y_(nsp-1)) -
  361. c vector of conservative variables.
  362. c-------------------------------------------------------------
  363. wl.jac(1,1)=1.0d0
  364. wl.jac(1,2)=0.0d0
  365. wl.jac(1,3)=0.0d0
  366. wl.jac(1,4)=0.0d0
  367. do 83 i=1,(nsp-1)
  368. wl.jac(1,4+i)=0.0d0
  369. 83 continue
  370. c------------------------------
  371. wl.jac(2,1)=-UXC/RC
  372. wl.jac(2,2)=1.0d0/RC
  373. wl.jac(2,3)=0.0d0
  374. wl.jac(2,4)=0.0d0
  375. do 84 i=1,(nsp-1)
  376. wl.jac(2,4+i)=0.0d0
  377. 84 continue
  378. c------------------------------
  379. wl.jac(3,1)=-UYC/RC
  380. wl.jac(3,2)=0.0d0
  381. wl.jac(3,3)=1.0d0/RC
  382. wl.jac(3,4)=0.0d0
  383. do 85 i=1,(nsp-1)
  384. wl.jac(3,4+i)=0.0d0
  385. 85 continue
  386. c------------------------------
  387. br1=0.0d0
  388. do 86 i=1,(nsp-1)
  389. br1=br1+dgdyc.vv(i)*yc.yet(i)
  390. 86 continue
  391. br1=br1*PC/(RC*(GAMC-1.0D0))
  392. wl.jac(4,1)=(GAMC-1.0D0)*(UXC*UXC+UYC*UYC)/2.0d0-br1
  393. wl.jac(4,2)=-UXC*(GAMC-1.0D0)
  394. wl.jac(4,3)=-UYC*(GAMC-1.0D0)
  395. wl.jac(4,4)=(GAMC-1.0D0)
  396. do 87 i=1,(nsp-1)
  397. wl.jac(4,4+i)=dgdyc.vv(i)*PC/(RC*(GAMC-1.0D0))
  398. 87 continue
  399. c------------------------------
  400. do 88 i=1,(nsp-1)
  401. do 89 j=1,4
  402. wl.jac(4+i,j)=0.0d0
  403. if(j.eq.1) wl.jac(4+i,j)=-yc.yet(i)/RC
  404. 89 continue
  405. c------------
  406. do 890 j=5,(4+nsp-1)
  407. wl.jac(4+i,j)=0.0d0
  408. if(4+i.eq.j) then
  409. wl.jac(4+i,j)=1.0d0/RC
  410. endif
  411. 890 continue
  412. 88 continue
  413. c------------------------------------------------
  414. C------------------------------------------------
  415. do 114 i=1,(3+nsp)
  416. do 115 j=1,(3+nsp)
  417. jtt.jac(i,j)=0.0d0
  418. do 116 k=1,(3+nsp)
  419. jtt.jac(i,j)=jtt.jac(i,j)+
  420. & jr.jac(i,k)*wl.jac(k,j)/COEF
  421. 116 continue
  422. 115 continue
  423. 114 continue
  424. C----------------------------------------------------------------
  425. C******* Jacobian with respect to conservative variables
  426. C----------------------------------------------------------------
  427. IF(IJAC.EQ.1)THEN
  428. DO 9 II = 1,(3+NSP)
  429. DO 15 JJ = 1,(3+NSP)
  430. KV = (II-1)*(3+NSP)
  431. C----------------------------------
  432. CELL = IMATRI.LIZAFM(1,KV+JJ)
  433. CELL.AM(IFAC,1,1) = JTT.JAC(II,JJ)*COEF
  434. 15 CONTINUE
  435. 9 CONTINUE
  436. ELSEIF(IJAC.EQ.2)THEN
  437. DO 20 II = 1,(3+NSP)
  438. DO 25 JJ = 1,(3+NSP)
  439. KV = (II-1)*(3+NSP)
  440. C----------------------------------
  441. CELL = IMATRI.LIZAFM(1,KV+JJ)
  442. CELL.AM(IFAC,1,1) = JR.JAC(II,JJ)
  443. 25 CONTINUE
  444. 20 CONTINUE
  445. ENDIF
  446. c--------------------------------------------------
  447. ENDDO
  448. C
  449. SEGDES MELEFC
  450. C
  451. SEGDES MLEMC
  452. SEGDES MLEMCB
  453. SEGDES MLEMF
  454. C
  455. SEGDES MPNORM
  456. SEGDES MPVOL
  457. SEGDES MPSURF
  458. SEGDES MPRC
  459. SEGDES MPPC
  460. SEGDES MPVC
  461. SEGDES MPYC
  462. SEGDES MPLIM
  463. SEGDES YC
  464. SEGDES YF
  465. SEGDES CP
  466. SEGDES CV
  467. SEGDES JLL
  468. SEGDES JRR
  469. SEGDES JTT
  470. SEGDES JFL
  471. SEGDES JFR
  472. SEGDES JR
  473. SEGDES WL
  474. SEGDES DGDYC
  475. SEGDES MATRIK
  476. DO 80 II=1,NBME
  477. CELL = IMATRI.LIZAFM(1,II)
  478. SEGDES CELL
  479. 80 CONTINUE
  480. SEGDES IMATRI
  481. C---------------------------------------------
  482. 9999 CONTINUE
  483. RETURN
  484. END
  485.  
  486.  
  487.  
  488.  
  489.  
  490.  
  491.  
  492.  
  493.  
  494.  
  495.  
  496.  
  497.  
  498.  
  499.  

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