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. -INC CCOPTIO
  109. -INC SMCHPOI
  110. -INC SMELEME
  111. -INC SMLMOTS
  112. -INC SMLENTI
  113. POINTEUR MPRN.MPOVAL, MPUN.MPOVAL, MPPN.MPOVAL, MPYN.MPOVAL,
  114. & MPNORM.MPOVAL, MPVOLU.MPOVAL, MPOVSU.MPOVAL
  115. POINTEUR MELEMC.MELEME, MELEMF.MELEME, MELEFE.MELEME,
  116. & MELEDU.MELEME, MELLIM.MELEME
  117. POINTEUR MLENTC.MLENTI, MLENTF.MLENTI, MLELIM.MLENTI
  118. POINTEUR CELL.IZAFM
  119. POINTEUR MLMINC.MLMOTS
  120. C----------------------------------------------------
  121. -INC SMLREEL
  122. POINTEUR MLRECP.MLREEL, MLRECV.MLREEL
  123. C-------------------------------------------------------
  124. C********* Les Jacobians ******************************
  125. C-------------------------------------------------------
  126. SEGMENT JACEL
  127. REAL*8 JAC(4+NSP,4+NSP)
  128. ENDSEGMENT
  129. POINTEUR JTL.JACEL, JTR.JACEL, JTT.JACEL
  130. C--------------------------------------------------------
  131. C KRIPAD pour la correspondance global/local des centres
  132. C--------------------------------------------------------
  133. CALL KRIPAD(MELLIM,MLELIM)
  134. CALL KRIPAD(MELEMC,MLENTC)
  135. C------------------
  136. SEGACT MELEMC
  137. C------------------
  138. SEGACT MELEFE
  139. C------------------
  140. CALL LICHT(ICHPSU,MPOVSU,TYPE,IGEOMF)
  141. CALL LICHT(INORM,MPNORM,TYPE,IGEOMF)
  142. CALL LICHT(ICHPVO,MPVOLU,TYPE,IGEOMC)
  143. C----------------------------------------------
  144. MELEMF = IGEOMF
  145. CALL KRIPAD(MELEMF,MLENTF)
  146. SEGACT MELEMF
  147. C-----------------------------------------------
  148. CALL LICHT(IRN,MPRN,TYPE,IGEOMC)
  149. CALL LICHT(IPN,MPPN,TYPE,IGEOMC)
  150. CALL LICHT(IUN,MPUN,TYPE,IGEOMC)
  151. CALL LICHT(IYN,MPYN,TYPE,IGEOMC)
  152. C-----------------------------------------------
  153. NFAC = MELEFE.NUM(/2)
  154. C-----------------------------------------------
  155. C**** Maillage des inconnues primales
  156. C-----------------------------------------------
  157. NBSOUS = 0
  158. NBREF = 0
  159. NBELEM = NFAC
  160. NBNN = 2
  161. C-----------------------------------------------
  162. SEGINI MELEDU
  163. C----------------------
  164. C** MELEDU = 'SEG2' **
  165. C----------------------
  166. MELEDU.ITYPEL = 2
  167. C--------------
  168. NRIGE = 7
  169. NMATRI = 1
  170. NKID = 9
  171. NKMT = 7
  172. C--------------
  173. SEGINI MATRIK
  174. IMAT = MATRIK
  175. MATRIK.IRIGEL(1,1) = MELEDU
  176. MATRIK.IRIGEL(2,1) = MELEDU
  177. C-------------------------
  178. C Matrice non symetrique
  179. C-------------------------
  180. MATRIK.IRIGEL(7,1) = 2
  181. C-------------------------
  182. NBME = (4+NSP)*(4+NSP)
  183. NBSOUS = 1
  184. SEGINI IMATRI
  185. MLMINC = ILIINC
  186. SEGACT MLMINC
  187. MATRIK.IRIGEL(4,1) = IMATRI
  188. C-----------------------------------------------
  189. DO 1 J=1,(NSP+4)
  190. KV=(J-1)*(4+NSP)
  191. IMATRI.LISPRI(KV+1) = MLMINC.MOTS(1)
  192. IMATRI.LISPRI(KV+2) = MLMINC.MOTS(2)
  193. IMATRI.LISPRI(KV+3) = MLMINC.MOTS(3)
  194. IMATRI.LISPRI(KV+4) = MLMINC.MOTS(4)
  195. IMATRI.LISPRI(KV+5) = MLMINC.MOTS(5)
  196. DO 2 I=1,(NSP-1)
  197. IMATRI.LISPRI(KV+5+I) = MLMINC.MOTS(5+I)
  198. 2 CONTINUE
  199. 1 CONTINUE
  200. C-----------------------------------------------
  201. DO 3 J=1,(NSP+4)
  202. KV=(J-1)*(4+NSP)
  203. IMATRI.LISDUA(KV+1) = MLMINC.MOTS(j)
  204. IMATRI.LISDUA(KV+2) = MLMINC.MOTS(j)
  205. IMATRI.LISDUA(KV+3) = MLMINC.MOTS(j)
  206. IMATRI.LISDUA(KV+4) = MLMINC.MOTS(j)
  207. IMATRI.LISDUA(KV+5) = MLMINC.MOTS(j)
  208. DO 4 I=1,(NSP-1)
  209. IMATRI.LISDUA(KV+5+I) = MLMINC.MOTS(j)
  210. 4 CONTINUE
  211. 3 CONTINUE
  212. C-----------------------------------------------
  213. C-----------------------------------------------
  214. NBEL = NBELEM
  215. NBSOUS = 1
  216. NP = 2
  217. MP = 2
  218. C-----------------------------------------------------------
  219. C-----------------------------------------------------------
  220. DO 5 I=1,NBME
  221. SEGINI CELL
  222. IMATRI.LIZAFM(1,I) = CELL
  223. 5 CONTINUE
  224. C---------------------------------
  225. C---------------------------------
  226. DO IFAC = 1, NFAC, 1
  227. NGCF = MELEFE.NUM(2,IFAC)
  228. NLCF = MLENTF.LECT(NGCF)
  229. IF(NLCF .NE. IFAC)THEN
  230. WRITE(IOIMP,*) 'Il ne faut pas jouer avec la table domaine'
  231. CALL ERREUR(5)
  232. GOTO 9999
  233. ENDIF
  234. NLFL = MLELIM.LECT(NGCF)
  235. NGCG = MELEFE.NUM(1,IFAC)
  236. NGCD = MELEFE.NUM(3,IFAC)
  237. IF(NLFL .NE. 0)THEN
  238. C--------------------------------------------------------
  239. C The point belongs on BC -> No contribution to jacobian!
  240. C--------------------------------------------------------
  241. MELEDU.NUM(1,IFAC) = NGCG
  242. MELEDU.NUM(2,IFAC) = NGCD
  243. ELSEIF(NGCG .NE. NGCD)THEN
  244. C-----------------------------------
  245. C********** Les MELEMEs
  246. C-----------------------------------
  247. MELEDU.NUM(1,IFAC) = NGCG
  248. MELEDU.NUM(2,IFAC) = NGCD
  249. C-----------------------------------
  250. C********** Les etats G et D
  251. C-----------------------------------
  252. NLCG = MLENTC.LECT(NGCG)
  253. NLCD = MLENTC.LECT(NGCD)
  254. C-----------------
  255. ROG = MPRN.VPOCHA(NLCG,1)
  256. PG = MPPN.VPOCHA(NLCG,1)
  257. UXG = MPUN.VPOCHA(NLCG,1)
  258. UYG = MPUN.VPOCHA(NLCG,2)
  259. UZG = MPUN.VPOCHA(NLCG,3)
  260. VOLG = MPVOLU.VPOCHA(NLCG,1)
  261. C-------------------------------------------------
  262. WVEC_L(1)=ROG
  263. WVEC_L(2)=UXG
  264. WVEC_L(3)=UYG
  265. WVEC_L(4)=UZG
  266. WVEC_L(5)=PG
  267. C-------------------------------------------------
  268. ROD = MPRN.VPOCHA(NLCD,1)
  269. PD = MPPN.VPOCHA(NLCD,1)
  270. UXD = MPUN.VPOCHA(NLCD,1)
  271. UYD = MPUN.VPOCHA(NLCD,2)
  272. UZD = MPUN.VPOCHA(NLCD,3)
  273. VOLD = MPVOLU.VPOCHA(NLCD,1)
  274. C------------------------------------------------
  275. WVEC_R(1)=ROD
  276. WVEC_R(2)=UXD
  277. WVEC_R(3)=UYD
  278. WVEC_R(4)=UZD
  279. WVEC_R(5)=PD
  280. C------------------------------------------------
  281. C La normale G->D
  282. C La tangente
  283. C------------------------------------------------
  284. SURF = MPOVSU.VPOCHA(NLCF,1)
  285. NVECT(1) = MPNORM.VPOCHA(NLCF,7)
  286. NVECT(2) = MPNORM.VPOCHA(NLCF,8)
  287. NVECT(3) = MPNORM.VPOCHA(NLCF,9)
  288. c-----------------------------------------------
  289. TVECT1(1) = MPNORM.VPOCHA(NLCF,1)
  290. TVECT1(2) = MPNORM.VPOCHA(NLCF,2)
  291. TVECT1(3) = MPNORM.VPOCHA(NLCF,3)
  292. c----------------------------------------------
  293. TVECT2(1) = MPNORM.VPOCHA(NLCF,4)
  294. TVECT2(2) = MPNORM.VPOCHA(NLCF,5)
  295. TVECT2(3) = MPNORM.VPOCHA(NLCF,6)
  296. C-------------------------------------------------
  297. CALL cms3d(NSP,JLL,JRR,WVEC_L,WVEC_R,NVECT,TVECT1,TVECT2,
  298. & MPYN,MLRECP,MLRECV,NLCG,NLCD)
  299. C-----------------------------------------------------
  300. C********** AB.AM(IFAC,IPRIM,IDUAL)
  301. C A = nom de l'inconnu duale (Ro,rUX,rUY,RET)
  302. C B = nom de l'inconnu primale (Ro,rUX,rUY,RET)
  303. C IPRIM = 1, 2 -> G, D
  304. C IDUAL = 1, 2 -> G, D
  305. C i.e.
  306. C A_IDUAL = AB.AM(IFAC,IPRIM,IDUAL) * B_IPRIM + ...
  307. C
  308. C
  309. C********** Dual RN
  310. C----------------------------------------------------
  311. JTL=JLL
  312. JTR=JRR
  313. SEGACT JTL
  314. SEGACT JTR
  315. DO 10 II = 1,(4+NSP)
  316. DO 15 JJ = 1,(4+NSP)
  317. KV = (II-1)*(4+NSP)
  318. C----------------------------------
  319. CELL = IMATRI.LIZAFM(1,KV+JJ)
  320. FUNCEL = SURF * JTL.JAC(II,JJ)
  321. CELL.AM(IFAC,1,1) = -1.0D0 * FUNCEL / VOLG
  322. CELL.AM(IFAC,1,2) = FUNCEL / VOLD
  323. 15 CONTINUE
  324. 10 CONTINUE
  325. C-----------------------------------------------------
  326. DO 20 II = 1,(4+NSP)
  327. DO 25 JJ = 1,(4+NSP)
  328. KV = (II-1)*(4+NSP)
  329. CELL = IMATRI.LIZAFM(1,KV+JJ)
  330. FUNCEL = SURF * JTR.JAC(II,JJ)
  331. CELL.AM(IFAC,2,2) = FUNCEL / VOLD
  332. CELL.AM(IFAC,2,1) = -FUNCEL / VOLG
  333. 25 CONTINUE
  334. 20 CONTINUE
  335. SEGDES JTR
  336. SEGDES JTL
  337. C-----------------------------------------------------
  338. ELSE
  339. C-----------------------------------------------------
  340. C************* Murs (NGCG = NGCD) ******************
  341. C-----------------------------------------------------
  342. MELEDU.NUM(1,IFAC) = NGCG
  343. MELEDU.NUM(2,IFAC) = NGCD
  344. NLCG = MLENTC.LECT(NGCG)
  345. C--------------------------------------
  346. ROG = MPRN.VPOCHA(NLCG,1)
  347. PG = MPPN.VPOCHA(NLCG,1)
  348. UXG = MPUN.VPOCHA(NLCG,1)
  349. UYG = MPUN.VPOCHA(NLCG,2)
  350. UZG = MPUN.VPOCHA(NLCG,3)
  351. VOLG = MPVOLU.VPOCHA(NLCG,1)
  352. C-------------------------------------------
  353. WVEC_L(1)=ROG
  354. WVEC_L(2)=UXG
  355. WVEC_L(3)=UYG
  356. WVEC_L(4)=UZG
  357. WVEC_L(5)=PG
  358. C-------------------------------------------------
  359. SURF = MPOVSU.VPOCHA(NLCF,1)
  360. NVECT(1) = MPNORM.VPOCHA(NLCF,7)
  361. NVECT(2) = MPNORM.VPOCHA(NLCF,8)
  362. NVECT(3) = MPNORM.VPOCHA(NLCF,9)
  363. c--------------------------------------------
  364. TVECT1(1) = MPNORM.VPOCHA(NLCF,1)
  365. TVECT1(2) = MPNORM.VPOCHA(NLCF,2)
  366. TVECT1(3) = MPNORM.VPOCHA(NLCF,3)
  367. c----------------------------------------------
  368. TVECT2(1) = MPNORM.VPOCHA(NLCF,4)
  369. TVECT2(2) = MPNORM.VPOCHA(NLCF,5)
  370. TVECT2(3) = MPNORM.VPOCHA(NLCF,6)
  371. C------- COEFFICIENTS ----------------------------
  372. C11=TVECT1(2)*TVECT2(3)-TVECT1(3)*TVECT2(2)
  373. C12=NVECT(2)*TVECT2(3)-TVECT2(2)*NVECT(3)
  374. C13=NVECT(2)*TVECT1(3)-TVECT1(2)*NVECT(3)
  375. C---------------------------------
  376. C21=TVECT1(1)*TVECT2(3)-TVECT1(3)*TVECT2(1)
  377. C22=NVECT(1)*TVECT2(3)-TVECT2(1)*NVECT(3)
  378. C23=NVECT(1)*TVECT1(3)-TVECT1(1)*NVECT(3)
  379. C---------------------------------
  380. C31=TVECT1(1)*TVECT2(2)-TVECT1(2)*TVECT2(1)
  381. C32=NVECT(1)*TVECT2(2)-TVECT2(1)*NVECT(2)
  382. C33=NVECT(1)*TVECT1(2)-TVECT1(1)*NVECT(2)
  383. DET=NVECT(1)*C11-NVECT(2)*C21+NVECT(3)*C31
  384. C---------------------------------
  385. ZC11=-NVECT(1)*C11-TVECT1(1)*C12+TVECT2(1)*C13
  386. ZC12=-NVECT(2)*C11-TVECT1(2)*C12+TVECT2(2)*C13
  387. ZC13=-NVECT(3)*C11-TVECT1(3)*C12+TVECT2(3)*C13
  388. C---------------------------------
  389. ZC21=NVECT(1)*C21+TVECT1(1)*C22-TVECT2(1)*C23
  390. ZC22=NVECT(2)*C21+TVECT1(2)*C22-TVECT2(2)*C23
  391. ZC23=NVECT(3)*C21+TVECT1(3)*C22-TVECT2(3)*C23
  392. C---------------------------------
  393. ZC31=-NVECT(1)*C31-TVECT1(1)*C32+TVECT2(1)*C33
  394. ZC32=-NVECT(2)*C31-TVECT1(2)*C32+TVECT2(2)*C33
  395. ZC33=-NVECT(3)*C31-TVECT1(3)*C32+TVECT2(3)*C33
  396. C-------------------------------------------------
  397. ROD = ROG
  398. PD = PG
  399. UXD = (ZC11*UXG+ZC12*UYG+ZC13*UZG)/DET
  400. UYD = (ZC21*UXG+ZC22*UYG+ZC23*UZG)/DET
  401. UZD = (ZC31*UXG+ZC32*UYG+ZC33*UZG)/DET
  402. VOLD = VOLG
  403. C------------------------------------------------
  404. WVEC_R(1)=ROD
  405. WVEC_R(2)=UXD
  406. WVEC_R(3)=UYD
  407. WVEC_R(4)=UZD
  408. WVEC_R(5)=PD
  409. C-------------------------------------------
  410. C********** La normale sortante
  411. C-------------------------------------------
  412. CALL cwms3d(NSP,JLL,WVEC_L,WVEC_R,NVECT,
  413. & TVECT1,TVECT2,MPYN,MLRECP,MLRECV,NLCG)
  414. C----------------------------------------------------
  415. JTL=JLL
  416. SEGACT JTL
  417. DO 70 II = 1,(4+NSP)
  418. DO 75 JJ = 1,(4+NSP)
  419. KV = (II-1)*(4+NSP)
  420. C---------------------------------
  421. CELL = IMATRI.LIZAFM(1,KV+JJ)
  422. FUNCEL = SURF * JTL.JAC(II,JJ)
  423. CELL.AM(IFAC,1,1) = -1.0D0 * FUNCEL / VOLG
  424. CELL.AM(IFAC,1,2) = 0.0D0
  425. 75 CONTINUE
  426. 70 CONTINUE
  427. SEGDES JTL
  428. C--------------------------------------------
  429. DO 40 II = 1,(4+NSP)
  430. DO 45 JJ = 1,(4+NSP)
  431. KV = (II-1)*(4+NSP)
  432. CELL = IMATRI.LIZAFM(1,KV+JJ)
  433. CELL.AM(IFAC,2,2) = 0.0D0
  434. CELL.AM(IFAC,2,1) = 0.0D0
  435. 45 CONTINUE
  436. 40 CONTINUE
  437. C--------------------------------------------
  438. ENDIF
  439. ENDDO
  440. C
  441. SEGDES MELEMC
  442. SEGDES MELEFE
  443. SEGDES MELEMF
  444. C
  445. SEGDES MPOVSU
  446. SEGDES MPVOLU
  447. SEGDES MPNORM
  448. C
  449. SEGDES MPRN
  450. SEGDES MPPN
  451. SEGDES MPUN
  452. SEGDES MPYN
  453. C
  454. SEGDES MELEDU
  455. SEGDES MATRIK
  456. DO 80 II=1,NBME
  457. CELL = IMATRI.LIZAFM(1,II)
  458. SEGDES CELL
  459. 80 CONTINUE
  460. SEGDES IMATRI
  461. C
  462. SEGSUP MLENTC
  463. SEGSUP MLENTF
  464. SEGSUP MLELIM
  465. SEGDES MLMINC
  466. IF(MELLIM .NE.0) SEGDES MELLIM
  467.  
  468. 9999 CONTINUE
  469. RETURN
  470. END
  471.  
  472.  
  473.  
  474.  
  475.  
  476.  
  477.  
  478.  
  479.  
  480.  
  481.  
  482.  
  483.  
  484.  
  485.  
  486.  
  487.  

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