1 : C PRIDEM SOURCE BECC 09/11/18 21:15:06 6542 2 : SUBROUTINE PRIDEM() 3 : C************************************************************************ 4 : C 5 : C PROJET : CASTEM 2000 6 : C 7 : C NOM : PRIDEM 8 : C 9 : C DESCRIPTION : Voir PRIMIT 10 : C 11 : C RDEM approach for combustion. 12 : C Computation of the primitive variables. 13 : C 14 : C LANGAGE : FORTRAN 77 + ESOPE 2000 (avec extensions CISI) 15 : C 16 : C AUTEUR : A. BECCANTINI, DEN/DM2S/SFME/LTMF 17 : C 18 : C************************************************************************ 19 : C 20 : C 21 : C APPELES (Calcul) : PRIDE1 22 : C 23 : C 24 : C************************************************************************ 25 : C 26 : C 26 : C CALL (GIBIANE) : 27 : C 28 : C RCHV1 RCHV2 RCHP1 RCHP2 RCHT1 RCHT2 = 'PRIM' 'DEM' TABPGAS 29 : C CHPAL1 CHPAL2 CHPARN1 CHPARN2 CHPAGN1 CHPAGN2 30 : C CHPARET1 CHPARET2 CHPTGUE1 CHPTGUE2 EPS ; 31 : C 32 : C 33 : C ENTREES : 34 : C 35 : C 36 : C TABPGAS : TABLE which contains 37 : C * 'SPECIES' 38 : C * 'CHEM_COEF' 39 : C * 'MASSFRA' initial and final mass fraction of 40 : C the first appearing in 'SPECIES', 41 : C final mass fractions of the species with 42 : C positive coefficients in 'CHEM_COEF', 43 : C initial mass fractions for the species with 44 : C negative coefficients in 'CHEM_COEF' 45 : C * 'RUNIV' = universal gas constant, 46 : C * ESPi = table containing the properties of 47 : C the species ESPi 48 : C * 'TMAX' maximum temperature for cv expansion; 49 : C for T>'TMAX', cv(T)=cv('TMAX') 50 : C * ESPI . 'A' 51 : C CV_i = \sum_{j=0,k} A_{i,j} T^j 52 : C * ESPI . 'W' (Kg/mole) 53 : C * ESPI . 'H0K' 54 : C e_{0,i} = h_{0,i} = h_{T_0,i} - {R_i * T_0 + 55 : C {\sum_{j=0,k} A_{i,j} / (j+1) T_0^(j+1)}}; 56 : C 57 : C CHPAL1 : CHPOINT which contains the volume fraction alpha_1 58 : C of 1 (one component, 'SCAL'). 59 : C 60 : C CHPAL2 : CHPOINT which contains the volume fraction alpha_2 61 : C of 2 (one component, 'SCAL'). 62 : C 63 : C CHPARN1 : CHPOINT which contains the alpha_1 * density of 1 64 : C (one component, 'SCAL'). 65 : C 66 : C CHPARN2 : CHPOINT which contains the alpha_2 * density of 2 67 : C (one component, 'SCAL'). 68 : C 69 : C CHPAGN1 : CHPOINT which contains the alpha_1 * momentom of 1 70 : C (one component, 'SCAL'). 71 : C 72 : C CHPAGN2 : CHPOINT which contains the alpha_2 * momentom of 2 73 : C (one component, 'SCAL'). 74 : C 75 : C CHPARET1: CHPOINT which contains the alpha_1 * total energy 76 : C of 1 (one component, 'SCAL'). 77 : C 78 : C CHPARET2: CHPOINT which contains the alpha_2 * total energy 79 : C of 2 (one component, 'SCAL'). 80 : C 81 : C CHPTGUE1: CHPOINT which contains the guess value 82 : C for the temperature of 1 (one component, 'SCAL'). 83 : C 84 : C CHPTGUE2: CHPOINT which contains the guess value 85 : C for the temperature of 2 (one component, 'SCAL'). 86 : C 87 : C K0 : FLOTTANT which contains the fundamental flame 88 : C speed 89 : C 90 : C EPS : FLOTTANT such that if ALPHA_i < EPS, we can say 91 : C that species i does not exists 92 : C 93 : C 94 : C SORTIES : 95 : C 96 : C RCHV1 : CHPOINT which contains the speed of 1 97 : C 98 : C RCHV2 : CHPOINT which contains the speed of 2 99 : C 100 : C RCHP1 : CHPOINT which contains the pressure of 1 101 : C 102 : C RCHP2 : CHPOINT which contains the pressure of 2 103 : C 104 : C RCHT1 : CHPOINT which contains the temperature of 1 105 : C 106 : C RCHT2 : CHPOINT which contains the temperature of 2 107 : C 108 : C 109 : C************************************************************************ 110 : C 111 : C HISTORIQUE : Crée le 06.09.09. 112 : C 113 : C************************************************************************ 114 : C 115 : C 116 : C**** Les variables 117 : C 118 : IMPLICIT INTEGER(I-N) 119 : INTEGER ICOND, IRETOU, INDIC, NBCOMP 120 : & , NESP, NESP1, ICEN, IALP1, IALP2 121 : & , IARN1, IARN2 122 : & , IAGN1, IAGN2, IARET1, IARET2, ITG1, ITG2 123 : & , IPGAS, IESP 124 : & , IVN1, IVN2, IPN1, IPN2, IRN1, IRN2, ITN1, ITN2 125 : & , I1, I2, JGM, JGN, NORD, NORDP1, NORD1 126 : REAL*8 VALER(2),VAL1,VAL2, TMAX, RUNIV, EPS 127 : CHARACTER*(40) MESERR(2),MESCEL 128 : CHARACTER*(8) TYPE 129 : CHARACTER*(4) MOT1(1) 130 : LOGICAL LOGNEG, LOGBOR, LOGAN, LOGNC, LOGIPG 131 : C 132 : C**** Variables en ACCTAB 133 : C 134 : INTEGER IVALI, IRETI,IVALR, IRETR 135 : REAL*8 XVALI, XVALR 136 : LOGICAL LOGII, LOGIR 137 : CHARACTER*(8) CHARR,MTYPI 138 : C 139 : C**** Segment des proprietes du gaz 140 : C 141 : SEGMENT PROPHY 142 : REAL*8 ACV(NORD+1,NESP), W(NESP), H0K(NESP) 143 : ENDSEGMENT 144 : C 145 : C**** Les Includes 146 : C 147 : -INC CCOPTIO 148 : -INC SMCHPOI 149 : -INC SMLMOTS 150 : -INC SMLREEL 151 : POINTEUR MLMOSC.MLMOTS, MLMESP.MLMOTS 152 : POINTEUR MLRMFR.MLREEL, MLRCHE.MLREEL 153 : C 154 : C**** Initialisation des parametres d'erreur 155 : C 156 : LOGAN = .FALSE. 157 : LOGNEG = .FALSE. 158 : LOGBOR = .FALSE. 159 : LOGNC = .FALSE. 160 : LOGIPG = .FALSE. 161 : MESCEL = ' ' 162 : MESERR(1) = MESCEL 163 : MESERR(2) = MESCEL 164 : MOTERR(1:40) = MESCEL(1:40) 165 : VALER(1) = 0.0D0 166 : VALER(2) = 0.0D0 167 : VAL1 = 0.0D0 168 : VAL2 = 0.0D0 169 : C 170 : C**** Initialisation des variables en ACCTAB 171 : C 172 : IVALI = 0 173 : IVALR = 0 174 : XVALI = 0.0D0 175 : XVALR = 0.0D0 176 : LOGII = .FALSE. 177 : LOGIR = .FALSE. 178 : IRETI = 0 179 : IRETR = 0 180 : CHARR = ' ' 181 : C 182 : C**** Initialisation des MOT1(1) (noms des composantes) 183 : C 184 : MOT1(1) = ' ' 185 : C 186 : C**** N.B. On veut lire les objets sequentiellement. 187 : C Donc on utilise QUETYP pour controler que 188 : C le type de l'objet soit le bon. 189 : C 190 : C**** Lecture de la table des proprietes du gaz 191 : C 192 : ICOND = 1 192 : CALL QUETYP(TYPE,ICOND,IRETOU) 193 : IF(IERR .NE. 0)GOTO 9999 194 : IF(TYPE .NE. 'TABLE ')THEN 195 : C 196 : C******* Message d'erreur standard 197 : C 37 2 198 : C On ne trouve pas d'objet de type %m1:8 199 : C 200 : MOTERR(1:8) = 'TABLE ' 200 : CALL ERREUR(37) 201 : GOTO 9999 202 : ELSE 203 : ICOND = 1 203 : CALL LIROBJ(TYPE,IPGAS,ICOND,IRETOU) 204 : IF(IERR .NE. 0)GOTO 9999 205 : ENDIF 206 : C 207 : C**** Ordre des polynoms pour les cv_i 208 : C 209 : MTYPI = 'MOT ' 210 : TYPE = ' ' 210 : CALL ACCTAB(IPGAS,MTYPI,IVALI,XVALI,'NORD',LOGII,IRETI, 211 : & TYPE,NORD,XVALR,CHARR,LOGIR,IESP) 212 : IF(TYPE .NE. 'ENTIER ')THEN 213 : C 214 : C******* Message d'erreur standard 215 : C -301 0 %m1:40 216 : C 217 : MOTERR(1:40) = 'TAB1 . NORD = ??? ' 218 : WRITE(IOIMP,*) MOTERR(1:40) 219 : C 220 : C******* Message d'erreur standard 221 : C 21 2 222 : C Données incompatibles 223 : C 223 : CALL ERREUR(21) 224 : GOTO 9999 225 : ENDIF 226 : NORDP1 = NORD + 1 227 : C 228 : C**** 'TMAX' 229 : C 230 : MTYPI = 'MOT ' 231 : TYPE = ' ' 231 : CALL ACCTAB(IPGAS,MTYPI,IVALI,XVALI,'TMAX',LOGII,IRETI, 232 : & TYPE,IVALR,XVALR,CHARR,LOGIR,IESP) 233 : IF(TYPE .NE. 'FLOTTANT')THEN 234 : C 235 : C******* Message d'erreur standard 236 : C -301 0 %m1:40 237 : C 238 : MOTERR(1:40) = 'TAB1 . TMAX = ??? ' 239 : WRITE(IOIMP,*) MOTERR(1:40) 240 : C 241 : C******* Message d'erreur standard 242 : C 21 2 243 : C Données incompatibles 244 : C 244 : CALL ERREUR(21) 245 : GOTO 9999 246 : ENDIF 247 : TMAX = XVALR 248 : C 249 : C**** 'RUNIV' 250 : C 251 : MTYPI = 'MOT ' 252 : TYPE = ' ' 252 : CALL ACCTAB(IPGAS,MTYPI,IVALI,XVALI,'RUNIV',LOGII,IRETI, 253 : & TYPE,IVALR,XVALR,CHARR,LOGIR,IESP) 254 : IF(TYPE .NE. 'FLOTTANT')THEN 255 : C 256 : C******* Message d'erreur standard 257 : C -301 0 %m1:40 258 : C 259 : MOTERR(1:40) = 'TAB1 . RUNIV = ??? ' 260 : WRITE(IOIMP,*) MOTERR(1:40) 261 : C 262 : C******* Message d'erreur standard 263 : C 21 2 264 : C Données incompatibles 265 : C 265 : CALL ERREUR(21) 266 : GOTO 9999 267 : ENDIF 268 : RUNIV = XVALR 269 : C 270 : C**** Les especes 271 : C 272 : TYPE = ' ' 272 : CALL ACMO(IPGAS,'SPECIES',TYPE,MLMESP) 273 : IF(TYPE .NE. 'LISTMOTS')THEN 274 : C 275 : C******* Message d'erreur standard 276 : C -301 0 %m1:40 277 : C 278 : MOTERR(1:40) = 'TAB1 . SPECIES = ??? ' 279 : WRITE(IOIMP,*) MOTERR(1:40) 280 : C 281 : C******* Message d'erreur standard 282 : C 21 2 283 : C Données incompatibles 284 : C 284 : CALL ERREUR(21) 285 : GOTO 9999 286 : ELSE 287 : SEGACT MLMESP 288 : NESP = MLMESP.MOTS(/2) 289 : SEGDES MLMESP 290 : ENDIF 291 : C 292 : C**** 'MASSFRA' 293 : C 294 : TYPE = ' ' 294 : CALL ACMO(IPGAS,'MASSFRA',TYPE,MLRMFR) 295 : IF(TYPE .NE. 'LISTREEL')THEN 296 : C 297 : C******* Message d'erreur standard 298 : C -301 0 %m1:40 299 : C 300 : MOTERR(1:40) = 'TAB1 . MASSFRA = ??? ' 301 : WRITE(IOIMP,*) MOTERR(1:40) 302 : C 303 : C******* Message d'erreur standard 304 : C 21 2 305 : C Données incompatibles 306 : C 306 : CALL ERREUR(21) 307 : GOTO 9999 308 : ELSE 309 : SEGACT MLRMFR 310 : NESP1 = MLRMFR.PROG(/1) 311 : IF (NESP1 .NE. NESP) THEN 312 : MOTERR(1:40) = 'TAB1 . MASSFRA = ??? ' 313 : WRITE(IOIMP,*) MOTERR(1:40) 314 : MOTERR(1:40) = 'TAB1 . SPECIES = ??? ' 315 : WRITE(IOIMP,*) MOTERR(1:40) 316 : C 317 : C******* Message d'erreur standard 318 : C 21 2 319 : C Données incompatibles 320 : C 320 : CALL ERREUR(21) 321 : GOTO 9999 322 : ENDIF 323 : SEGDES MLRMFR 324 : ENDIF 325 : C 326 : C**** 'CHEMCOEF' 327 : C 328 : TYPE = ' ' 328 : CALL ACMO(IPGAS,'CHEMCOEF',TYPE,MLRCHE) 329 : IF(TYPE .NE. 'LISTREEL')THEN 330 : C 331 : C******* Message d'erreur standard 332 : C -301 0 %m1:40 333 : C 334 : write(*,*) TYPE 335 : MOTERR(1:40) = 'TAB1 . CHEMCOEF = ??? ' 336 : WRITE(IOIMP,*) MOTERR(1:40) 337 : C 338 : C******* Message d'erreur standard 339 : C 21 2 340 : C Données incompatibles 341 : C 341 : CALL ERREUR(21) 342 : GOTO 9999 343 : ELSE 344 : SEGACT MLRCHE 345 : NESP1 = MLRCHE.PROG(/1) 346 : IF (NESP1 .NE. NESP) THEN 347 : MOTERR(1:40) = 'TAB1 . CHEMCOEF = ??? ' 348 : WRITE(IOIMP,*) MOTERR(1:40) 349 : MOTERR(1:40) = 'TAB1 . SPECIES = ??? ' 350 : WRITE(IOIMP,*) MOTERR(1:40) 351 : C 352 : C******* Message d'erreur standard 353 : C 21 2 354 : C Données incompatibles 355 : C 355 : CALL ERREUR(21) 356 : GOTO 9999 357 : ENDIF 358 : SEGDES MLRCHE 359 : ENDIF 360 : C 361 : C**** On rempli les segment PROPHY 362 : C Ordre: IPGAS . 'SPECIES' 363 : C 364 : SEGINI PROPHY 365 : SEGACT MLMESP 366 : C 367 : C**** N.B. MOT1 est un CHARACTER*(4) 368 : C 369 : DO I1 = 1, NESP 370 : MOT1(1) = MLMESP.MOTS(I1) 371 : C 371 : C******* CALL ACMF(...) ne marche pas parce que on a 372 : C des blanches dans nos composantes 373 : C 374 : MTYPI = 'MOT ' 375 : TYPE = ' ' 375 : CALL ACCTAB(IPGAS,MTYPI,IVALI,XVALI,MOT1(1), LOGII,IRETI, 376 : & TYPE,IVALR,XVALR,CHARR,LOGIR,IESP) 377 : C 378 : C******* En IESP a la table IPGAS.MOT1(1) 379 : C 380 : IF((IERR .NE. 0) .OR. (TYPE .NE. 'TABLE ')) THEN 381 : 382 : C 383 : C********** Message d'erreur standard 384 : C -301 0 %m1:40 385 : C 386 : MOTERR = ' ' 387 : MOTERR(1:7) = 'TAB1 . ' 388 : MOTERR(8:11) = MOT1(1) 389 : MOTERR(13:17) = '= ???' 390 : WRITE(IOIMP,*) MOTERR(1:40) 391 : C 392 : C********** Message d'erreur standard 393 : C 21 2 394 : C Données incompatibles 395 : C 395 : CALL ERREUR(21) 396 : GOTO 9999 397 : ENDIF 398 : C 399 : C******* W 400 : C 401 : MTYPI = 'MOT ' 402 : TYPE = ' ' 402 : CALL ACCTAB(IESP,MTYPI,IVALI,XVALI, 'W' , LOGII,IRETI, 403 : & TYPE,IVALR, XVALR ,CHARR,LOGIR,IRETR) 404 : IF((IERR .NE. 0) .OR. (TYPE .NE. 'FLOTTANT')) THEN 405 : C 406 : C********** Message d'erreur standard 407 : C -301 0 %m1:40 408 : C 409 : MOTERR = ' ' 410 : MOTERR(1:7) = 'TAB1 . ' 411 : MOTERR(8:11) = MOT1(1) 412 : MOTERR(13:23) = ' . W = ??? ' 413 : WRITE(IOIMP,*) MOTERR(1:40) 414 : C 415 : C********** Message d'erreur standard 416 : C 21 2 417 : C Données incompatibles 418 : C 418 : CALL ERREUR(21) 419 : GOTO 9999 420 : ENDIF 421 : PROPHY.W(I1)=XVALR 422 : C 423 : C******* H0K 424 : C 425 : MTYPI = 'MOT ' 426 : TYPE = ' ' 426 : CALL ACCTAB(IESP,MTYPI,IVALI,XVALI, 'H0K' , LOGII,IRETI, 427 : & TYPE,IVALR, XVALR ,CHARR,LOGIR,IRETR) 428 : IF((IERR .NE. 0) .OR. (TYPE .NE. 'FLOTTANT')) THEN 429 : 430 : C 431 : C********** Message d'erreur standard 432 : C -301 0 %m1:40 433 : C 434 : MOTERR = ' ' 435 : MOTERR(1:7) = 'TAB1 . ' 436 : MOTERR(8:11) = MOT1(1) 437 : MOTERR(13:25) = ' . H0K = ??? ' 438 : WRITE(IOIMP,*) MOTERR(1:40) 439 : C 440 : C********** Message d'erreur standard 441 : C 21 2 442 : C Données incompatibles 443 : C 443 : CALL ERREUR(21) 444 : GOTO 9999 445 : ENDIF 446 : PROPHY.H0K(I1)=XVALR 447 : C 448 : C******* A 449 : C 450 : MTYPI = 'MOT ' 451 : TYPE = ' ' 451 : CALL ACCTAB(IESP,MTYPI,IVALI,XVALI, 'A' , LOGII,IRETI, 452 : & TYPE,IVALR, XVALR ,CHARR,LOGIR,IRETR) 453 : IF((IERR .NE. 0) .OR. (TYPE .NE. 'LISTREEL')) THEN 454 : 455 : C 456 : C********** Message d'erreur standard 457 : C -301 0 %m1:40 458 : C 459 : MOTERR = ' ' 460 : MOTERR(1:7) = 'TAB1 . ' 461 : MOTERR(8:11) = MOT1(1) 462 : MOTERR(13:23) = ' . A = ??? ' 463 : WRITE(IOIMP,*) MOTERR(1:40) 464 : C 465 : C********** Message d'erreur standard 466 : C 21 2 467 : C Données incompatibles 468 : C 468 : CALL ERREUR(21) 469 : GOTO 9999 470 : ENDIF 471 : MLREEL = IRETR 472 : SEGACT MLREEL 473 : NORD1 = MLREEL.PROG(/1) 474 : IF(NORD1 .NE. NORDP1)THEN 475 : C 476 : C********** Message d'erreur standard 477 : C -301 0 %m1:40 478 : C 479 : MOTERR = ' ' 480 : MOTERR(1:10) = 'DIME(TAB1.' 481 : MOTERR(11:14) = MOT1(1) 482 : MOTERR(15:37) = '.A) != (TAB1.NORD) + 1' 483 : WRITE(IOIMP,*) MOTERR(1:40) 484 : C 485 : C********** Message d'erreur standard 486 : C 21 2 487 : C Données incompatibles 488 : C 488 : CALL ERREUR(21) 489 : GOTO 9999 490 : ENDIF 491 : 492 : C 493 : C******* Dans le calcul, c'est plus utile ACV dans la forme 494 : C ACV(exponente,espece) 495 : C 496 : DO I2 = 1, NORDP1 497 : PROPHY.ACV(I2,I1)= MLREEL.PROG(I2) 498 : ENDDO 499 : SEGDES MLREEL 500 : ENDDO 501 : SEGDES MLMESP 502 : C 503 : C**** La table IPGAS donc a ete controllee et PROPHY est rempli 504 : C 505 : C 506 : C**** Lecture du CHPOINT ALPHA1 507 : C 508 : TYPE='CHPOINT ' 509 : ICOND = 1 509 : CALL LIROBJ(TYPE,IALP1,ICOND,IRETOU) 510 : IF(IERR .NE. 0)GOTO 9999 511 : C 512 : C**** On cherche le pointeur de son maillage et on l'impose sur les 513 : C autres CHPOINTs 514 : C 515 : MCHPOI = IALP1 516 : SEGACT MCHPOI 517 : MSOUPO = MCHPOI.IPCHP(1) 518 : SEGACT MSOUPO 519 : ICEN = MSOUPO.IGEOC 520 : SEGDES MSOUPO 521 : SEGDES MCHPOI 522 : C 523 : C**** Control du CHPOINT: QUEPOI 524 : C 525 : C INDIC = 1 -> on impose le pointeur du support geometrique (ICEN) 526 : C N.B. Le CHPOINT peut changer de structure pour 527 : C avoir SPG = ICEN!!!! 528 : C INDIC = 0 -> on ne fait que verifier le support geometrique 529 : C (ICEN). Si le SPG sont differents INDIC = -4 en sortie 530 : C 531 : C NBCOMP > 0 -> numero des composantes 532 : C 533 : C MOT1(1) = ' ' obligatoire s'on connais pas les noms des composantes 534 : C 535 : INDIC = 1 536 : NBCOMP = 1 537 : MOT1(1) = 'SCAL' 537 : CALL QUEPOI(IALP1, ICEN, INDIC, NBCOMP, MOT1) 538 : IF(IERR .NE. 0)THEN 539 : C 540 : C******** Message d'erreur standard 541 : C -301 0 %m1:40 542 : C 543 : MOTERR = 'IALP1 = ??? ' 544 : WRITE(IOIMP,*) MOTERR(1:40) 545 : 546 : GOTO 9999 547 : ENDIF 548 : C 549 : C**** Lecture du CHPOINT ALPHA2 550 : C 551 : TYPE='CHPOINT ' 552 : ICOND = 1 552 : CALL LIROBJ(TYPE,IALP2,ICOND,IRETOU) 553 : C 554 : C**** Control du CHPOINT: QUEPOI 555 : C 556 : INDIC = 1 557 : NBCOMP = 1 558 : MOT1(1) = 'SCAL' 558 : CALL QUEPOI(IALP2, ICEN, INDIC, NBCOMP, MOT1) 559 : IF(IERR .NE. 0)THEN 560 : C 561 : C******** Message d'erreur standard 562 : C -301 0 %m1:40 563 : C 564 : MOTERR = 'IALP2 = ??? ' 565 : WRITE(IOIMP,*) MOTERR(1:40) 566 : 567 : GOTO 9999 568 : ENDIF 569 : C 570 : C**** Lecture du CHPOINT IARN1 571 : C 572 : ICOND = 1 573 : TYPE='CHPOINT ' 573 : CALL LIROBJ(TYPE,IARN1,ICOND,IRETOU) 574 : C 575 : C**** Control du CHPOINT: QUEPOI 576 : C 577 : INDIC = 1 578 : NBCOMP = 1 579 : MOT1(1) = 'SCAL' 579 : CALL QUEPOI(IARN1, ICEN, INDIC, NBCOMP, MOT1) 580 : IF(IERR .NE. 0)THEN 581 : C 582 : C******** Message d'erreur standard 583 : C -301 0 %m1:40 584 : C 585 : MOTERR = 'IARN1 = ??? ' 586 : WRITE(IOIMP,*) MOTERR(1:40) 587 : 588 : GOTO 9999 589 : ENDIF 590 : C 591 : C**** Lecture du CHPOINT IARN2 592 : C 593 : ICOND = 1 594 : TYPE='CHPOINT ' 594 : CALL LIROBJ(TYPE,IARN2,ICOND,IRETOU) 595 : IF(IERR .NE. 0)GOTO 9999 596 : C 597 : C**** Control du CHPOINT: QUEPOI 598 : C 599 : INDIC = 1 600 : NBCOMP = 1 601 : MOT1(1) = 'SCAL' 601 : CALL QUEPOI(IARN2, ICEN, INDIC, NBCOMP, MOT1) 602 : IF(IERR .NE. 0)THEN 603 : C 604 : C******** Message d'erreur standard 605 : C -301 0 %m1:40 606 : C 607 : MOTERR = 'IARN2 = ??? ' 608 : WRITE(IOIMP,*) MOTERR(1:40) 609 : 610 : GOTO 9999 611 : ENDIF 612 : C 613 : C**** Lecture du CHPOINT IAGN1 ( debits) 614 : C 615 : TYPE='CHPOINT ' 616 : ICOND = 1 616 : CALL LIROBJ(TYPE,IAGN1,ICOND,IRETOU) 617 : IF(IERR .NE. 0)GOTO 9999 618 : C 619 : C**** Control du CHPOINT 620 : C 621 : JGN = 4 622 : JGM = IDIM 623 : SEGINI MLMOTS 624 : MLMOTS.MOTS(1) = 'UX ' 625 : MLMOTS.MOTS(2) = 'UY ' 626 : IF(IDIM .EQ. 3) MLMOTS.MOTS(3) = 'UZ ' 627 : C 628 : C**** On controlle l'ordre de composantes de IAGN1 629 : C 629 : CALL QUEPO1(IAGN1, ICEN, MLMOTS) 630 : IF(IERR .NE. 0)THEN 631 : C 632 : C******** Message d'erreur standard 633 : C -301 0 %m1:40 634 : C 635 : MOTERR = 'IAGN1 = ??? ' 636 : WRITE(IOIMP,*) MOTERR(1:40) 637 : 638 : GOTO 9999 639 : ENDIF 640 : C 641 : C**** Lecture du CHPOINT IAGN2 ( debits) 642 : C 643 : ICOND = 1 644 : TYPE='CHPOINT ' 644 : CALL LIROBJ(TYPE,IAGN2,ICOND,IRETOU) 645 : IF(IERR .NE. 0)GOTO 9999 646 : C 647 : C**** Control du CHPOINT 648 : C 649 : JGN = 4 650 : JGM = IDIM 651 : SEGINI MLMOTS 652 : MLMOTS.MOTS(1) = 'UX ' 653 : MLMOTS.MOTS(2) = 'UY ' 654 : IF(IDIM .EQ. 3) MLMOTS.MOTS(3) = 'UZ ' 655 : C 656 : C**** On controlle l'ordre de composantes de IAGN2 657 : C 657 : CALL QUEPO1(IAGN2, ICEN, MLMOTS) 658 : IF(IERR .NE. 0)THEN 659 : C 660 : C******** Message d'erreur standard 661 : C -301 0 %m1:40 662 : C 663 : MOTERR = 'IAGN2 = ??? ' 664 : WRITE(IOIMP,*) MOTERR(1:40) 665 : 666 : GOTO 9999 667 : ENDIF 668 : C 669 : C**** Lecture du CHPOINT IARET1 670 : C 671 : ICOND = 1 672 : TYPE='CHPOINT ' 672 : CALL LIROBJ(TYPE,IARET1,ICOND,IRETOU) 673 : IF(IERR .NE. 0)GOTO 9999 674 : C 675 : C**** Control du CHPOINT: QUEPOI 676 : C 677 : INDIC = 1 678 : NBCOMP = 1 679 : MOT1(1) = 'SCAL' 679 : CALL QUEPOI(IARET1, ICEN, INDIC, NBCOMP, MOT1) 680 : IF(IERR .NE. 0)THEN 681 : C 682 : C******** Message d'erreur standard 683 : C -301 0 %m1:40 684 : C 685 : MOTERR = 'IARET1 = ??? ' 686 : WRITE(IOIMP,*) MOTERR(1:40) 687 : 688 : GOTO 9999 689 : ENDIF 690 : C 691 : C**** Lecture du CHPOINT IARET2 692 : C 693 : ICOND = 1 694 : TYPE='CHPOINT ' 694 : CALL LIROBJ(TYPE,IARET2,ICOND,IRETOU) 695 : IF(IERR .NE. 0)GOTO 9999 696 : C 697 : C**** Control du CHPOINT: QUEPOI 698 : C 699 : INDIC = 1 700 : NBCOMP = 1 701 : MOT1(1) = 'SCAL' 701 : CALL QUEPOI(IARET2, ICEN, INDIC, NBCOMP, MOT1) 702 : IF(IERR .NE. 0)THEN 703 : C 704 : C******** Message d'erreur standard 705 : C -301 0 %m1:40 706 : C 707 : MOTERR = 'IARET2 = ??? ' 708 : WRITE(IOIMP,*) MOTERR(1:40) 709 : 710 : GOTO 9999 711 : ENDIF 712 : C 713 : C**** Lecture du CHPOINT Tguess ITG1 714 : C 715 : ICOND = 1 716 : TYPE = 'CHPOINT ' 716 : CALL LIROBJ(TYPE,ITG1,ICOND,IRETOU) 717 : IF(IERR .NE. 0)GOTO 9999 718 : C 719 : C****** Control du CHPOINT 720 : C 721 : INDIC = 1 722 : NBCOMP = 1 723 : MOT1(1) = 'SCAL' 723 : CALL QUEPOI(ITG1, ICEN, INDIC, NBCOMP, MOT1) 724 : IF(IERR .NE. 0)THEN 725 : C 726 : C******* Message d'erreur standard 727 : C -301 0 %m1:40 728 : C 729 : MOTERR = 'ITG1 = ??? ' 730 : WRITE(IOIMP,*) MOTERR(1:40) 731 : GOTO 9999 732 : ENDIF 733 : C 734 : C**** Lecture du CHPOINT Tguess ITG2 735 : C 736 : ICOND = 1 737 : TYPE='CHPOINT ' 737 : CALL LIROBJ(TYPE,ITG2,ICOND,IRETOU) 738 : IF(IERR .NE. 0)GOTO 9999 739 : C 740 : C**** Control du CHPOINT: QUEPOI 741 : C 742 : INDIC = 1 743 : NBCOMP = 1 744 : MOT1(1) = 'SCAL' 744 : CALL QUEPOI(ITG2, ICEN, INDIC, NBCOMP, MOT1) 745 : IF(IERR .NE. 0)THEN 746 : C 747 : C******** Message d'erreur standard 748 : C -301 0 %m1:40 749 : C 750 : MOTERR = 'ITG2 = ??? ' 751 : WRITE(IOIMP,*) MOTERR(1:40) 752 : GOTO 9999 753 : ENDIF 754 : C 755 : ICOND = 1 755 : CALL LIRREE(EPS, ICOND, IRETOU) 756 : IF(IERR .NE. 0)GOTO 9999 757 : C 758 : C**** Creation of the CHAMPOIN for the results 759 : C 760 : TYPE = 'CHPOINT ' 760 : CALL KRCHP1(TYPE, ICEN, IVN1, MLMOTS) 760 : CALL KRCHP1(TYPE, ICEN, IVN2, MLMOTS) 761 : SEGSUP MLMOTS 762 : JGN = 4 763 : JGM = 1 764 : SEGINI MLMOTS 765 : MLMOTS.MOTS(1) = 'SCAL' 765 : CALL KRCHP1(TYPE, ICEN, IPN1, MLMOTS) 765 : CALL KRCHP1(TYPE, ICEN, IPN2, MLMOTS) 765 : CALL KRCHP1(TYPE, ICEN, IRN1, MLMOTS) 765 : CALL KRCHP1(TYPE, ICEN, IRN2, MLMOTS) 765 : CALL KRCHP1(TYPE, ICEN, ITN1, MLMOTS) 765 : CALL KRCHP1(TYPE, ICEN, ITN2, MLMOTS) 766 : SEGSUP MLMOTS 767 : C 768 : C**** Computation of the primitive variables 769 : C 769 : CALL PRIDE1(NESP,NORD,TMAX,RUNIV,PROPHY, 770 : & MLRCHE,MLRMFR, 771 : & ICEN,IALP1,IALP2,IARN1,IARN2,IAGN1,IAGN2,IARET1,IARET2, 772 : & ITG1,ITG2,IRN1,IRN2, 773 : & IVN1,IVN2,IPN1,IPN2,ITN1,ITN2, 774 : & EPS, 775 : & LOGAN,LOGIPG,LOGNEG,LOGBOR,LOGNC, 776 : & VALER,VAL1,VAL2) 777 : C 778 : IF(IERR .NE. 0)THEN 779 : WRITE(IOIMP,*) 'subroutine pride1' 780 : GOTO 9999 781 : ENDIF 782 : CCCC 783 : CCCC**** Calcul des sorties. 784 : CCCC 785 : CCCC Jusque a la NESP = nombre d'especes qui apparessent 786 : CCCC dans les equations d'Euler 787 : CCCC 788 : CCCC Maintenant NESP = nombre total d'espece 789 : CCCC 790 : CCC NESP = NESP + 1 790 : CCC CALL PRIMI2(NESP,NORDP1,NSCA,PROPHY, 791 : CCC & ICEN,IRO,IROVIT,IROET,IROY,IROSCA,LOGTEM,IT, 792 : CCC & IVIT,IPRES,ITEMP,IY,ISCA,IGAM, 793 : CCC & LOGAN,LOGNEG,LOGBOR,LOGIPG,LOGNC,MESERR, 794 : CCC & VALER,VAL1,VAL2) 795 : CCCC 796 : IF(LOGAN)THEN 797 : C 798 : C******* Message d'erreur standard 799 : C 5 3 800 : C Erreur anormale.contactez votre support 801 : C 801 : CALL ERREUR(5) 802 : GOTO 9999 803 : ELSE 804 : IF(LOGIPG)THEN 805 : C 806 : C********** CV(T) < 0 807 : C 808 : C 809 : C********** Message d'erreur standard 810 : C -301 0 %m1:40 811 : C 812 : MOTERR(1:40) = 'cv(T) < 0 ? R < 0 ? ' 813 : WRITE(IOIMP,*) MOTERR(1:40) 814 : MOTERR(1:40) = 'TAB1 = ??? ' 815 : WRITE(IOIMP,*) MOTERR(1:40) 816 : C 817 : C********** Message d'erreur standard 818 : C 21 2 819 : C Données incompatibles 820 : C 820 : CALL ERREUR(21) 821 : C IF(LOGTRI)THEN 822 : C IERR = 0 823 : C ELSE 824 : GOTO 9999 825 : C ENDIF 826 : ENDIF 827 : IF(LOGNC)THEN 828 : C 829 : C********** Newton - Raphson ne converge pas !!! 830 : C 831 : C 832 : C********** Message d'erreur standard 833 : C -301 0 %m1:40 834 : C 835 : MOTERR(1:40) = 'Newton - Raphson ' 836 : WRITE(IOIMP,*) MOTERR(1:40) 837 : C 838 : C********** Message d'erreur standard 839 : C 460 2 840 : C Pas de convergence dans les itérations internes 841 : C 841 : CALL ERREUR(460) 842 : GOTO 9999 843 : ENDIF 844 : IF(LOGNEG)THEN 845 : C 846 : C******* Pression (energie thermique) ou densité negative 847 : C 848 : C 849 : C******* Message d'erreur standard 850 : C 41 2 851 : C %m1:8 = %r1 inférieur à %r2 852 : C 853 : MESCEL = MESERR(1) 854 : MOTERR(1:8) = MESCEL(1:8) 855 : REAERR(1) = REAL(VALER(1)) 856 : REAERR(2) = 0.0 856 : CALL ERREUR(41) 857 : GOTO 9999 858 : ENDIF 859 : IF(LOGBOR)THEN 860 : C 861 : C ou Y !\in YMIN,YMAX 862 : C 863 : C******* Message d'erreur standard 864 : C 42 2 865 : C %m1:8 = %r1 non compris entre %r2 et %r3 866 : C 867 : MESCEL = MESERR(2) 868 : MOTERR(1:8) = MESCEL(1:8) 869 : REAERR(1) = REAL(VALER(2)) 870 : REAERR(2) = REAL(VAL1) 871 : REAERR(3) = REAL(VAL2) 871 : CALL ERREUR(42) 872 : C We artificially change the value of IERR in order to 873 : C continue the computation 874 : IERR = 0 875 : GOTO 9999 876 : ENDIF 877 : ENDIF 878 : C 879 : C*****Ecriture du CHPOINT contenant les "gamma". 880 : C 880 : CALL ECROBJ('CHPOINT',ITN2) 880 : CALL ECROBJ('CHPOINT',ITN1) 880 : CALL ECROBJ('CHPOINT',IPN2) 880 : CALL ECROBJ('CHPOINT',IPN1) 880 : CALL ECROBJ('CHPOINT',IVN2) 880 : CALL ECROBJ('CHPOINT',IVN1) 880 : CALL ECROBJ('CHPOINT',IRN2) 880 : CALL ECROBJ('CHPOINT',IRN1) 881 : 9999 CONTINUE 882 : C 883 : RETURN 884 : END 885 :
© Cast3M 2003 - All rights reserved.
Disclaimer