1 : C KBMMS2 SOURCE MAGN 09/07/30 21:17:30 6458
2 : SUBROUTINE KBMMS2(NSP,ILIINC,IRN,IUN,IPN,IYN,MLRECP,
3 : & MLRECV,INORM,ICHPVO,ICHPSU,IUINF,IUPRI,
4 : & MELEMC,MELEFE,MELLIM,IMAT)
5 : C
6 : C************************************************************************
7 : C
8 : C PROJET : CASTEM 2000
9 : C
10 : C NOM : KONMSP ("convection for multispecies")
11 : C
12 : C DESCRIPTION : Voir KON19 (appele par KON19)
13 : C Calcul du jacobien du résidu pour la méthode
14 : C AUSMplus
15 : C
15 : C Cas deux 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) : CONMSP, CWMSP
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 ICHPOVO : CHPOINT VOLUME contenant le volume
51 : C
52 : C ICHPOSU : 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 N.B.: Sur le bord on ne calcule rien!
69 : C
70 : C SORTIES
71 : C
72 : C IMAT : pointeur de la MATRIK du jacobien du residu
73 : C
74 : C************************************************************************
75 : C
76 : C HISTORIQUE (Anomalies et modifications éventuelles)
77 : C
78 : C HISTORIQUE :
79 : C
80 : C************************************************************************
81 : C
82 : C
83 : C N.B.: On suppose qu'on a déjà controllé RO, P > 0
84 : C GAMMA \in (1,3)
85 : C Si non il faut le faire!!!
86 : C
87 : C************************************************************************
88 : C
89 : C
90 : C**** Variables de COOPTIO
91 : C
92 : IMPLICIT INTEGER(I-N)
93 : INTEGER IPLLB, IERPER, IERMAX, IERR, INTERR
94 : & ,IOTER, IOLEC, IOIMP, IOCAR, IOACQ
95 : & ,IOPER, IOSGB, IOGRA, IOSAU, IORES
96 : & ,IECHO, IIMPI, IOSPI
97 : & ,IDIM
98 : & ,MCOORD
99 : & ,IFOMOD, NIFOUR, IFOUR, NSDPGE, IONIVE
100 : & ,NGMAXY, IZROSF, ISOTYP, IOSCR,LTEXLU
101 : & ,NORINC,NORVAL,NORIND,NORVAD
102 : & ,NUCROU, IPSAUV
103 : C---------------------------------------------------
104 : INTEGER NSP,II,JJ,KV,I,J,IFICLE,IPREFI,JLL,JRR
105 : INTEGER ILIINC, IRN,IUN,IPN,IYN,INORM,ICHPVO,ICHPSU
106 : & , IMAT, IGEOMC, IGEOMF, IUPRI, IUINF
107 : & , NFAC, NBSOUS, NBREF, NBELEM, NBNN, NRIGE, NMATRI, NKID
108 : & , NKMT, NBME, NBEL, MP, NP
109 : & , IFAC, NGCF, NLCF, NGCG, NGCD, NLCG, NLCD, NLFL
110 : REAL*8 ROG, PG, UXG, UYG, VOLG
111 : & , ROD, PD, UXD, UYD, VOLD
112 : & , SURF, FUNCEL, UPR_L, UPR_R, UPR_M, V_INF
113 : REAL*8 WVEC_L(4), WVEC_R(4), NVECT(2), TVECT(2)
114 : REAL*8 ZC1, ZC2, ZC3, ZC4
115 : CHARACTER*8 TYPE
116 : C
117 : C**** LES INCLUDES
118 : C
119 : -INC CCOPTIO
120 : -INC SMCHPOI
121 : -INC SMMATRIK
122 : -INC SMELEME
123 : -INC SMLMOTS
124 : -INC SMLENTI
125 : POINTEUR MPRN.MPOVAL, MPUN.MPOVAL, MPPN.MPOVAL, MPYN.MPOVAL,
126 : & MPNORM.MPOVAL, MPVOLU.MPOVAL, MPOVSU.MPOVAL
127 : POINTEUR MELEMC.MELEME, MELEMF.MELEME, MELEFE.MELEME,
128 : & MELEDU.MELEME, MELLIM.MELEME
129 : POINTEUR MLENTC.MLENTI, MLENTF.MLENTI, MLELIM.MLENTI
130 : POINTEUR CELL.IZAFM
131 : POINTEUR MLMINC.MLMOTS, MPUPRI.MPOVAL, MPUINF.MPOVAL
132 : C----------------------------------------------------
133 : -INC SMLREEL
134 : POINTEUR MLRECP.MLREEL, MLRECV.MLREEL
135 : C-------------------------------------------------------
136 : C********* Les Jacobians ******************************
137 : C-------------------------------------------------------
138 : SEGMENT JACEL
139 : REAL*8 JAC(3+NSP,3+NSP)
140 : ENDSEGMENT
141 : POINTEUR JTL.JACEL, JTR.JACEL, JTT.JACEL
142 : C--------------------------------------------------------
143 : C KRIPAD pour la correspondance global/local des centres
144 : C--------------------------------------------------------
144 : CALL KRIPAD(MELLIM,MLELIM)
144 : CALL KRIPAD(MELEMC,MLENTC)
145 : C------------------
146 : SEGACT MELEMC
147 : C------------------
148 : SEGACT MELEFE
149 : C------------------
149 : CALL LICHT(ICHPSU,MPOVSU,TYPE,IGEOMF)
149 : CALL LICHT(INORM,MPNORM,TYPE,IGEOMF)
149 : CALL LICHT(ICHPVO,MPVOLU,TYPE,IGEOMC)
150 : C----------------------------------------------
151 : MELEMF = IGEOMF
151 : CALL KRIPAD(MELEMF,MLENTF)
152 : SEGACT MELEMF
153 : C-----------------------------------------------
153 : CALL LICHT(IRN,MPRN,TYPE,IGEOMC)
153 : CALL LICHT(IPN,MPPN,TYPE,IGEOMC)
153 : 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 = (3+NSP)*(3+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+3)
192 : KV=(J-1)*(3+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 : DO 2 I=1,(NSP-1)
198 : IMATRI.LISPRI(KV+4+I) = MLMINC.MOTS(4+I)
199 : 2 CONTINUE
200 : 1 CONTINUE
201 : C-----------------------------------------------
202 : DO 3 J=1,(NSP+3)
203 : KV=(J-1)*(3+NSP)
204 : IMATRI.LISDUA(KV+1) = MLMINC.MOTS(j)
205 : IMATRI.LISDUA(KV+2) = MLMINC.MOTS(j)
206 : IMATRI.LISDUA(KV+3) = MLMINC.MOTS(j)
207 : IMATRI.LISDUA(KV+4) = MLMINC.MOTS(j)
208 : DO 4 I=1,(NSP-1)
209 : IMATRI.LISDUA(KV+4+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 Bas Mach
226 : C**************************************************************
226 : CALL LICHT(IUPRI,MPUPRI,TYPE,IGEOMC)
226 : CALL LICHT(IUINF,MPUINF,TYPE,IGEOMC)
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'
232 : CALL ERREUR(5)
233 : GOTO 9999
234 : ENDIF
235 : NLFL = MLELIM.LECT(NGCF)
236 : NGCG = MELEFE.NUM(1,IFAC)
237 : NGCD = MELEFE.NUM(3,IFAC)
238 : IF(NLFL .NE. 0)THEN
239 : C---------------------------------------------------------
240 : C The point belongs to BC -> No contribution to jacobian!
241 : C---------------------------------------------------------
242 : MELEDU.NUM(1,IFAC) = NGCG
243 : MELEDU.NUM(2,IFAC) = NGCD
244 : ELSEIF(NGCG .NE. NGCD)THEN
245 : C-----------------------------------
246 : C********** Les MELEMEs
247 : C-----------------------------------
248 : MELEDU.NUM(1,IFAC) = NGCG
249 : MELEDU.NUM(2,IFAC) = NGCD
250 : C-----------------------------------
251 : C********** Les etats G et D
252 : C-----------------------------------
253 : NLCG = MLENTC.LECT(NGCG)
254 : NLCD = MLENTC.LECT(NGCD)
255 : C-----------------
256 : ROG = MPRN.VPOCHA(NLCG,1)
257 : PG = MPPN.VPOCHA(NLCG,1)
258 : UXG = MPUN.VPOCHA(NLCG,1)
259 : UYG = MPUN.VPOCHA(NLCG,2)
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)=PG
266 : C-------------------------------------------------
267 : ROD = MPRN.VPOCHA(NLCD,1)
268 : PD = MPPN.VPOCHA(NLCD,1)
269 : UXD = MPUN.VPOCHA(NLCD,1)
270 : UYD = MPUN.VPOCHA(NLCD,2)
271 : VOLD = MPVOLU.VPOCHA(NLCD,1)
272 : C------------------------------------------------
273 : UPR_L=MPUPRI.VPOCHA(NLCG,1)
274 : UPR_R=MPUPRI.VPOCHA(NLCD,1)
275 : UPR_M=MAX(UPR_L,UPR_R)
276 : V_INF=MAX(MPUINF.VPOCHA(NLCG,1),MPUINF.VPOCHA(NLCD,1))
277 : V_INF=MAX(UPR_M,V_INF)
278 : C------------------------------------------------
279 : WVEC_R(1)=ROD
280 : WVEC_R(2)=UXD
281 : WVEC_R(3)=UYD
282 : WVEC_R(4)=PD
283 : C------------------------------------------------
284 : C La normale G->D
285 : C La tangente
286 : C------------------------------------------------
287 : SURF = MPOVSU.VPOCHA(NLCF,1)
288 : NVECT(1) = MPNORM.VPOCHA(NLCF,1)
289 : NVECT(2) = MPNORM.VPOCHA(NLCF,2)
290 : TVECT(1) = -1.0D0 * NVECT(2)
291 : TVECT(2) = NVECT(1)
292 : C-------------------------------------------------
292 : CALL jbmms2(NSP,JLL,JRR,WVEC_L,WVEC_R,NVECT,TVECT,
293 : & MPYN,MLRECP,MLRECV,NLCG,NLCD,V_INF)
294 : C-----------------------------------------------------
295 : C********** AB.AM(IFAC,IPRIM,IDUAL)
296 : C A = nom de l'inconnu duale (Ro,rUX,rUY,RET)
297 : C B = nom de l'inconnu primale (Ro,rUX,rUY,RET)
298 : C IPRIM = 1, 2 -> G, D
299 : C IDUAL = 1, 2 -> G, D
300 : C i.e.
301 : C A_IDUAL = AB.AM(IFAC,IPRIM,IDUAL) * B_IPRIM + ...
302 : C
303 : C
304 : C********** Dual RN
305 : C----------------------------------------------------
306 : JTL=JLL
307 : JTR=JRR
308 : SEGACT JTL
309 : SEGACT JTR
310 : DO 10 II = 1,(3+NSP)
311 : DO 15 JJ = 1,(3+NSP)
312 : KV = (II-1)*(3+NSP)
313 : C----------------------------------
314 : CELL = IMATRI.LIZAFM(1,KV+JJ)
315 : FUNCEL = SURF * JTL.JAC(II,JJ)
316 : CELL.AM(IFAC,1,1) = -1.0D0 * FUNCEL / VOLG
317 : CELL.AM(IFAC,1,2) = FUNCEL / VOLD
318 : 15 CONTINUE
319 : 10 CONTINUE
320 : C-----------------------------------------------------
321 : DO 20 II = 1,(3+NSP)
322 : DO 25 JJ = 1,(3+NSP)
323 : KV = (II-1)*(3+NSP)
324 : CELL = IMATRI.LIZAFM(1,KV+JJ)
325 : FUNCEL = SURF * JTR.JAC(II,JJ)
326 : CELL.AM(IFAC,2,2) = FUNCEL / VOLD
327 : CELL.AM(IFAC,2,1) = -FUNCEL / VOLG
328 : 25 CONTINUE
329 : 20 CONTINUE
330 : SEGDES JTR
331 : SEGDES JTL
332 : C-----------------------------------------------------
333 : ELSE
334 : C-----------------------------------------------------
335 : C************* Murs (NGCG = NGCD) ******************
336 : C-----------------------------------------------------
337 : MELEDU.NUM(1,IFAC) = NGCG
338 : MELEDU.NUM(2,IFAC) = NGCD
339 : NLCG = MLENTC.LECT(NGCG)
340 : C--------------------------------------
341 : ROG = MPRN.VPOCHA(NLCG,1)
342 : PG = MPPN.VPOCHA(NLCG,1)
343 : UXG = MPUN.VPOCHA(NLCG,1)
344 : UYG = MPUN.VPOCHA(NLCG,2)
345 : VOLG = MPVOLU.VPOCHA(NLCG,1)
346 : C-------------------------------------------
347 : WVEC_L(1)=ROG
348 : WVEC_L(2)=UXG
349 : WVEC_L(3)=UYG
350 : WVEC_L(4)=PG
351 : C-------------------------------------------------
352 : SURF = MPOVSU.VPOCHA(NLCF,1)
353 : NVECT(1) = MPNORM.VPOCHA(NLCF,1)
354 : NVECT(2) = MPNORM.VPOCHA(NLCF,2)
355 : TVECT(1) =-NVECT(2)
356 : TVECT(2) = NVECT(1)
357 : C------- COEFFICIENTS ----------------------------
358 : ZC1=NVECT(1)*TVECT(2)+TVECT(1)*NVECT(2)
359 : ZC2=NVECT(1)*TVECT(2)-TVECT(1)*NVECT(2)
360 : ZC3=2.0D0*NVECT(1)*TVECT(1)
361 : ZC4=2.0D0*NVECT(2)*TVECT(2)
362 : C-------------------------------------------------
363 : ROD = ROG
364 : PD = PG
365 : UXD = -ZC1*UXG/ZC2-ZC4*UYG/ZC2
366 : UYD = ZC3*UXG/ZC2+ZC1*UYG/ZC2
367 : VOLD = VOLG
368 : C------------------------------------------------
369 : UPR_L=MPUPRI.VPOCHA(NLCG,1)
370 : UPR_R=UPR_L
371 : V_INF=MAX(MPUINF.VPOCHA(NLCG,1),UPR_L)
372 : C------------------------------------------------
373 : WVEC_R(1)=ROD
374 : WVEC_R(2)=UXD
375 : WVEC_R(3)=UYD
376 : WVEC_R(4)=PD
377 : C-------------------------------------------
378 : C********** La normale sortante
379 : C-------------------------------------------
379 : CALL jbmmw2(NSP,JLL,WVEC_L,WVEC_R,NVECT,
380 : & TVECT,MPYN,MLRECP,MLRECV,NLCG,V_INF)
381 : C----------------------------------------------------
382 : JTL=JLL
383 : SEGACT JTL
384 : DO 70 II = 1,(3+NSP)
385 : DO 75 JJ = 1,(3+NSP)
386 : KV = (II-1)*(3+NSP)
387 : C---------------------------------
388 : CELL = IMATRI.LIZAFM(1,KV+JJ)
389 : FUNCEL = SURF * JTL.JAC(II,JJ)
390 : CELL.AM(IFAC,1,1) = -1.0D0 * FUNCEL / VOLG
391 : CELL.AM(IFAC,1,2) = 0.0D0
392 : 75 CONTINUE
393 : 70 CONTINUE
394 : SEGDES JTL
395 : C--------------------------------------------
396 : DO 40 II = 1,(3+NSP)
397 : DO 45 JJ = 1,(3+NSP)
398 : KV = (II-1)*(3+NSP)
399 : CELL = IMATRI.LIZAFM(1,KV+JJ)
400 : CELL.AM(IFAC,2,2) = 0.0D0
401 : CELL.AM(IFAC,2,1) = 0.0D0
402 : 45 CONTINUE
403 : 40 CONTINUE
404 : C--------------------------------------------
405 : ENDIF
406 : ENDDO
407 : C
408 : SEGDES MELEMC
409 : SEGDES MELEFE
410 : SEGDES MELEMF
411 : C
412 : SEGDES MPOVSU
413 : SEGDES MPVOLU
414 : SEGDES MPNORM
415 : C
416 : SEGDES MPRN
417 : SEGDES MPPN
418 : SEGDES MPUN
419 : SEGDES MPYN
420 : C
421 : SEGDES MELEDU
422 : SEGDES MATRIK
423 : DO 80 II=1,NBME
424 : CELL = IMATRI.LIZAFM(1,II)
425 : SEGDES CELL
426 : 80 CONTINUE
427 : SEGDES IMATRI
428 : C
429 : SEGSUP MLENTC
430 : SEGSUP MLENTF
431 : SEGSUP MLELIM
432 : SEGDES MLMINC
433 : C
434 : SEGDES MPUPRI
435 : SEGDES MPUINF
436 : IF(MELLIM .NE.0) SEGDES MELLIM
437 :
438 : 9999 CONTINUE
439 : RETURN
440 : END
441 :
442 :
443 :
444 :
445 :
446 :
447 :
448 :
449 :
450 :
451 :
452 :
453 :
454 :
455 :
© Cast3M 2003 - All rights reserved.
Disclaimer