Télécharger stationary_discontinuity.dgibi
* fichier : stationary_discontinuity.dgibi *********************************************************************** *********************************************************************** *********************************************************************** ***** CALCUL DU TUBE A CHOC; CAS MULTIESPECE ***** ***** GAZ multi-especes "THERMALLY PERFECT" ***** ***** ***** ***** FORMULATION VF COMPRESSIBLE EXPLICITE ***** ***** Colella-Glaz, discontinuité de contact stationnaire ***** ***** ***** ***** A. BECCANTINI DRN/DMT/SEMT/LTMF FEVRIER 2000 ***** *********************************************************************** *********************************************************************** *********************************************************************** 'MESSAGE' 'A mettre a jours' ; 'FIN' ; 'OPTION' 'ELEM' 'QUA4' ; 'OPTION' 'ISOV' 'SULI' ; 'OPTION' 'ECHO' 0 ; 'OPTION' 'TRAC' 'X'; GRAPH = FAUX ; * GRAPH = VRAI ; * *** Methodes possibles : * * 'VLH' * 'CG' METO = 'CG' ; ******************************************************* ***** PROCEDURE POUR TESTER CONVERGENCE ***** ******************************************************* 'DEBPROC' CALCUL ; 'ARGUMENT' RVX*'TABLE' ; *KDOMA = RV . 'DOMAINE' ; KDOMA = RV . 'MODTOT' ; RNi = RV . 'INCO' . 'RNI' ; RN0 = RV . 'INCO' . 'RN0' ; * *** RN0 = solution à t = tn_M1; * DD = RV . 'PASDETPS' . 'NUPASDT' ; NINTE = 20 ; NN = DD/NINTE ; LO = (DD-(NINTE*NN)) EGA 0 ; 'SI' ( LO ) ; ELI = 'MAXIMUM' (RNi '-' RN0) 'ABS' ; ELI = ('LOG' (ELI + 1.0E-20))/('LOG' 10.) ; 'MESSAGE' ('CHAINE' 'ITER ' DD ' ERREUR LINF ' ELI ); RV . 'INCO' . 'IT' = (RV . 'INCO' . 'IT') 'ET' IT ; RV . 'INCO' . 'ER' = (RV . 'INCO' . 'ER') 'ET' ER ; 'FINSI' ; RV . 'INCO' . 'RN0' = 'COPIER' RNi; 'FINPROC' ; ***************************************************** ***************************************************** ** PROCEDURE EXEX POUR FORMULATION VF COMPRESSIBLE ** ** CAS MULTIESPECES "THERMALLY PERFECT" ** ***************************************************** ***************************************************** 'DEBPROC' EXEX ; 'ARGUMENT' RV*TABLE ; ******************************************* **** Recherche de RV . *KONV . IDEUL **** ******************************************* * **** Nom de la table RV . *'KONV' -> NOMT * 'REPETER' BL1 NBOP; MCEL = 'EXTRAIRE' &BL1 RV . 'LISTOPER'; 'QUITTER' BL1; 'FINSI' ; 'FIN' BL1; IEUL = RV . NOMT . 'KOPT' . 'IDEUL'; 'SI' ('NON' (IEUL 'EGA' 3)); 'MESSAGE' 'EULERMST ???'; 'ERREUR' 21; 'FINSI' ; * Mono-espece ou multi-especes 'SI' ('EXISTE' (RV . NOMT) 'ARG5') ; LOGMUL = VRAI ; 'SINON' ; LOGMUL = FAUX ; 'FINSI' ; * **** CL * LOGLIM = RV . 'INCO' . 'CLIM' ; ****************************************** **** Ordre en espace, ordre en temps **** ****************************************** ORD_ESP = RV . 'ORDREESP' ; ORD_TPS = RV . 'ORDRETPS' ; 'MESSAGE' '--------------------------'; 'MESSAGE' 'Ordre en Espace :' ord_esp; 'MESSAGE' 'Ordre en Temps :' ord_tps; 'MESSAGE' '--------------------------'; 'SI' ((ORD_ESP 'EGA' 1) 'ET' (ORD_TPS 'EGA' 2)); 'MESSAGE' ; 'MESSAGE' (CHAINE 'Ordre en Espace doit etre 2'); 'MESSAGE' (CHAINE 'On impose ça.'); 'MESSAGE' ; RV . 'ORDREESP' = 2 ; 'MESSAGE' ; 'MESSAGE' '--------------------------'; 'MESSAGE' 'Ordre en Espace :' ord_esp; 'MESSAGE' 'Ordre en Temps :' ord_tps; 'MESSAGE' '--------------------------'; 'FINSI' ; ****************************** **** La table 'PASDETPS' **** ****************************** TPSI = RV . 'TPSI' ; TFIN = RV . 'TFINAL'; RV . 'PASDETPS' . 'TPS' = TPSI; * **** DELTAT-1 est un argument de PRET (prediction) * Donc on doit l'initialiser. * RV . 'PASDETPS' . 'DELTAT-1' = 0.0D0; ********************* **** Les TABLES ***** ********************* * **** RV . 'INCO' * RV . 'DOMAINE' * RV . 'KIZD' * RV . 'KIZG' * **** RV . 'INCO' -> KINCO * KINCO = (RV . 'INCO') ; * **** RV . 'DOMAINE' -> KDOMA * *KDOMA = (RV . 'DOMAINE') ; KDOMA = (RV . 'MODTOT') ; KDOMA2 = (RV . 'DOMAINE') ; * **** RV . 'KIZD' contient les volumes des ELTs * 'SI' ('NON' ('EXISTE' RV 'KIZD')) ; 'KDIA' RV ; 'FINSI' ; * ***** RV . 'KIZG' contient les flux aux interfaces. * 'SI' ('NON' ('EXISTE' RV 'KIZG')) ; RV . 'KIZG' = 'TABLE' 'KIZG' ; 'FINSI' ; 'SI' LOGMUL ; ********************************************************* **** Multi-especes, boucle Sur les Pas de Temps **** ********************************************************* * **** Evaluation de coeff pour le calcule des pentes * KINCO . 'V' KINCO . 'P' KINCO . 'T' KINCO . 'Y' (KINCO . 'RNI') (KINCO . 'GNI') (KINCO . 'ENI') (KINCO. 'RYNI'); GRADR ALR COEFR = 'PENT' KDOMA 'CENTRE' 'EULESCAL' 'LIMITEUR' (KINCO . 'RNI'); GRADP ALP COEFP = 'PENT' KDOMA 'CENTRE' 'EULESCAL' 'LIMITEUR' (KINCO . 'P'); GRADV ALV COEFV = 'PENT' KDOMA 'CENTRE' 'EULEVECT' 'LIMITEUR' (KINCO . 'V'); GRADY ALY COEFY = 'PENT' KDOMA 'CENTRE' 'EULESCAL' 'LIMITEUR' (KINCO . 'Y'); I = 0 ; 'REPETER' BLOC1 (RV . 'ITMA') ; I = I + 1 ; * ***** Les variables primitives * KINCO . 'V' KINCO . 'P' KINCO . 'T' KINCO . 'Y' (KINCO . 'RNI') (KINCO . 'GNI') (KINCO . 'ENI') (KINCO. 'RYNI'); 'SI' (ORD_ESP 'EGA' 1) ; ROF VITF PF YF = 'PRET' 'PERFTEMP' ORD_ESP ORD_TPS KDOMA (KINCO . 'IPGAZ') (KINCO . 'RNI') (KINCO . 'V') (KINCO . 'P') (KINCO . 'Y') ; 'SINON'; * ***** Ordre 2 en espace => calcul des pentes * GRADR ALR = 'PENT' KDOMA 'CENTRE' 'EULESCAL' 'LIMITEUR' (KINCO . 'RNI') 'GRADGEO' COEFR ; GRADP ALP = 'PENT' KDOMA 'CENTRE' 'EULESCAL' 'LIMITEUR' (KINCO . 'P') 'GRADGEO' COEFP ; GRADV ALV = 'PENT' KDOMA 'CENTRE' 'EULEVECT' 'LIMITEUR' (KINCO . 'V') 'GRADGEO' COEFV ; GRADY ALY = 'PENT' KDOMA 'CENTRE' 'EULESCAL' 'LIMITEUR' (KINCO . 'Y') 'GRADGEO' COEFY ; 'SI' (ORD_TPS 'EGA' 1); ROF VITF PF YF = 'PRET' 'PERFTEMP' ORD_ESP ORD_TPS KDOMA (KINCO . 'IPGAZ') (KINCO . 'RNI') GRADR ALR (KINCO . 'V') GRADV ALV (KINCO . 'P') GRADP ALP (KINCO . 'Y') GRADY ALY ; 'SINON' ; * ********* Ordre 2 en temps * ROF VITF PF YF = 'PRET' 'PERFTEMP' ORD_ESP ORD_TPS KDOMA (KINCO . 'IPGAZ') (KINCO . 'RNI') GRADR ALR (KINCO . 'V') GRADV ALV (KINCO . 'P') GRADP ALP (KINCO . 'Y') GRADY ALY (KINCO . 'GAMMA') ((RV . 'PASDETPS' . 'DELTAT-1')/2.0); 'FINSI' ; 'FINSI' ; * *********** Creation de MCHAML de type 'FACEL' pour les * calcul de flux aux interfaces KINCO . 'RNF' = ROF ; KINCO . 'VITNF' = VITF ; KINCO . 'PNF' = PF ; KINCO . 'YF' = YF ; * ********* Boucle sur les operateurs * 'REPETER' BLOC2 NBOP ; NOMPER = 'EXTRAIRE' &BLOC2 (RV . 'LISTOPER'); ('TEXTE' NOMPER) (RV . NOTABLE) ; 'FIN' BLOC2 ; * ********* Mise a jour de la table RV . 'PASDETPS' * 'SI' ('EXISTE' RV 'DTI'); DTI = 'MINIMUM' (RV . 'PASDETPS' . 'DTCONV') ); 'SINON'; DTI = (RV . 'PASDETPS' . 'DTCONV'); 'FINSI'; RV . 'PASDETPS' . 'DELTAT' = DTI ; TMPS = RV . 'PASDETPS' . 'TPS'; DTI0 = TFIN '-' TMPS; 'SI' (DTI0 '<EG' DTI); DTI = DTI0; RV . 'PASDETPS' . 'DELTAT' = DTI ; LOGQUIT = VRAI; 'SINON' ; LOGQUIT = FAUX ; 'FINSI' ; * ******** On avance au pas de temps suivant * * N.B. Tn+1 = Tn + (CFL * RV . 'PASDETPS' . 'DELTAT'); * * ******** On detrui les choses qui ne servent plus * * * Les variables primitives * 'OUBL' KINCO 'V'; 'OUBL' KINCO 'P'; 'OUBL' KINCO 'T'; 'OUBL' KINCO 'Y'; 'OUBL' KINCO 'GAMMA'; * * Les MCHAML faces * 'DETR' ROF ; 'DETR' VITF ; 'DETR' PF ; 'DETR' YF; 'OUBL' KINCO 'RNF'; 'OUBL' KINCO 'VITNF'; 'OUBL' KINCO 'PNF'; 'OUBL' KINCO 'YF'; * * Les pentes 'ET' les limiteurs * 'SI' (ORD_ESP 'EGA' 2); 'DETR' GRADR; 'DETR' GRADP; 'DETR' GRADV; 'DETR' GRADY; 'DETR' ALR; 'DETR' ALP; 'DETR' ALV; 'DETR' ALY; 'FINSI' ; * ******** On impose le CL * * 'SI' LOGLIM; PROLIM RV ; 'FINSI' ; 'SI' LOGQUIT; 'QUITTER' BLOC1; 'FINSI'; 'MENAGE' ; 'FIN' BLOC1 ; ******************************************** **** FIN de Boucle Sur les Pas de Temps *** ******************************************** 'SINON' ; ********************************************************* **** Mono-espece, boucle sur les pas de temps **** ********************************************************* KINCO . 'V' KINCO . 'P' KINCO . 'T' (KINCO . 'RNI') (KINCO . 'GNI') (KINCO . 'ENI') ; GRADR ALR COEFR = 'PENT' KDOMA 'CENTRE' 'EULESCAL' 'LIMITEUR' (KINCO . 'RNI'); GRADP ALP COEFP = 'PENT' KDOMA 'CENTRE' 'EULESCAL' 'LIMITEUR' (KINCO . 'P'); GRADV ALV COEFV = 'PENT' KDOMA 'CENTRE' 'EULEVECT' 'LIMITEUR' (KINCO . 'V'); I = 0 ; 'REPETER' BLOC1 (RV . 'ITMA') ; I = I + 1 ; * ***** Les variables primitives * KINCO . 'V' KINCO . 'P' KINCO . 'T' (KINCO . 'RNI') (KINCO . 'GNI') (KINCO . 'ENI') ; 'SI' (ORD_ESP 'EGA' 1) ; ROF VITF PF = 'PRET' 'PERFTEMP' ORD_ESP ORD_TPS KDOMA (KINCO . 'IPGAZ') (KINCO . 'RNI') (KINCO . 'V') (KINCO . 'P') ; 'SINON'; * ***** Ordre 2 en espace => calcul des pentes * GRADR ALR = 'PENT' KDOMA 'CENTRE' 'EULESCAL' 'LIMITEUR' (KINCO . 'RNI') 'GRADGEO' COEFR ; GRADP ALP = 'PENT' KDOMA 'CENTRE' 'EULESCAL' 'LIMITEUR' (KINCO . 'P') 'GRADGEO' COEFP ; GRADV ALV = 'PENT' KDOMA 'CENTRE' 'EULEVECT' 'LIMITEUR' (KINCO . 'V') 'GRADGEO' COEFV ; 'SI' (ORD_TPS 'EGA' 1); ROF VITF PF = 'PRET' 'PERFTEMP' ORD_ESP ORD_TPS KDOMA (KINCO . 'IPGAZ') (KINCO . 'RNI') GRADR ALR (KINCO . 'V') GRADV ALV (KINCO . 'P') GRADP ALP ; 'SINON' ; * ********* Ordre 2 en temps * ROF VITF PF = 'PRET' 'PERFTEMP' ORD_ESP ORD_TPS KDOMA (KINCO . 'IPGAZ') (KINCO . 'RNI') GRADR ALR (KINCO . 'V') GRADV ALV (KINCO . 'P') GRADP ALP (KINCO . 'GAMMA') ((RV . 'PASDETPS' . 'DELTAT-1')/2.0); 'FINSI' ; 'FINSI' ; * *********** Creation de MCHAML de type 'FACEL' pour les * calcul de flux aux interfaces KINCO . 'RNF' = ROF ; KINCO . 'VITNF' = VITF ; KINCO . 'PNF' = PF ; * ********* Boucle sur les operateurs * 'REPETER' BLOC2 NBOP ; NOMPER = 'EXTRAIRE' &BLOC2 (RV . 'LISTOPER'); ('TEXTE' NOMPER) (RV . NOTABLE) ; 'FIN' BLOC2 ; * ********* Mise a jour de la table RV . 'PASDETPS' * 'SI' ('EXISTE' RV 'DTI'); DTI = 'MINIMUM' (RV . 'PASDETPS' . 'DTCONV') ); 'SINON'; DTI = (RV . 'PASDETPS' . 'DTCONV'); 'FINSI'; RV . 'PASDETPS' . 'DELTAT' = DTI ; TMPS = RV . 'PASDETPS' . 'TPS'; DTI0 = TFIN '-' TMPS; 'SI' (DTI0 '<EG' DTI); DTI = DTI0; RV . 'PASDETPS' . 'DELTAT' = DTI ; LOGQUIT = VRAI; 'SINON' ; LOGQUIT = FAUX ; 'FINSI' ; * ******** On avance au pas de temps suivant * * N.B. Tn+1 = Tn + (CFL * RV . 'PASDETPS' . 'DELTAT'); * * ******** On detrui les choses qui ne servent plus * * * Les variables primitives * 'OUBL' KINCO 'V'; 'OUBL' KINCO 'P'; 'OUBL' KINCO 'T'; 'OUBL' KINCO 'GAMMA'; * * Les MCHAML faces * 'DETR' ROF ; 'DETR' VITF ; 'DETR' PF ; 'OUBL' KINCO 'RNF'; 'OUBL' KINCO 'VITNF'; 'OUBL' KINCO 'PNF'; * * Les pentes 'ET' les limiteurs * 'SI' (ORD_ESP 'EGA' 2); 'DETR' GRADR; 'DETR' GRADP; 'DETR' GRADV; 'DETR' ALR; 'DETR' ALP; 'DETR' ALV; 'FINSI' ; * ******** On impose le CL * * 'SI' LOGLIM; PROLIM RV ; 'FINSI' ; 'SI' LOGQUIT; 'QUITTER' BLOC1; 'FINSI'; 'MENAGE' ; 'FIN' BLOC1 ; ******************************************** **** FIN de Boucle Sur les Pas de Temps *** ******************************************** 'FINSI' ; 'FINPROC' ; ***************************************************** ***************************************************** ** FIN PROCEDURE EXEX ** ***************************************************** ***************************************************** ***************************************************** ***************************************************** ***** PROCEDURE PROLIM ***** ***************************************************** ***************************************************** * * **** Cas Euler * 'DEBPROC' PROLIM ; 'ARGUMENT' RV*TABLE ; * **** RV . 'DOMAINE' -> KDOMA * *KDOMA = RV . 'DOMAINE' ; KDOMA = RV . 'MODTOT' ; KDOMA2 = RV . 'DOMAINE' ; KINCO = RV . 'INCO' ; KINCO . 'RNI' = 'KCHT' KDOMA 'SCAL' 'CENTRE' (KINCO . 'RNI') (KINCO . 'RNLIM'); (KINCO . 'GNI') (KINCO . 'GNLIM'); KINCO . 'ENI' = 'KCHT' KDOMA 'SCAL' 'CENTRE' (KINCO . 'ENI') (KINCO . 'ENLIM'); KINCO . 'RYNI' = 'KCHT' KDOMA 'SCAL' 'CENTRE' (KINCO . 'RYNLIM') ; 'FINPROC' ; ***************************************************** ***************************************************** ***** FIN PROCEDURE PROLIM ***** ***************************************************** ***************************************************** ************ * MAILLAGE * ************ RAF = 2 ; NY = 4 ; NX = 3 '*' RAF ; L = 1.0D0 ; DX = L '/' NX '/' 2.0D0; H = NY '*' DX ; xD = 0.5D0 '*' L ; xG = -1.0D0 '*' xD ; yH = 0.5D0 '*' H ; yB = -1.0D0 '*' yH ; A1 = xG yB ; A2 = 0.0D0 yB ; A3 = xD yB ; A4 = xD yH ; A5 = 0.0D0 yH ; A6 = xG yH ; VECTG = XG 0.0D0 ; VECTD = XD 0.0D0 ; xBG = xG '-' DX; XBD = xD '+' DX; B1 = xBG yB; B2 = xBG yH; B3 = xBD yB; B4 = xBD yH; DOM1 = LAT12 'TRANSLATION' NX VECTG ; DOM2 = LAT12 'TRANSLATION' NX VECTD ; VECTFG = (-1.0D0 '*' DX) 0.0D0; VECTFD = DX 0.0D0; FRONTG = LAT1 'TRANSLATION' 1 VECTFG; FRONTD = LAT2 'TRANSLATION' 1 VECTFD; * *** Rotation * ANGLE = 30.0D0; ORIG = 0.0D0 0.0D0; 'MESSAGE'; 'MESSAGE' (CHAIN 'ANGLE = ' ANGLE); 'MESSAGE'; DOM1 = DOM1 'TOURNER' ANGLE ORIG; DOM2 = DOM2 'TOURNER' ANGLE ORIG; FRONTG = FRONTG 'TOURNER' ANGLE ORIG; FRONTD = FRONTD 'TOURNER' ANGLE ORIG; DOMINT = DOM1 'ET' DOM2 ; 'ELIMINATION' DOMINT 1D-6; FRONT = FRONTG 'ET' FRONTD ; 'ELIMINATION' FRONT 1D-6; DOMTOT = DOMINT 'ET' FRONT; 'ELIMINATION' DOMTOT 1D-6; ********************** *** OBJETS MODELES *** ********************** MDOMTOT = 'CHANGER' DOMTOT 'QUAF' ; MDOMINT = 'CHANGER' DOMINT 'QUAF' ; MDOM1 = 'CHANGER' DOM1 'QUAF' ; MDOM2 = 'CHANGER' DOM2 'QUAF' ; MFRONTG = 'CHANGER' FRONTG 'QUAF' ; MFRONTD = 'CHANGER' FRONTD 'QUAF' ; MFRONT = 'CHANGER' FRONT 'QUAF' ; MFRONTG 'ET' MFRONTD 'ET' MFRONT) 1.E-5; MDNS = 'EULER' ; * ******* Creation de la ligne Utilisée pour le Post-Traitement * reliant les points centres * XINIT = XG '+' (0.5D0 '*' DX) ; YINIT = YB '+' (0.5D0 '*' DX) ; XFIN = XD '-' (0.5D0 '*' DX) ; YFIN = YINIT ; PINI = XINIT YINIT; PFIN = XFIN YFIN ; IAUX = (2 '*' NX) '-' 1 ; COURB = PINI 'DROIT' IAUX PFIN; COURB = COURB 'TOURNER' ANGLE ORIG; COURB = COURB 'COULEUR' 'VERT'; 'SI' GRAPH ; 'TITRE' ('CHAINE' 'Maillage '); 'FINSI' ; *********************** **** LA TABLE PGAZ **** *********************** PGAZ = 'TABLE' ; * **** Ordre des polynoms cv_i * PGAZ . 'NORD' = 4 ; * **** Especes qui sont dans les equations d'Euler * * **** Espece qui n'y est pas * PGAZ . 'ESPNEULE' = 'O2 '; * PGAZ . 'N2 ' = 'TABLE' ; PGAZ . 'O2 ' = 'TABLE' ; * **** R (J/Kg/K) * PGAZ . 'N2 ' . 'R' = 296.8 ; PGAZ . 'O2 ' . 'R' = 259.8 ; * **** Regressions polynomials * 2.33636971E-08 -1.53304905E-12; 8.78233606E-09 -3.05514485E-13 ; * **** "Enthalpies" (ou energies) de formations a OK (J/Kg) * Note: ce sont des entites numeriques * PGAZ . 'O2 ' . 'H0K' = -2.634D5 ; PGAZ . 'N2 ' . 'H0K' = -2.953D5 ; ********************************************************** *** Fin PGAZ ********************************************* ********************************************************** * *** Etatg * * s = 0 ; YN2g = 0.77785 ; UMYN2g = 1.0 '-' YN2g ; RTOTg = ((PGAZ . 'O2 ' . 'R') * UMYN2g) '+' ((PGAZ . 'N2 ' . 'R') * YN2G ) ; pg = 1.0D6; rog = 1.0D2 ; Tg = pg '/' (rog '*' RTOTg) ; ung = s ; utg = 0.0 ; * A1 = (('EXTRAIRE' (PGAZ . 'O2 ' . 'A') 1) '*' UMYN2g) '+' (('EXTRAIRE' (PGAZ . 'N2 ' . 'A') 1) '*' YN2g); A2 = (('EXTRAIRE' (PGAZ . 'O2 ' . 'A') 2) '*' UMYN2g) '+' (('EXTRAIRE' (PGAZ . 'N2 ' . 'A') 2) '*' YN2g); A3 = (('EXTRAIRE' (PGAZ . 'O2 ' . 'A') 3) '*' UMYN2g) '+' (('EXTRAIRE' (PGAZ . 'N2 ' . 'A') 3) '*' YN2g); A4 = (('EXTRAIRE' (PGAZ . 'O2 ' . 'A') 4) '*' UMYN2g) '+' (('EXTRAIRE' (PGAZ . 'N2 ' . 'A') 4) '*' YN2g); A5 = (('EXTRAIRE' (PGAZ . 'O2 ' . 'A') 5) '*' UMYN2g) '+' (('EXTRAIRE' (PGAZ . 'N2 ' . 'A') 5) '*' YN2g); T2 = Tg '*' Tg ; T3 = T2 '*' Tg ; T4 = T3 '*' Tg ; T5 = T4 '*' Tg ; ETHERg = ((A1 * Tg) '+' ((A2 '/' 2.0) * T2) '+' ((A3 '/' 3.0) * T3) '+' ((A4 '/' 4.0) * T4) '+' ((A5 '/' 5.0) * T5)) ; retg = rog '*' (ETHERg '+' (0.5 '*' ung '*' ung)) ; rouxg = ((ung '*' ('COS' ANGLE)) '-' (utg '*' ('SIN' ANGLE))) '*' rog ; rouyg = ((ung '*' ('SIN' ANGLE)) '+' (utg '*' ('COS' ANGLE))) '*' rog; * *** Etatd * YN2d = 0.77785 ; UMYN2d = 1.0 '-' YN2d ; RTOTd = ((PGAZ . 'O2 ' . 'R') * UMYN2d) '+' ((PGAZ . 'N2 ' . 'R') * YN2D ) ; pd = pg ; rod = rog '/' 5.0 ; Td = pd '/' (rod '*' RTOTd) ; und = ung ; utd = utg ; * A1 = (('EXTRAIRE' (PGAZ . 'O2 ' . 'A') 1) '*' UMYN2d) '+' (('EXTRAIRE' (PGAZ . 'N2 ' . 'A') 1) '*' YN2d); A2 = (('EXTRAIRE' (PGAZ . 'O2 ' . 'A') 2) '*' UMYN2d) '+' (('EXTRAIRE' (PGAZ . 'N2 ' . 'A') 2) '*' YN2d); A3 = (('EXTRAIRE' (PGAZ . 'O2 ' . 'A') 3) '*' UMYN2d) '+' (('EXTRAIRE' (PGAZ . 'N2 ' . 'A') 3) '*' YN2d); A4 = (('EXTRAIRE' (PGAZ . 'O2 ' . 'A') 4) '*' UMYN2d) '+' (('EXTRAIRE' (PGAZ . 'N2 ' . 'A') 4) '*' YN2d); A5 = (('EXTRAIRE' (PGAZ . 'O2 ' . 'A') 5) '*' UMYN2d) '+' (('EXTRAIRE' (PGAZ . 'N2 ' . 'A') 5) '*' YN2d); T2 = Td '*' Td ; T3 = T2 '*' Td ; T4 = T3 '*' Td ; T5 = T4 '*' Td ; ETHERd = ((A1 * Td) '+' ((A2 '/' 2.0) * T2) '+' ((A3 '/' 3.0) * T3) '+' ((A4 '/' 4.0) * T4) '+' ((A5 '/' 5.0) * T5)) ; rouxd = ((und '*' ('COS' ANGLE)) '-' (utd '*' ('SIN' ANGLE))) '*' rod; rouyd = ((und '*' ('SIN' ANGLE)) '+' (utd '*' ('COS' ANGLE))) '*' rod; retd = rod '*' (ETHERd '+' (0.5 '*' und '*' und)) ; * *** ro * * *** ro u, ro v * * *** ro e * * *** ro y * ******************************************************** *** CREATION DE 'MODE' POUR GRAPHIQUER LE CHAMELEM *** ******************************************************** * **** Les debits dans le repaire tube * GNN = (GNX * ('COS' ANGLE)) '+' (GNY * ('SIN' ANGLE)); GNT = (GNY * ('COS' ANGLE)) '-' (GNX * ('SIN' ANGLE)); * *** Les variables primitives * VN PN TN YN GAMN = 'PRIM' 'PERFTEMP' PGAZ RN GN REN RYN ; RN0 = 'COPIER' RN ; * *** GRAPHIQUE DES C.I. * 'SI' GRAPH ; * *** CREATION DE CHAMELEM * 'FINSI' ; *********************** **** La table EQEX **** *********************** NITER = 100 ; * TFIN = 0.0005 ; TFIN = 100 ; 'TPSI' 0. 'TFINAL' TFIN * ***** Option VF = volumes finis * CONS = conservative * EULER = euler monoespece * VANLEER = methode * 'OPTI' 'VF' 'CONS' 'EULERMST' METO * ***** Procedure CALCUL * 'ZONE' $DOMTOT 'OPER' 'CALCUL' * ***** Operateur 'KONV' * 'ZONE' $DOMTOT 'OPER' 'KONV' * ***** Arguments de 'KONV' (maximum 8 chatacters) * 'IPGAZ' 'RNF' 'VITNF' 'PNF' 'YF' * ***** LIST des inconnues * 'INCO' 'RNI' 'GNI' 'ENI' 'RYNI' ; * *** La table RV . INCO (de soustype INCO); * RV . 'INCO' = TABLE 'INCO' ; * *** Stokage des inconnues des equations d'Euler. * RV . 'INCO' . 'RNI' = 'COPIER' RN ; RV . 'INCO' . 'GNI' = 'COPIER' GN ; RV . 'INCO' . 'ENI' = 'COPIER' REN ; RV . 'INCO' . 'RYNI' = 'COPIER' RYN ; * *** Le gaz * RV . 'INCO' . 'IPGAZ' = PGAZ ; * *** CONDITIONS LIMITS * RV . 'INCO' . 'CLIM' = VRAI ; RV . 'INCO' . 'RYNLIM' = 'KCHT' $FRONT 'SCAL' 'CENTRE' 'COMP' 'N2 ' RYN ; * **** Ordre en espace * Ordre en temps * IE = 2 ; IT = 2 ; RV . 'ORDREESP' = IE ; RV . 'ORDRETPS' = IT ; * **** Test de la convergence * RV . 'INCO' . 'RN0' = 'COPIER' RN ; *??? RV . 'MODTOT' = $DOMTOT ; * *********************************** *** Execution EXEX *** *********************************** * 'MESSAGE' ; 'MESSAGE' '**************************'; 'MESSAGE' ('CHAINE' 'METHODE : ' METO) ; 'MESSAGE' '**************************'; 'MESSAGE' ; 'TEMPS' ; EXEX RV ; 'TEMPS' ; TFIN = RV. 'PASDETPS'. 'TPS'; * **** SOLUTIONS * * **** Les variables conservatives * RN = 'COPIER' (RV . 'INCO' . 'RNI'); GN = 'COPIER' (RV . 'INCO'. 'GNI'); REN = 'COPIER' (RV . 'INCO'. 'ENI'); RYN = 'COPIER' (RV . 'INCO'. 'RYNI'); * **** Les variables primitives * VN PN TN YN GAMN = 'PRIM' 'PERFTEMP' PGAZ RN GN REN RYN ; * **** La vitesse dans le repaire tube * VNN = (VNX * ('COS' ANGLE)) '+' (VNY * ('SIN' ANGLE)); VNT = (VNY * ('COS' ANGLE)) '-' (VNX * ('SIN' ANGLE)); GNN = VNN * RN; GNT = VNT * RN; * *** GRAPHIQUE DES SOLUTIONS * 'SI' GRAPH ; * *** CREATION DE CHAMELEM * 'FINSI' ; * *** Objects evolutions * xx yy = 'COORDONNEE' Courb; x0 = 'MINIMUM' lxx; x1 = 'MAXIMUM' lxx ; * ro tro = CHAINE '1D ' METO ' : RO IT ' IT ' IE ' IE ' tmps ' TFIN ' elem ' 'QUA4' ; * u tu = CHAINE '1D ' METO ' : U IT ' IT ' IE ' IE ' tmps ' TFIN ' elem ' 'QUA4' ; * v tv = CHAINE '1D ' METO ' : V IT ' IT ' IE ' IE ' tmps ' TFIN ' elem ' 'QUA4' ; 'SI' GRAPH ; 'DESSIN' evv 'TITRE' tv 'XBOR' x0 x1; 'FINSI' ; * p tp = CHAINE '1D ' METO ' : P IT ' IT ' IE ' IE ' tmps ' TFIN ' elem ' 'QUA4' ; * gamma tgam = CHAINE '1D ' METO ' : GAMMA IT ' IT ' IE ' IE ' tmps ' TFIN ' elem ' 'QUA4' ; * T tT = 'CHAINE' '1D ' METO ' : T IT ' IT ' IE ' IE ' tmps ' TFIN ' elem ' 'QUA4' ; * Fractions massiques EVY = 'TABLE' ; TY = CHAINE '1D ' METO ' : Y IT ' IT ' IE ' IE ' tmps ' TFIN ' elem ' 'QUA4' ; YOLD = 'COPIER' YN ; 'REPETER' BLOCCO NUMECO; NOMCEL = 'EXTRAIRE' NOMECO &BLOCCO ; Y0 = 'KCHT' $DOMTOT 'SCAL' 'CENTRE' 'COMP' NOMCEL N0 = NOMCEL ; EVY . &BLOCCO = EVY0; 'FIN' BLOCCO; 'LOG|E|inf' (RV . 'INCO' . 'ER') ; 'SI' GRAPH ; 'DESSIN' evro 'TITRE' tro 'XBOR' x0 x1 ; 'DESSIN' evT 'TITRE' tT 'XBOR' x0 x1 ; 'DESSIN' evp 'TITRE' tp 'XBOR' x0 x1 ; 'DESSIN' evu 'TITRE' tu 'XBOR' x0 x1 ; 'DESSIN' evgam 'TITRE' tgam 'XBOR' x0 x1 ; 'DESSIN' EVERR 'TITRE' ('CHAINE' METO ' IE =' IE ' IT =' IT) ; 'FINSI' ; * TESTS ERRO = ('MAXIMUM' (RN '-' RN0) 'ABS') '/' ('MAXIMUM' RN0); 'SI' (ERRO > 1.0D-4) ; 'MESSAGE' 'CG perde la discontinuite de contact stationnaire' ; 'ERREUR' 5; 'FINSI'; 'FIN' ;
© Cast3M 2003 - Tous droits réservés.
Mentions légales