1 : C PRIMI2 SOURCE CHAT 05/01/13 02:27:10 5004 2 : SUBROUTINE PRIMI2(NESP,NORDP1,NSCA,PROPHY, 3 : & ICEN,IRO,IROVIT,IROET,IROY,IROSCA,LOGTEM,IT, 4 : & IVIT,IPRES,ITEMP,IY,ISCA,IGAM, 5 : & LOGAN,LOGNEG,LOGBOR,LOGIPG,LOGNC,MESERR, 6 : & VALER,VAL1,VAL2) 7 : C 8 : C************************************************************************ 9 : C 10 : C PROJET : CASTEM 2000 11 : C 12 : C NOM : PRIMI2 13 : C 14 : C DESCRIPTION : VOIR PRIMI1 15 : C 16 : C Melange des gaz "thermally perfect". 17 : C 18 : C Calcul de la vitesse, de la pression, de la 19 : C temperature, des fractions massiques, de "gamma" 20 : C 21 : C LANGAGE : FORTRAN 77 + ESOPE 2000 (avec extensions CISI) 22 : C 23 : C AUTEUR : A. BECCANTINI, DRN/DMT/SEMT/TTMF 24 : C 25 : C************************************************************************ 26 : C 27 : C APPELES : LICHT, ACCTAB 28 : C 29 : C************************************************************************ 30 : C 31 : C ENTREES : NESP : nombre d'especes dans le melange. 32 : C 33 : C NORDP1 : ordre des polynoms du cv_i + 1 34 : C 35 : C NSCA : nombre scalaires passif splittés 36 : C 37 : C PROPHY : pointeur sur un segment qui contient les 38 : C proprietes du melange 39 : C 40 : C ICEN : MELEME contenant les centres des ELTs 41 : C 42 : C IRO : CHPOINT contenant la masse volumique. 43 : C 44 : C IROVIT : CHPOINT contenant les dèbits 45 : C ( NDIM composantes); 46 : C 47 : C IROET : CHPOINT contenat l'énergie totale per 48 : C unité de volume (RHO Et); 49 : C 50 : C IROY : CHPOINT contenant la masse de 51 : C differentes especes; 52 : C 53 : C IROSCA : CHPOINT contenant les scalaires passifs splittes 54 : C fois la masse; 55 : C 56 : C LOGTEM : si .TRUE., la temperature de premier tentative 57 : C pour Newton-Rapson est donne; 58 : C 59 : C IT : pointeur sur ls CHPOINT de la temperature de 60 : C premier tentative; 61 : C 62 : C SORTIES : 63 : C 64 : C IVIT : CHPOINT contenant la vitesse 65 : C 66 : C IPRES : CHPOINT contenant la pression du gaz; 67 : C 68 : C ITEMP : CHPOINT contenant la temperature du 69 : C gaz; 70 : C 71 : C IY : CHPOINT contenant les fractions 72 : C massiques du gaz; 73 : C 74 : C ISCA : CHPOINT contenant les scalaires passifs splittes; 75 : C 76 : C IGAM : CHPOINT contenant les "gamma" du gaz; 77 : C 78 : C LOGAN : anomalie detectée 79 : C 80 : C LOGNEG : (LOGICAL): si .TRUE. une densité ou un temperature 81 : C negative a été detectée -> le programme s'arrete 82 : C et on peut visualiser le champe de densite ou de 83 : C temperature (sa valeur stockée en MESERR(1) et 84 : C VALER(1)) 85 : C 86 : C LOGBOR : (LOGICAL) 87 : C si .TRUE. la fraction massique a ete 88 : C detecté dehor YMIN et YMAX 89 : C (sa valeur stockée en MESERR(2) et VALER(2) 90 : C ,VAL1,VAL2) 91 : C ou 92 : C gamma a été detecté dehor GAMMIN et GAMMAX 93 : C (sa valeur stockée en MESERR(2) et 94 : C VALER(2),VAL1,VAL2) 95 : C 96 : C LOGIPG : si .TRUE., cv(T) < 0 97 : C 98 : C LOGNC : si .TRUE., Newton Rapson pour le calcul de T 99 : C n'a pas converge 100 : C 101 : C MESERR(2), 102 : C VALER(2), 103 : C VAL1, 104 : C VAL2 : pour message d'erreur 105 : C 106 : C 107 : C************************************************************************ 108 : C 109 : C HISTORIQUE (Anomalies et modifications éventuelles) 110 : C 111 : C HISTORIQUE : Créée le 16.12.98. 112 : C 113 : C************************************************************************ 114 : C 115 : C 116 : C**** Variables de COOPTIO 117 : C 118 : C INTEGER IPLLB, IERPER, IERMAX, IERR, INTERR 119 : C & ,IOTER, IOLEC, IOIMP, IOCAR, IOACQ 120 : C & ,IOPER, IOSGB, IOGRA, IOSAU, IORES 121 : C & ,IECHO, IIMPI, IOSPI 122 : C & ,IDIM 123 : C & ,MCOORD 124 : C & ,IFOMOD, NIFOUR, IFOUR, NSDPGE, IONIVE 125 : C & ,NGMAXY, IZROSF, ISOTYP, IOSCR,LTEXLU 126 : C & ,NORINC,NORVAL,NORIND,NORVAD 127 : C & ,NUCROU, IPSAUV 128 : C 129 : C**** Les variables 130 : C 131 : IMPLICIT INTEGER(I-N) 132 : INTEGER NESP,NORDP1,NSCA 133 : & ,ICEN,IRO,IROVIT,IROET,IROY,IROSCA 134 : & ,IVIT,IPRES,IY,ISCA,IGAM,IGEOMC 135 : & ,N1,NLCE,I1,ITEMP,IT,I2 136 : C 137 : C**** NESP = Nombre d'especes dans le gaz. 138 : C 139 : REAL*8 EPSI,VALER(2),VAL1,VAL2 140 : & ,YMIN,YMAX 141 : & ,RO,UX,UY,UZ,RETOT,ETHER,CELL 142 : & ,RTOT,CVTOT,T0 143 : & ,GAMMA,GAMMIN,GAMMAX 144 : & ,T 145 : CHARACTER*(8) TYPE 146 : CHARACTER*(40) MESERR(2) 147 : LOGICAL LOGNEG, LOGBOR, LOGAN, LOGNC, LOGIPG, LOGTEM 148 : C 149 : C**** Valeur minimum di Y et T(°C) 150 : C 151 : C N.B.: il doit etre le meme dans gamma.eso. 152 : C 153 : PARAMETER(EPSI=1.0D-4, 154 : & YMIN=-EPSI,YMAX=1.0D0+EPSI, 155 : & GAMMIN=1.0D0,GAMMAX=3.0D0) 156 : C 157 : C**** Les includes 158 : C 159 : -INC CCOPTIO 160 : -INC SMCHPOI 161 : POINTEUR MPORO.MPOVAL, MPOROV.MPOVAL, MPOROE.MPOVAL, 162 : & MPOROY.MPOVAL, MPOVIT.MPOVAL, MPOTEM.MPOVAL, 163 : & MPOPRE.MPOVAL, MPOY.MPOVAL, MPOGAM.MPOVAL, 164 : & MSOUPY.MSOUPO, MPOTET.MPOVAL, MPOROS.MPOVAL, 165 : & MPOSCA.MPOVAL, MSOUPS.MSOUPO 166 : -INC SMELEME 167 : C 168 : C**** Y(NESP) 169 : C 170 : SEGMENT FRAMAS 171 : REAL*8 Y(NESP) 172 : ENDSEGMENT 173 : SEGINI FRAMAS 174 : C 175 : C**** Segment du propriete du gaz 176 : C 177 : SEGMENT PROPHY 178 : REAL*8 ACV(NORD+1,NESP), R(NESP), HOK(NESP) 179 : ENDSEGMENT 180 : C 181 : C**** Initialisation des variables pour la gestion des erreurs pas ici, 182 : C mais avant 183 : C 184 : SEGMENT SUMCV 185 : REAL*8 ACVTOT(NORDP1) 186 : ENDSEGMENT 187 : SEGINI SUMCV 188 : C 189 : C**** Activation du MELEME "CENTRE" 190 : C 191 : IPT1 = ICEN 192 : SEGACT IPT1 193 : N1 = IPT1.NUM(/2) 194 : SEGDES IPT1 195 : C 196 : C*** IVIT ressemble CHPOINT IROVIT (DEBIT). 197 : C 198 : MCHPO1 = IROVIT 199 : SEGACT MCHPO1 200 : MSOUP1 = MCHPO1.IPCHP(1) 201 : SEGDES MCHPO1 202 : SEGACT MSOUP1 203 : MPOROV = MSOUP1.IPOVAL 204 : SEGDES MSOUP1 205 : SEGACT MPOROV 206 : C 207 : C**** IVIT 208 : C 209 : SEGINI, MPOVIT = MPOROV 210 : SEGINI, MSOUP2 = MSOUP1 211 : MSOUP2.IPOVAL = MPOVIT 212 : SEGINI, MCHPO2 = MCHPO1 213 : MCHPO2.IPCHP(1)= MSOUP2 214 : SEGDES MSOUP2 215 : SEGDES MCHPO2 216 : IVIT = MCHPO2 217 : C 218 : C**** Creation des CHPOINTs IPRES et IGAM 219 : C 220 : C Ce CHPOINT ressemble à IRO 221 : C Donc on lit IRO 222 : C 223 : MCHPO1 = IRO 224 : SEGACT MCHPO1 225 : MSOUP1 = MCHPO1.IPCHP(1) 226 : SEGDES MCHPO1 227 : SEGACT MSOUP1 228 : MPORO = MSOUP1.IPOVAL 229 : SEGDES MSOUP1 230 : SEGACT MPORO 231 : C 232 : C*** IPRES 233 : C 234 : SEGINI, MPOPRE = MPORO 235 : SEGINI, MSOUP2 = MSOUP1 236 : MSOUP2.IPOVAL = MPOPRE 237 : SEGINI, MCHPO2 = MCHPO1 238 : MCHPO2.IPCHP(1)= MSOUP2 239 : SEGDES MSOUP2 240 : SEGDES MCHPO2 241 : IPRES = MCHPO2 242 : C 243 : C*** ITEMP 244 : C 245 : SEGINI, MPOTEM = MPORO 246 : SEGINI, MSOUP2 = MSOUP1 247 : MSOUP2.IPOVAL = MPOTEM 248 : SEGINI, MCHPO2 = MCHPO1 249 : MCHPO2.IPCHP(1)= MSOUP2 250 : SEGDES MSOUP2 251 : SEGDES MCHPO2 252 : ITEMP = MCHPO2 253 : C 254 : C*** IGAM 255 : C 256 : SEGINI, MPOGAM = MPORO 257 : SEGINI, MSOUP2 = MSOUP1 258 : MSOUP2.IPOVAL = MPOGAM 259 : SEGINI, MCHPO2 = MCHPO1 260 : MCHPO2.IPCHP(1)= MSOUP2 261 : SEGDES MSOUP2 262 : SEGDES MCHPO2 263 : IGAM = MCHPO2 264 : IF(IROY .GT. 0)THEN 265 : C 266 : C*** IY 267 : C 268 : C Ce CHPOINT ressemble à IROY 269 : C Donc on lit IROY 270 : C 271 : MCHPO1 = IROY 272 : SEGACT MCHPO1 273 : MSOUPY = MCHPO1.IPCHP(1) 274 : SEGDES MCHPO1 275 : SEGACT MSOUPY 276 : MPOROY = MSOUPY.IPOVAL 277 : SEGACT MPOROY 278 : C 279 : SEGINI, MPOY = MPOROY 280 : SEGINI, MSOUP2 = MSOUPY 281 : MSOUP2.IPOVAL = MPOY 282 : SEGINI, MCHPO2 = MCHPO1 283 : MCHPO2.IPCHP(1)= MSOUP2 284 : SEGDES MSOUP2 285 : SEGDES MCHPO2 286 : SEGDES MSOUPY 287 : IY = MCHPO2 288 : ELSE 289 : IY=0 290 : ENDIF 291 : C 292 : C*** ISCA 293 : C 294 : C Ce CHPOINT ressemble à IROSCA 295 : C Donc on lit IROSCA 296 : C 297 : IF(IROSCA .GT. 0)THEN 298 : MCHPO1 = IROSCA 299 : SEGACT MCHPO1 300 : MSOUPS = MCHPO1.IPCHP(1) 301 : SEGDES MCHPO1 302 : SEGACT MSOUPS 303 : MPOROS = MSOUPS.IPOVAL 304 : SEGACT MPOROS 305 : C 306 : SEGINI, MPOSCA = MPOROS 307 : SEGINI, MSOUP2 = MSOUPS 308 : MSOUP2.IPOVAL = MPOSCA 309 : SEGINI, MCHPO2 = MCHPO1 310 : MCHPO2.IPCHP(1)= MSOUP2 311 : SEGDES MSOUP2 312 : SEGDES MCHPO2 313 : SEGDES MSOUPS 314 : ISCA = MCHPO2 315 : ELSE 316 : ISCA = 0 317 : ENDIF 318 : C 319 : C**** Lecture des MPOVALs des autres MCHPOIs 320 : C 320 : CALL LICHT(IROET,MPOROE,TYPE,IGEOMC) 321 : C 322 : C**** LICHT active les MPOVALs en *MOD 323 : C 324 : C i.e. 325 : C 326 : C SEGACT MPOROE*MOD 327 : C 328 : C 329 : C**** Temperature de premier tentative 330 : C 331 : IF(LOGTEM)THEN 331 : CALL LICHT(IT,MPOTET,TYPE,IGEOMC) 332 : ELSE 333 : SEGINI, MPOTET = MPOROE 334 : DO I1=1,N1 335 : MPOTET.VPOCHA(I1,1) = 600 336 : ENDDO 337 : ENDIF 338 : C 339 : C**** RICAPITOLATIF 340 : C 341 : C On a activé que les CHPOINTs suivants 342 : C 343 : C MPORO < -- > RO 344 : C MPOROV < -- > DEBITS 345 : C MPOROE < -- > ROET 346 : C MPOROY < -- > ROY (si NESP >1) 347 : C MPOROS < -- > ROS (si NSCA >1) 348 : C MPOVIT < -- > IVIT 349 : C MPOPRE < -- > IPRES 350 : C MPOTEM < -- > ITEM 351 : C MPOY < -- > IY (si NESP > 1) 352 : C MPOSCA < -- > ISCA (si NSCA > 1) 353 : C MPOGAM < -- > IGAM 354 : C MPOTET < -- > Temperature de premier tentative 355 : C 356 : C C 357 : C 358 : C**** BOUCLE SUR LES CENTRES pour le calcul des MPOVAL 359 : C 360 : DO NLCE = 1, N1 361 : C 362 : C******* Les differents variables a chaque centre 363 : C 364 : RO = MPORO.VPOCHA(NLCE,1) 365 : IF(RO .LE. 0.0D0)THEN 366 : VALER(1) = RO 367 : MESERR(1) = 'RO ' 368 : LOGNEG = .TRUE. 369 : C 370 : C********** RO < 0: le programme s'arrete mais apres le calcul des 371 : C CHPOINTs 372 : C 373 : ENDIF 374 : UX = MPOROV.VPOCHA(NLCE,1)/RO 375 : UY = MPOROV.VPOCHA(NLCE,2)/RO 376 : C 377 : C******* VITESSE -> MPOVIT 378 : C 379 : MPOVIT.VPOCHA(NLCE,1) = UX 380 : MPOVIT.VPOCHA(NLCE,2) = UY 381 : CELL = UX * UX + UY * UY 382 : IF(IDIM .EQ. 3)THEN 383 : UZ = MPOROV.VPOCHA(NLCE,3)/RO 384 : MPOVIT.VPOCHA(NLCE,3) = UZ 385 : CELL = CELL + UZ * UZ 386 : ENDIF 387 : CELL = 0.5D0 * RO * CELL 388 : C 389 : C******* RETOT, ETHER 390 : C 391 : RETOT = MPOROE.VPOCHA(NLCE,1) 392 : ETHER = RETOT - CELL 393 : ETHER = ETHER / RO 394 : C 395 : C******* Y 396 : C 397 : C On a utilisé le ponteur par defaut, 398 : C i.e. Y(I1) = FRAMAS.Y(I1) 399 : C 400 : C 401 : Y(NESP) = 1.0D0 402 : C 403 : C******* SUMCV.ACVTOT a ZERO 404 : C 405 : DO I1 = 1, NORDP1 406 : SUMCV.ACVTOT(I1) = 0.0D0 407 : ENDDO 408 : C 409 : C******* N.B.: NESP => 1 410 : C 411 : RTOT = 0.0D0 412 : IF(NESP .GT. 1)THEN 413 : DO I1 = 1, NESP-1 414 : Y(I1) = MPOROY.VPOCHA(NLCE,I1)/RO 415 : MPOY.VPOCHA(NLCE,I1)=Y(I1) 416 : Y(NESP) = Y(NESP) - Y(I1) 417 : IF((Y(I1) .LT. YMIN) .OR. (Y(I1) .GT. YMAX))THEN 418 : MESERR(2) = 'Y ' 419 : VALER(2) = Y(I1) 420 : LOGBOR = .TRUE. 421 : VAL1 = YMIN 422 : VAL2 = YMAX 423 : C 424 : C************* Y !\in (YMIN,YMAX) : le programme s'arrete 425 : C mais apres le calcul des CHPOINTs 426 : C 427 : ENDIF 428 : DO I2 = 1, NORDP1 429 : SUMCV.ACVTOT(I2) = SUMCV.ACVTOT(I2) + 430 : & (Y(I1) * PROPHY.ACV(I2,I1)) 431 : ENDDO 432 : RTOT = RTOT + Y(I1) * PROPHY.R(I1) 433 : ENDDO 434 : IF((Y(NESP) .LT. YMIN) .OR. (Y(NESP) .GT. YMAX))THEN 435 : MESERR(2) = 'Y ' 436 : VALER(2) = Y(NESP) 437 : LOGBOR = .TRUE. 438 : VAL1 = YMIN 439 : VAL2 = YMAX 440 : C 441 : C********** Y !\in (YMIN,YMAX) : le programme s'arrete 442 : C mais apres le calcul des CHPOINTs 443 : C 444 : ENDIF 445 : ENDIF 446 : DO I2 = 1, NORDP1 447 : SUMCV.ACVTOT(I2) = SUMCV.ACVTOT(I2) 448 : & + (Y(NESP) * PROPHY.ACV(I2,NESP)) 449 : ENDDO 450 : RTOT = RTOT + Y(NESP) * PROPHY.R(NESP) 451 : IF(RTOT .LE. 0.0D0)THEN 452 : LOGIPG = .TRUE. 453 : ENDIF 454 : C 455 : C******* Calcul de la temperature et de la pression 456 : C 457 : T0 = MPOTET.VPOCHA(NLCE,1) 457 : CALL PRITEM(NORDP1,SUMCV.ACVTOT,ETHER,T0,T,LOGNC,LOGIPG) 458 : MPOTEM.VPOCHA(NLCE,1) = T 459 : MPOPRE.VPOCHA(NLCE,1) = RTOT * T * RO 460 : IF(T .LE. 0.0D0)THEN 461 : VALER(1) = T 462 : MESERR(1) = 'T ' 463 : LOGNEG = .TRUE. 464 : C 465 : C********** T < 0: le programme s'arrete mais apres le calcul des 466 : C CHPOINTs 467 : C 468 : ENDIF 469 : 470 : C 471 : C******* Calcul du GAMMA 472 : C 473 : C Le CV 474 : C 475 : CVTOT = SUMCV.ACVTOT(1) 476 : T0 = 1.0D0 477 : DO I1 = 2, NORDP1 478 : T0 = T0 * T 479 : CELL = SUMCV.ACVTOT(I1) * T0 480 : CVTOT = CVTOT + CELL 481 : ENDDO 482 : GAMMA = (CVTOT + RTOT) / CVTOT 483 : IF((GAMMA .LT. GAMMIN) .OR. (GAMMA .GT. GAMMAX))THEN 484 : MESERR(2) = 'GAMMA ' 485 : VALER(2) = GAMMA 486 : VAL1 = GAMMIN 487 : VAL2 = GAMMAX 488 : LOGBOR = .TRUE. 489 : C 490 : C********** GAMMA !\in (GAMMIN,GAMMAX) : le programme s'arrete 491 : C mais apres le calcul des CHPOINTs 492 : C 493 : ENDIF 494 : MPOGAM.VPOCHA(NLCE,1) = GAMMA 495 : C 496 : C******* Calcul du vecteur de scalaire 497 : C 498 : DO I1 = 1, NSCA, 1 499 : MPOSCA.VPOCHA(NLCE,I1) = MPOROS.VPOCHA(NLCE,I1) / RO 500 : ENDDO 501 : C 502 : C**** Fin boucle sur les points 503 : C 504 : ENDDO 505 : C 506 : 9999 CONTINUE 507 : C 508 : SEGDES MPORO 509 : SEGDES MPOROV 510 : SEGDES MPOROE 511 : IF(NESP .GT. 1)THEN 512 : SEGDES MPOROY 513 : SEGDES MPOY 514 : ENDIF 515 : IF(NSCA .GE. 1)THEN 516 : SEGDES MPOROS 517 : SEGDES MPOSCA 518 : ENDIF 519 : SEGDES MPOVIT 520 : SEGDES MPOPRE 521 : SEGDES MPOGAM 522 : SEGDES MPOTEM 523 : IF(LOGTEM)THEN 524 : SEGDES MPOTET 525 : ELSE 526 : SEGSUP MPOTET 527 : ENDIF 528 : SEGSUP FRAMAS 529 : SEGSUP SUMCV 530 : C 531 : RETURN 532 : END 533 : 534 : 535 : 536 : 537 : 538 :
© Cast3M 2003 - All rights reserved.
Disclaimer