Télécharger konms3.eso

Retour à la liste

Numérotation des lignes :

  1. C KONMS3 SOURCE PV 16/11/17 22:00:11 9180
  2. SUBROUTINE KONMS3(NSP,ILIINC,IRN,IUN,IPN,IYN,MLRECP,
  3. & MLRECV,INORM,ICHPVO,ICHPSU,MELEMC,MELEFE,MELLIM,IMAT)
  4. C
  5. C************************************************************************
  6. C
  7. C PROJET : CASTEM 2000
  8. C
  9. C NOM : KONMS3 ("convection for multispecies")
  10. C
  11. C DESCRIPTION : Voir KON19 (appele par KON19)
  12. C Calcul du jacobien du résidu pour la méthode
  13. C AUSM+
  14. C
  15. C Cas trois dimensions, gaz "calorically perfect"
  16. C MULTISPECIES!!!!!!!!
  17. C
  18. C LANGAGE : FORTRAN 77 + ESOPE 2000 (avec estensions CISI)
  19. C
  20. C AUTEUR : S. KUDRIAKOV, DM2S/SFME/LTMF
  21. C
  22. C************************************************************************
  23. C
  24. C
  25. C APPELES (Outils
  26. C CASTEM) : KRIPAD, LICHT, ERREUR
  27. C
  28. C APPELES (Calcul) : CMS3D, CWMS3D
  29. C
  30. C************************************************************************
  31. C
  32. C ENTREES
  33. C
  34. C NSP : number of species (total) ;
  35. C
  36. C ILINC : liste des inconnues (pointeur d'un objet de type LISTMOTS)
  37. C
  38. C 1) Pointeurs des CHPOINT
  39. C
  40. C IRN : CHPOINT CENTRE contenant la masse volumique ;
  41. C
  42. C IUN : CHPOINT CENTRE contenant la vitesse ;
  43. C
  44. C IPN : CHPOINT CENTRE contenant la pression ;
  45. C
  46. C IYN : CHPOINT CENTRE contenant les fractions massiques ;
  47. C
  48. C INORM : CHPOINT FACE contenant les normales aux faces ;
  49. C
  50. C ICHPVO : CHPOINT VOLUME contenant le volume
  51. C
  52. C ICHPSU : CHPOINT FACE contenant la surface des faces
  53. C
  54. C
  55. C 2) Pointeurs des LIST REELS
  56. C
  57. C MLRECP : LIST REELS contenant les CP's des gases differents ;
  58. C
  59. C MLRECV : LIST REELS contenant les CV's des gases differents ;
  60. C
  61. C 3) Pointeurs de MELEME de la table DOMAINE
  62. C
  63. C MELEMC : MELEME 'CENTRE' du SPG des CENTRES
  64. C
  65. C MELEFE : MELEME 'FACEL' du connectivité Faces -> Elts
  66. C
  67. C MELLIM : MELEME SPG des conditions aux bords
  68. C
  69. C SORTIES
  70. C
  71. C IMAT : pointeur de la MATRIK du jacobien du residu
  72. C
  73. C************************************************************************
  74. C
  75. C HISTORIQUE (Anomalies et modifications éventuelles)
  76. C
  77. C HISTORIQUE :
  78. C
  79. C************************************************************************
  80. C
  81. C
  82. C N.B.: On suppose qu'on a déjà controllé RO, P > 0
  83. C GAMMA \in (1,3)
  84. C Si non il faut le faire!!!
  85. C
  86. C************************************************************************
  87. C
  88. IMPLICIT INTEGER(I-N)
  89. IMPLICIT REAL*8(A-H,O-Z)
  90.  
  91. C---------------------------------------------------
  92. INTEGER NSP,II,JJ,KV,I,J,JLL,JRR
  93. INTEGER ILIINC, IRN,IUN,IPN,IYN,INORM,ICHPVO,ICHPSU
  94. & , IMAT, IGEOMC, IGEOMF
  95. & , NFAC, NBSOUS, NBREF, NBELEM, NBNN, NRIGE, NMATRI, NKID
  96. & , NKMT, NBME, NBEL, MP, NP
  97. & , IFAC, NGCF, NLCF, NGCG, NGCD, NLCG, NLCD, NLFL
  98. REAL*8 ROG, PG, UXG, UYG, UZG, VOLG
  99. & , ROD, PD, UXD, UYD, UZD, VOLD
  100. & , SURF, FUNCEL
  101. REAL*8 WVEC_L(5), WVEC_R(5), NVECT(3), TVECT1(3), TVECT2(3)
  102. REAL*8 C11,C12,C13,C21,C22,C23,C31,C32,C33,DET
  103. REAL*8 ZC11,ZC12,ZC13,ZC21,ZC22,ZC23,ZC31,ZC32,ZC33
  104. CHARACTER*8 TYPE
  105. C
  106. C**** LES INCLUDES
  107. C
  108.  
  109. -INC PPARAM
  110. -INC CCOPTIO
  111. -INC SMCHPOI
  112. -INC SMELEME
  113. -INC SMLMOTS
  114. -INC SMLENTI
  115. POINTEUR MPRN.MPOVAL, MPUN.MPOVAL, MPPN.MPOVAL, MPYN.MPOVAL,
  116. & MPNORM.MPOVAL, MPVOLU.MPOVAL, MPOVSU.MPOVAL
  117. POINTEUR MELEMC.MELEME, MELEMF.MELEME, MELEFE.MELEME,
  118. & MELEDU.MELEME, MELLIM.MELEME
  119. POINTEUR MLENTC.MLENTI, MLENTF.MLENTI, MLELIM.MLENTI
  120. POINTEUR CELL.IZAFM
  121. POINTEUR MLMINC.MLMOTS
  122. C----------------------------------------------------
  123. -INC SMLREEL
  124. POINTEUR MLRECP.MLREEL, MLRECV.MLREEL
  125. C-------------------------------------------------------
  126. C********* Les Jacobians ******************************
  127. C-------------------------------------------------------
  128. SEGMENT JACEL
  129. REAL*8 JAC(4+NSP,4+NSP)
  130. ENDSEGMENT
  131. POINTEUR JTL.JACEL, JTR.JACEL, JTT.JACEL
  132. C--------------------------------------------------------
  133. C KRIPAD pour la correspondance global/local des centres
  134. C--------------------------------------------------------
  135. CALL KRIPAD(MELLIM,MLELIM)
  136. CALL KRIPAD(MELEMC,MLENTC)
  137. C------------------
  138. SEGACT MELEMC
  139. C------------------
  140. SEGACT MELEFE
  141. C------------------
  142. CALL LICHT(ICHPSU,MPOVSU,TYPE,IGEOMF)
  143. CALL LICHT(INORM,MPNORM,TYPE,IGEOMF)
  144. CALL LICHT(ICHPVO,MPVOLU,TYPE,IGEOMC)
  145. C----------------------------------------------
  146. MELEMF = IGEOMF
  147. CALL KRIPAD(MELEMF,MLENTF)
  148. SEGACT MELEMF
  149. C-----------------------------------------------
  150. CALL LICHT(IRN,MPRN,TYPE,IGEOMC)
  151. CALL LICHT(IPN,MPPN,TYPE,IGEOMC)
  152. CALL LICHT(IUN,MPUN,TYPE,IGEOMC)
  153. CALL LICHT(IYN,MPYN,TYPE,IGEOMC)
  154. C-----------------------------------------------
  155. NFAC = MELEFE.NUM(/2)
  156. C-----------------------------------------------
  157. C**** Maillage des inconnues primales
  158. C-----------------------------------------------
  159. NBSOUS = 0
  160. NBREF = 0
  161. NBELEM = NFAC
  162. NBNN = 2
  163. C-----------------------------------------------
  164. SEGINI MELEDU
  165. C----------------------
  166. C** MELEDU = 'SEG2' **
  167. C----------------------
  168. MELEDU.ITYPEL = 2
  169. C--------------
  170. NRIGE = 7
  171. NMATRI = 1
  172. NKID = 9
  173. NKMT = 7
  174. C--------------
  175. SEGINI MATRIK
  176. IMAT = MATRIK
  177. MATRIK.IRIGEL(1,1) = MELEDU
  178. MATRIK.IRIGEL(2,1) = MELEDU
  179. C-------------------------
  180. C Matrice non symetrique
  181. C-------------------------
  182. MATRIK.IRIGEL(7,1) = 2
  183. C-------------------------
  184. NBME = (4+NSP)*(4+NSP)
  185. NBSOUS = 1
  186. SEGINI IMATRI
  187. MLMINC = ILIINC
  188. SEGACT MLMINC
  189. MATRIK.IRIGEL(4,1) = IMATRI
  190. C-----------------------------------------------
  191. DO 1 J=1,(NSP+4)
  192. KV=(J-1)*(4+NSP)
  193. IMATRI.LISPRI(KV+1) = MLMINC.MOTS(1)
  194. IMATRI.LISPRI(KV+2) = MLMINC.MOTS(2)
  195. IMATRI.LISPRI(KV+3) = MLMINC.MOTS(3)
  196. IMATRI.LISPRI(KV+4) = MLMINC.MOTS(4)
  197. IMATRI.LISPRI(KV+5) = MLMINC.MOTS(5)
  198. DO 2 I=1,(NSP-1)
  199. IMATRI.LISPRI(KV+5+I) = MLMINC.MOTS(5+I)
  200. 2 CONTINUE
  201. 1 CONTINUE
  202. C-----------------------------------------------
  203. DO 3 J=1,(NSP+4)
  204. KV=(J-1)*(4+NSP)
  205. IMATRI.LISDUA(KV+1) = MLMINC.MOTS(j)
  206. IMATRI.LISDUA(KV+2) = MLMINC.MOTS(j)
  207. IMATRI.LISDUA(KV+3) = MLMINC.MOTS(j)
  208. IMATRI.LISDUA(KV+4) = MLMINC.MOTS(j)
  209. IMATRI.LISDUA(KV+5) = MLMINC.MOTS(j)
  210. DO 4 I=1,(NSP-1)
  211. IMATRI.LISDUA(KV+5+I) = MLMINC.MOTS(j)
  212. 4 CONTINUE
  213. 3 CONTINUE
  214. C-----------------------------------------------
  215. C-----------------------------------------------
  216. NBEL = NBELEM
  217. NBSOUS = 1
  218. NP = 2
  219. MP = 2
  220. C-----------------------------------------------------------
  221. C-----------------------------------------------------------
  222. DO 5 I=1,NBME
  223. SEGINI CELL
  224. IMATRI.LIZAFM(1,I) = CELL
  225. 5 CONTINUE
  226. C---------------------------------
  227. C---------------------------------
  228. DO IFAC = 1, NFAC, 1
  229. NGCF = MELEFE.NUM(2,IFAC)
  230. NLCF = MLENTF.LECT(NGCF)
  231. IF(NLCF .NE. IFAC)THEN
  232. WRITE(IOIMP,*) 'Il ne faut pas jouer avec la table domaine'
  233. CALL ERREUR(5)
  234. GOTO 9999
  235. ENDIF
  236. NLFL = MLELIM.LECT(NGCF)
  237. NGCG = MELEFE.NUM(1,IFAC)
  238. NGCD = MELEFE.NUM(3,IFAC)
  239. IF(NLFL .NE. 0)THEN
  240. C--------------------------------------------------------
  241. C The point belongs on BC -> No contribution to jacobian!
  242. C--------------------------------------------------------
  243. MELEDU.NUM(1,IFAC) = NGCG
  244. MELEDU.NUM(2,IFAC) = NGCD
  245. ELSEIF(NGCG .NE. NGCD)THEN
  246. C-----------------------------------
  247. C********** Les MELEMEs
  248. C-----------------------------------
  249. MELEDU.NUM(1,IFAC) = NGCG
  250. MELEDU.NUM(2,IFAC) = NGCD
  251. C-----------------------------------
  252. C********** Les etats G et D
  253. C-----------------------------------
  254. NLCG = MLENTC.LECT(NGCG)
  255. NLCD = MLENTC.LECT(NGCD)
  256. C-----------------
  257. ROG = MPRN.VPOCHA(NLCG,1)
  258. PG = MPPN.VPOCHA(NLCG,1)
  259. UXG = MPUN.VPOCHA(NLCG,1)
  260. UYG = MPUN.VPOCHA(NLCG,2)
  261. UZG = MPUN.VPOCHA(NLCG,3)
  262. VOLG = MPVOLU.VPOCHA(NLCG,1)
  263. C-------------------------------------------------
  264. WVEC_L(1)=ROG
  265. WVEC_L(2)=UXG
  266. WVEC_L(3)=UYG
  267. WVEC_L(4)=UZG
  268. WVEC_L(5)=PG
  269. C-------------------------------------------------
  270. ROD = MPRN.VPOCHA(NLCD,1)
  271. PD = MPPN.VPOCHA(NLCD,1)
  272. UXD = MPUN.VPOCHA(NLCD,1)
  273. UYD = MPUN.VPOCHA(NLCD,2)
  274. UZD = MPUN.VPOCHA(NLCD,3)
  275. VOLD = MPVOLU.VPOCHA(NLCD,1)
  276. C------------------------------------------------
  277. WVEC_R(1)=ROD
  278. WVEC_R(2)=UXD
  279. WVEC_R(3)=UYD
  280. WVEC_R(4)=UZD
  281. WVEC_R(5)=PD
  282. C------------------------------------------------
  283. C La normale G->D
  284. C La tangente
  285. C------------------------------------------------
  286. SURF = MPOVSU.VPOCHA(NLCF,1)
  287. NVECT(1) = MPNORM.VPOCHA(NLCF,7)
  288. NVECT(2) = MPNORM.VPOCHA(NLCF,8)
  289. NVECT(3) = MPNORM.VPOCHA(NLCF,9)
  290. c-----------------------------------------------
  291. TVECT1(1) = MPNORM.VPOCHA(NLCF,1)
  292. TVECT1(2) = MPNORM.VPOCHA(NLCF,2)
  293. TVECT1(3) = MPNORM.VPOCHA(NLCF,3)
  294. c----------------------------------------------
  295. TVECT2(1) = MPNORM.VPOCHA(NLCF,4)
  296. TVECT2(2) = MPNORM.VPOCHA(NLCF,5)
  297. TVECT2(3) = MPNORM.VPOCHA(NLCF,6)
  298. C-------------------------------------------------
  299. CALL cms3d(NSP,JLL,JRR,WVEC_L,WVEC_R,NVECT,TVECT1,TVECT2,
  300. & MPYN,MLRECP,MLRECV,NLCG,NLCD)
  301. C-----------------------------------------------------
  302. C********** AB.AM(IFAC,IPRIM,IDUAL)
  303. C A = nom de l'inconnu duale (Ro,rUX,rUY,RET)
  304. C B = nom de l'inconnu primale (Ro,rUX,rUY,RET)
  305. C IPRIM = 1, 2 -> G, D
  306. C IDUAL = 1, 2 -> G, D
  307. C i.e.
  308. C A_IDUAL = AB.AM(IFAC,IPRIM,IDUAL) * B_IPRIM + ...
  309. C
  310. C
  311. C********** Dual RN
  312. C----------------------------------------------------
  313. JTL=JLL
  314. JTR=JRR
  315. SEGACT JTL
  316. SEGACT JTR
  317. DO 10 II = 1,(4+NSP)
  318. DO 15 JJ = 1,(4+NSP)
  319. KV = (II-1)*(4+NSP)
  320. C----------------------------------
  321. CELL = IMATRI.LIZAFM(1,KV+JJ)
  322. FUNCEL = SURF * JTL.JAC(II,JJ)
  323. CELL.AM(IFAC,1,1) = -1.0D0 * FUNCEL / VOLG
  324. CELL.AM(IFAC,1,2) = FUNCEL / VOLD
  325. 15 CONTINUE
  326. 10 CONTINUE
  327. C-----------------------------------------------------
  328. DO 20 II = 1,(4+NSP)
  329. DO 25 JJ = 1,(4+NSP)
  330. KV = (II-1)*(4+NSP)
  331. CELL = IMATRI.LIZAFM(1,KV+JJ)
  332. FUNCEL = SURF * JTR.JAC(II,JJ)
  333. CELL.AM(IFAC,2,2) = FUNCEL / VOLD
  334. CELL.AM(IFAC,2,1) = -FUNCEL / VOLG
  335. 25 CONTINUE
  336. 20 CONTINUE
  337. SEGDES JTR
  338. SEGDES JTL
  339. C-----------------------------------------------------
  340. ELSE
  341. C-----------------------------------------------------
  342. C************* Murs (NGCG = NGCD) ******************
  343. C-----------------------------------------------------
  344. MELEDU.NUM(1,IFAC) = NGCG
  345. MELEDU.NUM(2,IFAC) = NGCD
  346. NLCG = MLENTC.LECT(NGCG)
  347. C--------------------------------------
  348. ROG = MPRN.VPOCHA(NLCG,1)
  349. PG = MPPN.VPOCHA(NLCG,1)
  350. UXG = MPUN.VPOCHA(NLCG,1)
  351. UYG = MPUN.VPOCHA(NLCG,2)
  352. UZG = MPUN.VPOCHA(NLCG,3)
  353. VOLG = MPVOLU.VPOCHA(NLCG,1)
  354. C-------------------------------------------
  355. WVEC_L(1)=ROG
  356. WVEC_L(2)=UXG
  357. WVEC_L(3)=UYG
  358. WVEC_L(4)=UZG
  359. WVEC_L(5)=PG
  360. C-------------------------------------------------
  361. SURF = MPOVSU.VPOCHA(NLCF,1)
  362. NVECT(1) = MPNORM.VPOCHA(NLCF,7)
  363. NVECT(2) = MPNORM.VPOCHA(NLCF,8)
  364. NVECT(3) = MPNORM.VPOCHA(NLCF,9)
  365. c--------------------------------------------
  366. TVECT1(1) = MPNORM.VPOCHA(NLCF,1)
  367. TVECT1(2) = MPNORM.VPOCHA(NLCF,2)
  368. TVECT1(3) = MPNORM.VPOCHA(NLCF,3)
  369. c----------------------------------------------
  370. TVECT2(1) = MPNORM.VPOCHA(NLCF,4)
  371. TVECT2(2) = MPNORM.VPOCHA(NLCF,5)
  372. TVECT2(3) = MPNORM.VPOCHA(NLCF,6)
  373. C------- COEFFICIENTS ----------------------------
  374. C11=TVECT1(2)*TVECT2(3)-TVECT1(3)*TVECT2(2)
  375. C12=NVECT(2)*TVECT2(3)-TVECT2(2)*NVECT(3)
  376. C13=NVECT(2)*TVECT1(3)-TVECT1(2)*NVECT(3)
  377. C---------------------------------
  378. C21=TVECT1(1)*TVECT2(3)-TVECT1(3)*TVECT2(1)
  379. C22=NVECT(1)*TVECT2(3)-TVECT2(1)*NVECT(3)
  380. C23=NVECT(1)*TVECT1(3)-TVECT1(1)*NVECT(3)
  381. C---------------------------------
  382. C31=TVECT1(1)*TVECT2(2)-TVECT1(2)*TVECT2(1)
  383. C32=NVECT(1)*TVECT2(2)-TVECT2(1)*NVECT(2)
  384. C33=NVECT(1)*TVECT1(2)-TVECT1(1)*NVECT(2)
  385. DET=NVECT(1)*C11-NVECT(2)*C21+NVECT(3)*C31
  386. C---------------------------------
  387. ZC11=-NVECT(1)*C11-TVECT1(1)*C12+TVECT2(1)*C13
  388. ZC12=-NVECT(2)*C11-TVECT1(2)*C12+TVECT2(2)*C13
  389. ZC13=-NVECT(3)*C11-TVECT1(3)*C12+TVECT2(3)*C13
  390. C---------------------------------
  391. ZC21=NVECT(1)*C21+TVECT1(1)*C22-TVECT2(1)*C23
  392. ZC22=NVECT(2)*C21+TVECT1(2)*C22-TVECT2(2)*C23
  393. ZC23=NVECT(3)*C21+TVECT1(3)*C22-TVECT2(3)*C23
  394. C---------------------------------
  395. ZC31=-NVECT(1)*C31-TVECT1(1)*C32+TVECT2(1)*C33
  396. ZC32=-NVECT(2)*C31-TVECT1(2)*C32+TVECT2(2)*C33
  397. ZC33=-NVECT(3)*C31-TVECT1(3)*C32+TVECT2(3)*C33
  398. C-------------------------------------------------
  399. ROD = ROG
  400. PD = PG
  401. UXD = (ZC11*UXG+ZC12*UYG+ZC13*UZG)/DET
  402. UYD = (ZC21*UXG+ZC22*UYG+ZC23*UZG)/DET
  403. UZD = (ZC31*UXG+ZC32*UYG+ZC33*UZG)/DET
  404. VOLD = VOLG
  405. C------------------------------------------------
  406. WVEC_R(1)=ROD
  407. WVEC_R(2)=UXD
  408. WVEC_R(3)=UYD
  409. WVEC_R(4)=UZD
  410. WVEC_R(5)=PD
  411. C-------------------------------------------
  412. C********** La normale sortante
  413. C-------------------------------------------
  414. CALL cwms3d(NSP,JLL,WVEC_L,WVEC_R,NVECT,
  415. & TVECT1,TVECT2,MPYN,MLRECP,MLRECV,NLCG)
  416. C----------------------------------------------------
  417. JTL=JLL
  418. SEGACT JTL
  419. DO 70 II = 1,(4+NSP)
  420. DO 75 JJ = 1,(4+NSP)
  421. KV = (II-1)*(4+NSP)
  422. C---------------------------------
  423. CELL = IMATRI.LIZAFM(1,KV+JJ)
  424. FUNCEL = SURF * JTL.JAC(II,JJ)
  425. CELL.AM(IFAC,1,1) = -1.0D0 * FUNCEL / VOLG
  426. CELL.AM(IFAC,1,2) = 0.0D0
  427. 75 CONTINUE
  428. 70 CONTINUE
  429. SEGDES JTL
  430. C--------------------------------------------
  431. DO 40 II = 1,(4+NSP)
  432. DO 45 JJ = 1,(4+NSP)
  433. KV = (II-1)*(4+NSP)
  434. CELL = IMATRI.LIZAFM(1,KV+JJ)
  435. CELL.AM(IFAC,2,2) = 0.0D0
  436. CELL.AM(IFAC,2,1) = 0.0D0
  437. 45 CONTINUE
  438. 40 CONTINUE
  439. C--------------------------------------------
  440. ENDIF
  441. ENDDO
  442. C
  443. SEGDES MELEMC
  444. SEGDES MELEFE
  445. SEGDES MELEMF
  446. C
  447. SEGDES MPOVSU
  448. SEGDES MPVOLU
  449. SEGDES MPNORM
  450. C
  451. SEGDES MPRN
  452. SEGDES MPPN
  453. SEGDES MPUN
  454. SEGDES MPYN
  455. C
  456. SEGDES MELEDU
  457. SEGDES MATRIK
  458. DO 80 II=1,NBME
  459. CELL = IMATRI.LIZAFM(1,II)
  460. SEGDES CELL
  461. 80 CONTINUE
  462. SEGDES IMATRI
  463. C
  464. SEGSUP MLENTC
  465. SEGSUP MLENTF
  466. SEGSUP MLELIM
  467. SEGDES MLMINC
  468. IF(MELLIM .NE.0) SEGDES MELLIM
  469.  
  470. 9999 CONTINUE
  471. RETURN
  472. END
  473.  
  474.  
  475.  
  476.  
  477.  
  478.  
  479.  
  480.  
  481.  
  482.  
  483.  
  484.  
  485.  
  486.  
  487.  
  488.  
  489.  

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