* fichier : crebe12.dgibi *-------------------------------------------- * GIBIANE file for the HDR test E12.3.2 * * in 2D * * We make 100 iterations on the coarse mesh * * and compare the overpressure in the room * * R1904 with earlier computed values * * CREBCOM combustion model is used * *-------------------------------------------- *-------------------------------------------- * Here we construct the mesh for the * * HDR test E12.3.2 * *-------------------------------------------- 'OPTION' 'ECHO' 1 'DIME' 2 'ELEM' 'QUA4' 'TRAC' 'X' ; *------------------------------------------- * Given a line and a mesh, this procedure * computes a new line which containes the * points of the mesh closest to the points * of the given line *------------------------------------------- 'DEBPROC' LIGPROC ; ******* inputs 'ARGUMENT' LIGG*'MAILLAGE' MAI1*'MAILLAGE' ; NPOI = 'NBNO' LIGG ; MAI2 = 'CHANGER' MAI1 'POI1' ; * First element of the line P1 = LIGG 'POIN' 1 ; P1P = MAI2 'POIN' 'PROC' P1 ; P2 = LIGG 'POIN' 2 ; P2P = MAI2 'POIN' 'PROC' P2 ; 'REPETER' BL1 (NPOI '-' 1) ; 'SI' ('NEG' P1P P2P) ; NPOI = (NPOI '-' &BL1) '+' 1 ; 'QUITTER' BL1 ; 'FINSI' ; P2 = LIGG 'POIN' (&BL1 '+' 2) ; P2P = LIGG 'POIN' 'PROC' P2 ; 'FIN' BL1 ; LIGRES = 'MANUEL' 'SEG2' P1P P2P ; * The other element of the line 'REPETER' BL1 (NPOI '-' 2) ; P1P = P2P ; P2 = LIGG 'POIN' (&BL1 '+' 2) ; P2P = MAI2 'POIN' 'PROC' P2 ; 'SI' ('NEG' P1P P2P) ; LIGRES = LIGRES 'ET' ('MANUEL' 'SEG2' P1P P2P) ; 'FINSI' ; 'FIN' BL1 ; 'RESPRO' LIGRES ; 'FINPROC' ; *-------------------------------------------------------- *---refinement of the mesh (10 or 16) RAF = 10 ; *-- Length of the thirst and second room L1 = 7.32D0 ; *-- Height of the first window H1 = 0.5D0 ; *-- Height of the first and second room H2 = 4.0D0 ; *--------------------------------- * Number of cells in each region *--------------------------------- N1 = RAF ; *---- number of points in window1 ** It is BETTER to have NY1 even! NY1 = 'ENTIER' (RAF '/' 4) ; *--------------------------------- DX = L1 '/' N1 ; *---- cell size in window DX1 = H1 '/' NY1 ; *---- number of points in the outer window NY2 = 'ENTIER' ((H2 '-' H1) '/' DX) ; ***---- wall thickness DW = H1 ; *----------------------------------------- * First and Second Room (R1904 and R1905) *----------------------------------------- LB = L1 '-' (7.5 '*' H1) ; A1 = 0.0 0.0 ; A12 = LB 0.0 ; A13 = LB (H1 '/' 2) ; A14 = 0.0 (H1 '/' 2) ; *------------------------------------ * "Lower-left" domain *------------------------------------- NY1M = NY1 '/' 2 ; A1A12 = A1 'DROIT' A12 'DINI' DX1 'DFIN' (DW '/' NY1) ; A12A13 = A12 'DROIT' NY1M A13 ; A13A14 = A13 'DROIT' A14 'DINI' (DW '/' NY1) 'DFIN' DX1 ; A14A1 = A14 'DROIT' NY1M A1 ; DOM1L = 'DALLER' A1A12 A12A13 A13A14 A14A1 'PLAN' ; *-------------------------------------- * "Lower-right" domain *-------------------------------------- LB1 = L1 '-' (2.5 '*' H1) ; AI1 = LB1 0.0 ; A2 = L1 0.0 ; AI3 = L1 (H1 '/' 2) ; AI4 = LB1 (H1 '/' 2) ; *------------------------------------ DENF = (DW '/' NY1) ; AI1A2 = AI1 'DROIT' A2 'DINI' DENF 'DFIN' DENF ; A2AI3 = A2 'DROIT' NY1M AI3 ; AI3AI4 = AI3 'DROIT' AI4 'DINI' DENF 'DFIN' DENF ; AI4AI1 = AI4 'DROIT' NY1M AI1 ; DOM1R = 'DALLER' AI1A2 A2AI3 AI3AI4 AI4AI1 'PLAN' ; DOM1 = DOM1L 'ET' DOM1R ; 'ELIMINATION' DOM1 1.0D-6 ; *---------------------------------------------------- * Line for postreatment *---------------------------------------------------- *** NNN is the number of points in the disc NNN = NY1 '*' 4 ; xcel = LB '+' (H1 '/' 2) ; ycel = H1 '/' 2 ; ACEL = xcel ycel ; LIGB = A13 'DROIT' NY1M ACEL ; *------------------ xmm = LB '+' (2.5 '*' H1) ; 'REPETER' BL1 NNN ; ACEL0 = ACEL ; xcel = xcel '+' DX1 ; dop = 18.0625 '*' H1 '*' H1 ; dop = dop '-' ((xcel-xmm) '*' (xcel-xmm)) ; dop = dop '**' 0.5 ; ycel = dop '-' (3.25 '*' H1) ; ACEL = xcel ycel ; LIGB = LIGB 'ET' (ACEL0 'DROIT' 1 ACEL) ; 'FIN' BL1 ; LIGB = LIGB 'ET' (ACEL 'DROIT' AI4 'DINI' DX1 'DFIN' DX1 ) ; *-------------------------------------------- LIGBG = ('INVERSE' A13A14) 'ET' LIGB 'ET' ('INVERSE' AI3AI4) ; *------------------------------------- * "HIGH" domain *------------------------------------- A3 = L1 H1 ; A4 = 0.0 H1 ; A5 = L1 H2 ; A6 = 0.0 H2 ; AI5 = LB H2 ; *--------------------------------- AI3A5 = AI3 'DROIT' A5 'DINI' DX1 'DFIN' DX1 ; A5AI5 = A5 'DROIT' AI5 'DINI' DX1 'DFIN' DX1 ; AI5A6 = AI5 'DROIT' A6 'DINI' DX1 'DFIN' DX1 ; LIGH = A5AI5 'ET' AI5A6 ; A6A14 = A6 'DROIT' A14 'DINI' DX1 'DFIN' DX1 ; *--------------------------- DOM2 = 'DALLER' LIGBG AI3A5 LIGH A6A14 'PLAN' ; *--------------------------------------------- * The Channel (junction between rooms) *--------------------------------------------- A21 = A2 'PLUS' (DW 0.0) ; A31 = A3 'PLUS' (DW 0.0) ; A2A21 = A2 'DROIT' NY1 A21 ; A21A31 = A21 'DROIT' NY1 A31 ; A31A3 = A31 'DROIT' NY1 A3 ; A3A2 = A3 'DROIT' NY1 A2 ; *------------------------- DOMC = 'DALLER' A2A21 A21A31 A31A3 A3A2 'PLAN' ; *-------------------------------------------- * Third room (R1801) *-------------------------------------------- *---- Length of the third room L3 = 5.0D0 ; *--- length before the window LB = 1.8D0 ; *--- width of window LW = 0.668D0 ; *--- length after the window LA = L3 '-' LB '-' LW ; *---------------------------------------- * Low Domain *---------------------------------------- H3 = (9.5 '-' H2) '-' H1 ; MH3 = (-1.0D0) '*' H3 ; A70 = (L1 '+' DW) MH3 ; A71 = (L1 '+' LB '+' DW) MH3 ; A72 = (L1 '+' LB '+' DW) 0.0 ; *---------------------------- DW1 = LW '/' NY1 ; A70A71 = A70 'DROIT' A71 'DINI' (DW '/' NY1) 'DFIN' DW1 ; A71A72 = A71 'DROIT' A72 'DINI' DX1 'DFIN' DX1 ; A72A21 = A72 'DROIT' A21 'DINI' DW1 'DFIN' (DW '/' NY1) ; A21A70 = A21 'DROIT' A70 'DINI' DX1 'DFIN' DX1 ; *------------------------------ DOM3A = 'DALLER' A70A71 A71A72 A72A21 A21A70 'PLAN' ; *------------------------------ A7W = A71 'PLUS' (LW 0.0) ; A7U = A72 'PLUS' (LW 0.0) ; *-------------------------- A71A7W = A71 'DROIT' NY1 A7W ; A7WA7U = A7W 'DROIT' A7U 'DINI' DX1 'DFIN' DX1 ; A7UA72 = A7U 'DROIT' NY1 A72 ; A72A71 = 'INVERSE' (A71A72) ; *------------------------- DOM3B = 'DALLER' A71A7W A7WA7U A7UA72 A72A71 'PLAN' ; *------------------------- A81 = A7W 'PLUS' (LA 0.0) ; A82 = A7U 'PLUS' (LA 0.0) ; *------------------------- A7WA81 = A7W 'DROIT' A81 'DINI' DW1 'DFIN' DX1 ; A81A82 = A81 'DROIT' A82 'DINI' DX1 'DFIN' DX1 ; A82A7U = A82 'DROIT' A7U 'DINI' DX1 'DFIN' DW1 ; A7UA7W = 'INVERSE' (A7WA7U) ; *------------------------- DOM3C = 'DALLER' A7WA81 A81A82 A82A7U A7UA7W 'PLAN' ; *------------------------- DOM3 = DOM3A 'ET' DOM3B 'ET' DOM3C ; 'ELIMINATION' DOM3 1.0D-6 ; *----------------------------------------- * Middle domain *----------------------------------------- A73 = A72 'PLUS' (0.0 H1) ; ***--- first layer A21A72 = A21 'DROIT' A72 'DINI' (DW '/' NY1) 'DFIN' DW1 ; A72A73 = A72 'DROIT' NY1 A73 ; A73A31 = A73 'DROIT' A31 'DINI' DW1 'DFIN' (DW '/' NY1) ; A31A21 = A31 'DROIT' NY1 A21 ; *----------------------------- DOM4A = 'DALLER' A21A72 A72A73 A73A31 A31A21 'PLAN' ; *------------------------------------------- A7M = A7U 'PLUS' (0.0 H1) ; *--------------- A72A7U = A72 'DROIT' NY1 A7U ; A7UA7M = A7U 'DROIT' NY1 A7M ; A7MA73 = A7M 'DROIT' NY1 A73 ; A73A72 = 'INVERSE' (A72A73) ; *--------------------- DOM4B = 'DALLER' A72A7U A7UA7M A7MA73 A73A72 'PLAN' ; *--------------------------- A83 = A82 'PLUS' (0.0 H1) ; *--------------------------- A7UA82 = A7U 'DROIT' A82 'DINI' DW1 'DFIN' DX1 ; A82A83 = A82 'DROIT' NY1 A83 ; A83A7M = A83 'DROIT' A7M 'DINI' DX1 'DFIN' DW1 ; A7MA7U = 'INVERSE' (A7UA7M) ; *--------------------------- DOM4C = 'DALLER' A7UA82 A82A83 A83A7M A7MA7U 'PLAN' ; *--------------------------- DOM4 = DOM4A 'ET' DOM4B 'ET' DOM4C ; 'ELIMINATION' DOM4 1.0D-6 ; *----------------------------------------- * High domain *----------------------------------------- A32 = A31 'PLUS' (0.0 H2) ; A74 = A73 'PLUS' (0.0 H2) ; *--------------------------- A31A73 = 'INVERSE' (A73A31) ; A73A74 = A73 'DROIT' A74 'DINI' DX1 'DFIN' DW1 ; A74A32 = A74 'DROIT' A32 'DINI' DW1 'DFIN' (DW '/' NY1) ; A32A31 = A32 'DROIT' A31 'DINI' DW1 'DFIN' DX1 ; *--------------------------- DOM5A = 'DALLER' A31A73 A73A74 A74A32 A32A31 'PLAN' ; *--------------------------- A7H = A7M 'PLUS' (0.0 H2) ; *--------------------------- A73A7M = 'INVERSE' (A7MA73) ; A7MA7H = A7M 'DROIT' A7H 'DINI' DX1 'DFIN' DW1 ; A7HA74 = A7H 'DROIT' NY1 A74 ; A74A73 = 'INVERSE' (A73A74) ; *--------------------------- DOM5B = 'DALLER' A73A7M A7MA7H A7HA74 A74A73 'PLAN' ; *--------------------------- A84 = A83 'PLUS' (0.0 H2) ; *--------------------------- A7MA83 = A7M 'DROIT' A83 'DINI' DW1 'DFIN' DX1 ; A83A84 = A83 'DROIT' A84 'DINI' DX1 'DFIN' DW1 ; A84A7H = A84 'DROIT' A7H 'DINI' DX1 'DFIN' DW1 ; A7HA7M = 'INVERSE' (A7MA7H) ; *----------------------------- DOM5C = 'DALLER' A7MA83 A83A84 A84A7H A7HA7M 'PLAN' ; *------------------------------------- DOM5 = DOM5A 'ET' DOM5B 'ET' DOM5C ; 'ELIMINATION' DOM5 1.0D-6 ; ****----------*****----------------------******------- * Upper vent *---**********-----**********************------******* A75 = A74 'PLUS' (0.0 LW) ; A7F = A7H 'PLUS' (0.0 LW) ; *--------------------------- A74A7H = A74 'DROIT' NY1 A7H ; A7HA7F = A7H 'DROIT' NY1 A7F ; A7FA75 = A7F 'DROIT' NY1 A75 ; A75A74 = A75 'DROIT' NY1 A74 ; *--------------------------- DOMW = 'DALLER' A74A7H A7HA7F A7FA75 A75A74 'PLAN' ; FRONTH = A7FA75 ; *----------------------------------* * Ignition Zone * *----------------------------------* ** height at which we ignite HIGN = H2 '/' 2.0D0 ; *------------------------- PP1 = 0.0 HIGN ; PP2 = PP1 'PLUS' (0.0 (DX1)) ; PP3 = PP2 'PLUS' (0.0 (DX1)) ; gg = (-1.0D0) '*' DX '*' 0.6 ; PP5 = PP1 'PLUS' (0.0 gg) ; PP51 = PP5 'PLUS' (0.0 gg) ; PO1 = A6A14 'POIN' 'PROCHE' PP1 ; PO2 = A6A14 'POIN' 'PROCHE' PP2 ; PO31 = A6A14 'POIN' 'PROCHE' PP3 ; PO5 = A6A14 'POIN' 'PROCHE' PP5 ; PO51 = A6A14 'POIN' 'PROCHE' PP51 ; **-------------------------------- PPP = DX1 0.0 ; PA1 = A1A12 'POIN' 'PROCHE' PPP ; PPP = (2.0 '*' DX) 0.0 ; PA2 = A1A12 'POIN' 'PROCHE' PPP ; XXX1 = 'COORDONNEE' 1 PA1 ; XXX2 = 'COORDONNEE' 1 PA2 ; YYY2 = 'COORDONNEE' 2 PO2 ; YYY3 = 'COORDONNEE' 2 PO31 ; PO4 = XXX1 YYY2 ; PO7 = XXX1 YYY3 ; *----------------------------- * Just eight cells *----------------------------- PO1PO2 = PO1 'DROIT' 1 PO2 ; PO2PO31 = PO2 'DROIT' 1 PO31 ; PO4PO7 = PO4 'DROIT' 1 PO7 ; PO5PO51 = PO5 'DROIT' 1 PO51 ; **------------------------------ 'SI' (RAF 'EGA' 10) ; NZD = 1 ; KEF = -1.0D0 ; 'SINON' ; KEF = -2.0D0 ; NZD = 2 ; 'FINSI' ; DIG1 = 'TRANSLATION' PO1PO2 NZD ((KEF '*' XXX1) 0.0) ; DIG5 = 'TRANSLATION' PO2PO31 NZD ((KEF '*' XXX1) 0.0) ; DIG6 = 'TRANSLATION' PO4PO7 1 ((XXX2-XXX1) 0.0) ; DIG7 = 'TRANSLATION' PO5PO51 1 (XXX1 0.0) ; 'SI' (RAF '>' 10) ; YO31 = 'COORDONNEE' 2 PO31 ; YO2 = 'COORDONNEE' 2 PO1 ; DIFY = YO31 '-' YO2 ; VECUP = (0.0 DIFY) ; DIG6 = DIG1 'PLUS' VECUP ; DIG7 = DIG5 'PLUS' VECUP ; 'FINSI' ; *----------------------------- 'SI' (RAF '>' 10) ; DOMI = DIG1 'ET' DIG5 'ET' DIG6 'ET' DIG7 ; 'SINON' ; DOMI = DIG1 'ET' DIG5 ; 'FINSI' ; 'ELIMINATION' DOMI 1.0D-6; * 'TRACER' DOMI ; *---------------------------------------------------- * Total Domain *---------------------------------------------------- DOMTOT = DOM1 'ET' DOM2 'ET' DOM3 'ET' DOM4 'ET' DOM5 'ET' DOMC 'ET' DOMW 'ET' DOMI ; 'ELIMINATION' DOMTOT 1D-6 ; *------------------------------------- * Total Domain without ignition region *------------------------------------- DREST = 'DIFF' DOMTOT DOMI ; *---------------------------------- * The fird room (R1801) *---------------------------------- DTHIRD = DOM3 'ET' DOM4 'ET' DOM5 'ET' DOMC 'ET' DOMW ; 'ELIMINATION' DTHIRD 1.0D-6 ; *---------------------------------------- * First and second room (R1904 and R1905) *---------------------------------------- DSEC = DOM1 'ET' DOM2 ; *------------------------------- * Lignes de postreatment *------------------------------- LIGA = LIGBG 'PLUS' (0.0 (DX1 '/' 2.0)) ; A3A31 = 'INVERSE' A31A3 ; A31A73 = 'INVERSE' A73A31 ; A73A7M = 'INVERSE' A7MA73 ; A7MA83 = 'INVERSE' A83A7M ; LIGT = A3A31 'ET' A31A73 'ET' A73A7M 'ET' A7MA83 ; 'ELIMINATION' LIGT 1.0D-6 ; dis1 = (DX1 '/' 2) '-' (DW '/' 2.0) ; LIGB2 = LIGT 'PLUS' (0.0 dis1) ; LIGC = LIGA 'ET' LIGB2 ; 'ELIMINATION' LIGC 1.0D-6; LIG1 = LIGC ; * 'TRACER' (DOMTOT 'ET' LIG1) ; *------------------------------- * Points for capturing variables *------------------------------- *--------- inside the first two rooms dish = (H2 '/' 2) '-' (DX1 '/' 2) ; disl = (H1 '/' 2) '+' (DX1 '/' 2) ; PCEL1 = (L1 '/' 6) dish ; PCEL2 = (L1 '/' 6) disl ; PCEL3 = (L1 '/' 2) dish ; PCEL4 = (L1 '/' 2) disl ; PCEL5 = (L1 '-' (L1 '/' 6)) dish ; PCEL6 = (L1 '-' (L1 '/' 6)) disl ; *-------- inside the third room dis1 = L1 '+' (2.0 '*' DW) ; PCEL7 = dis1 (2.0 '*' DW) ; dis2 = (LB '+' (LW '/' 2) '+' (DX1 '/' 2.0)) '+' L1 '+' DW ; PCEL8 = dis2 (DW '-' 5.0) ; PCEL9 = dis2 (4.5 '-' DW) ; * 'TRACER' (DOMTOT 'ET' PCEL1 'ET' PCEL2 'ET' PCEL3 'ET' PCEL4 * 'ET' PCEL5 'ET' PCEL6 'ET' PCEL7 'ET' PCEL8 'ET' PCEL9) ; * * 'TRACER' (('CONTOUR' DREST) * 'ET' PCEL1 'ET' PCEL2 'ET' PCEL3 'ET' PCEL4 * 'ET' PCEL5 'ET' PCEL6 'ET' PCEL7 'ET' PCEL8 'ET' PCEL9) ; *---------------------------------- * Creating the models *---------------------------------- MDNS = 'EULER' ; *--------------------------------- MDOMI = 'MODELISER' DOMI MDNS ; MDREST = 'MODELISER' DREST MDNS ; MFRONTH = 'MODELISER' FRONTH MDNS ; MDOMTOT = 'MODELISER' DOMTOT MDNS ; MDTHIRD = 'MODELISER' DTHIRD MDNS ; MDSEC = 'MODELISER' DSEC MDNS ; *---------------------------------- $DOMI = 'DOMA' MDOMI 'VF' ; $DREST = 'DOMA' MDREST 'VF' ; $FRONTH = 'DOMA' MFRONTH 'VF' ; $DOMTOT = 'DOMA' MDOMTOT 'VF' ; $DTHIRD = 'DOMA' MDTHIRD 'VF' ; $DSEC = 'DOMA' MDSEC 'VF' ; *---------------------------------- QDOMI = $DOMI . 'QUAF' ; QDREST = $DREST . 'QUAF' ; QFRONTH = $FRONTH . 'QUAF' ; QDOMTOT = $DOMTOT . 'QUAF' ; QDTHIRD = $DTHIRD . 'QUAF' ; QDSEC = $DSEC . 'QUAF' ; *---------------------------------- 'ELIMINATION' QDOMTOT (1.0D-6 '/' RAF) QDOMI ; 'ELIMINATION' QDOMTOT (1.0D-6 '/' RAF) QDREST ; 'ELIMINATION' QDOMTOT (1.0D-6 '/' RAF) QFRONTH ; 'ELIMINATION' QDOMTOT (1.0D-6 '/' RAF) QDTHIRD ; 'ELIMINATION' QDOMTOT (1.0D-6 '/' RAF) QDSEC ; *------------------------------- LIN1 = LIGPROC LIG1 ('DOMA' MDOMTOT 'CENTRE') ; * 'TRACER' (DOMTOT 'ET' LIN1) ; *---------------------------------- * Ligne at the "coridor" *---------------------------------- VECL = ((-0.5 '*' DX1) (0.5 '*' DX1)) ; PL2 = A21 'PLUS' VECL ; PL2 = ('DOMA' MDOMTOT 'CENTRE') 'POIN' 'PROCHE' PL2 ; VECL = ((-0.5 '*' DX1) (-0.5 '*' DX1)) ; PL3 = A31 'PLUS' VECL ; PL3 = ('DOMA' MDOMTOT 'CENTRE') 'POIN' 'PROCHE' PL3 ; LIG2 = PL2 'DROIT' (NY1 '-' 1) PL3 ; 'SI' (RAF 'EGA' 10) ; LIN2 = LIG2 ; 'SINON' ; LIN2 = LIGPROC LIG2 ('DOMA' MDOMTOT 'CENTRE') ; 'FINSI' ; * 'TRACER' (DOMTOT 'ET' LIG2) ; * 'TRACER' (DOMTOT 'ET' LIN2) ; *-------------------------------- MDLIN2 = 'MODELISER' LIN2 MDNS ; $DLIN2 = 'DOMA' MDLIN2 'VF' ; QDLIN2 = $DLIN2 . 'QUAF' ; 'ELIMINATION' QDOMTOT (1.0D-6 '/' RAF) QDLIN2 ; *---------------------------------- * "Centralising the capteurs" *---------------------------------- PCE1 = ('DOMA' MDOMTOT 'CENTRE') 'POIN' 'PROCHE' PCEL1 ; PCE2 = ('DOMA' MDOMTOT 'CENTRE') 'POIN' 'PROCHE' PCEL2 ; PCE3 = ('DOMA' MDOMTOT 'CENTRE') 'POIN' 'PROCHE' PCEL3 ; PCE4 = ('DOMA' MDOMTOT 'CENTRE') 'POIN' 'PROCHE' PCEL4 ; PCE5 = ('DOMA' MDOMTOT 'CENTRE') 'POIN' 'PROCHE' PCEL5 ; PCE6 = ('DOMA' MDOMTOT 'CENTRE') 'POIN' 'PROCHE' PCEL6 ; PCE7 = ('DOMA' MDOMTOT 'CENTRE') 'POIN' 'PROCHE' PCEL7 ; PCE8 = ('DOMA' MDOMTOT 'CENTRE') 'POIN' 'PROCHE' PCEL8 ; PCE9 = ('DOMA' MDOMTOT 'CENTRE') 'POIN' 'PROCHE' PCEL9 ; * 'TRACER' (DOMTOT 'ET' PCE1 'ET' PCE2 'ET' PCE3 'ET' PCE4 * 'ET' PCE5 'ET' PCE6 'ET' PCE7 'ET' PCE8 'ET' PCE9) ; * 'TRACER' (DOMTOT 'ET' ('DOMA' MDOMI 'CENTRE')) ; numbel = 'NBEL' DOMTOT ; *------------------------------------------- * Cells sizes: min, max and averaged *------------------------------------------- VOLC = 'DOMA' MDOMTOT 'VOLUME' ; SIZE = VOLC '**' 0.5 ; *------- minimal SMIN = 'MINIMUM' SIZE ; *------- maximal SMAX = 'MAXIMUM' SIZE ; *------- averaged CHUN = 'MANUEL' 'CHPO' ('DOMA' MDOMTOT 'CENTRE') 1 'SCAL' 1.0 'NATU' 'DISCRET' ; SIZEAV = 'XTY' SIZE CHUN ('MOTS' 'SCAL') ('MOTS' 'SCAL') ; SIZEAV = SIZEAV '/' numbel ; *------------------------------- * Saving the results *------------------------------- *'OPTION' 'SAUVER' ('CHAINE' '/serg/HDR12.3.2_dir/' * 'mail' RAF '.sauv') ; * 'SAUV' ; *----------------------------------------------- * Initial conditions '+' gas properties *----------------------------------------------- * 'OPTION' 'ECHO' 1 'DIME' 2 'ELEM' 'QUA4' * 'TRAC' 'X' ; * * RAF = 10 ; * * 'OPTION' 'REST' ('CHAINE' '/serg/HDR12.3.2_dir/mail' RAF '.sauv') ; * 'RESTITUER' ; zero = 1.0D-9 ; GRAPH = FAUX ; *************************************************************** ***** PROCEDURE POUR CALCULER LES 'CP' ET 'CV' ***** *************************************************************** 'DEBPROC' CPCV TN*'FLOTTANT' PGAZ*'TABLE' ; * CPP = 'TABLE' ; CVV = 'TABLE' ; TT = TN ; T2 = TT '*' TT ; T3 = T2 '*' TT ; T4 = T3 '*' TT ; 'REPETER' BL1 ('DIME' (PGAZ . 'ESPEULE')) ; NOMCOM = 'EXTRAIRE' &BL1 (PGAZ . 'ESPEULE'); A0 = ('EXTRAIRE' 1 (PGAZ . NOMCOM . 'A')) ; A1 = ('EXTRAIRE' 2 (PGAZ . NOMCOM . 'A')) '/' 2.0 ; A2 = ('EXTRAIRE' 3 (PGAZ . NOMCOM . 'A')) '/' 3.0 ; A3 = ('EXTRAIRE' 4 (PGAZ . NOMCOM . 'A')) '/' 4.0 ; A4 = ('EXTRAIRE' 5 (PGAZ . NOMCOM . 'A')) '/' 5.0 ; CVV . NOMCOM = A0 '+' (A1*TT) '+' (A2 '*' T2) '+' (A3 '*' T3) '+' (A4 '*' T4) ; CPP . NOMCOM = (CVV . NOMCOM) '+' (PGAZ . NOMCOM . 'R') ; 'FIN' BL1 ; NOMCOM = PGAZ . 'ESPNEULE' ; A0 = ('EXTRAIRE' 1 (PGAZ . NOMCOM . 'A')) ; A1 = ('EXTRAIRE' 2 (PGAZ . NOMCOM . 'A')) '/' 2.0 ; A2 = ('EXTRAIRE' 3 (PGAZ . NOMCOM . 'A')) '/' 3.0 ; A3 = ('EXTRAIRE' 4 (PGAZ . NOMCOM . 'A')) '/' 4.0 ; A4 = ('EXTRAIRE' 5 (PGAZ . NOMCOM . 'A')) '/' 5.0 ; CVV . NOMCOM = A0 '+' (A1*TT) '+' (A2 '*' T2) '+' (A3 '*' T3) '+' (A4 '*' T4) ; CPP . NOMCOM = (CVV . NOMCOM) '+' (PGAZ . NOMCOM . 'R') ; 'FINPROC' CPP CVV ; *************************************************************** ***** PROCEDURE POUR CALCULER LES VARIABLES CONSERVATIVES ***** *************************************************************** 'DEBPROC' CONS RN*'CHPOINT' VN*'CHPOINT' PN*'CHPOINT' YN*'CHPOINT' $DM*'MMODEL' PGAZ*'TABLE' ; * A0 = 'KCHT' $DM 'SCAL' 'CENTRE' 0.0 ; A1 = 'KCHT' $DM 'SCAL' 'CENTRE' 0.0 ; A2 = 'KCHT' $DM 'SCAL' 'CENTRE' 0.0 ; A3 = 'KCHT' $DM 'SCAL' 'CENTRE' 0.0 ; A4 = 'KCHT' $DM 'SCAL' 'CENTRE' 0.0 ; RTOT = 'KCHT' $DM 'SCAL' 'CENTRE' 0.0 ; YTOT = 'KCHT' $DM 'SCAL' 'CENTRE' 0.0 ; 'REPETER' BL1 ('DIME' (PGAZ . 'ESPEULE')) ; NOMCOM = 'EXTRAIRE' &BL1 (PGAZ . 'ESPEULE'); YCEL = 'EXCO' NOMCOM YN ; RYCEL = 'NOMC' (YCEL * RN) NOMCOM 'NATU' 'DISCRET' ; YTOT = YTOT '+' YCEL ; RTOT = RTOT '+' ((PGAZ . NOMCOM . 'R') * YCEL) ; A0 = A0 '+' (('EXTRAIRE' 1 (PGAZ . NOMCOM . 'A')) '*' YCEL) ; A1 = A1 '+' (('EXTRAIRE' 2 (PGAZ . NOMCOM . 'A')) '*' YCEL) ; A2 = A2 '+' (('EXTRAIRE' 3 (PGAZ . NOMCOM . 'A')) '*' YCEL) ; A3 = A3 '+' (('EXTRAIRE' 4 (PGAZ . NOMCOM . 'A')) '*' YCEL) ; A4 = A4 '+' (('EXTRAIRE' 5 (PGAZ . NOMCOM . 'A')) '*' YCEL) ; 'SI' ('EGA' &BL1 1) ; RYN = 'KCHT' $DM 'SCAL' 'CENTRE' 'COMP' NOMCOM RYCEL ; 'SINON' ; RYN = RYN 'ET' ('KCHT' $DM 'SCAL' 'CENTRE' 'COMP' NOMCOM RYCEL) ; 'FINSI' ; 'FIN' BL1 ; YCEL = ('KCHT' $DM 'SCAL' 'CENTRE' 1.0) '-' YTOT ; NOMCOM = PGAZ . 'ESPNEULE' ; A0 = A0 '+' (('EXTRAIRE' 1 (PGAZ . NOMCOM . 'A')) '*' YCEL) ; A1 = A1 '+' (('EXTRAIRE' 2 (PGAZ . NOMCOM . 'A')) '*' YCEL) ; A2 = A2 '+' (('EXTRAIRE' 3 (PGAZ . NOMCOM . 'A')) '*' YCEL) ; A3 = A3 '+' (('EXTRAIRE' 4 (PGAZ . NOMCOM . 'A')) '*' YCEL) ; A4 = A4 '+' (('EXTRAIRE' 5 (PGAZ . NOMCOM . 'A')) '*' YCEL) ; RTOT = RTOT '+' ((PGAZ . NOMCOM . 'R') '*' YCEL); TN = PN '/' (RN '*' RTOT) ; T2 = TN '*' TN ; T3 = T2 '*' TN ; T4 = T3 '*' TN ; T5 = T4 '*' TN ; ETHER = (A0 * TN) '+' ((A1 '/' 2.0) * T2) '+' ((A2 '/' 3.0) * T3) '+' ((A3 '/' 4.0) * T4) '+' ((A4 '/' 5.0) * T5) ; GN = RN '*' VN ('MOTS' 'SCAL' 'SCAL') ('MOTS' 'UX ' 'UY ') ('MOTS' 'UX ' 'UY ') ; GN = 'KCHT' $DM 'VECT' 'CENTRE' GN ; ECIN = 0.5D0 '*' ('PSCAL' GN VN ('MOTS' 'UX ' 'UY ') ('MOTS' 'UX ' 'UY ')) ; REN = (RN '*' ETHER) '+' ECIN ; REN = 'KCHT' $DM 'SCAL' 'CENTRE' REN ; 'FINPROC' GN REN RYN ; ************************************************************ ************************************************************ * * ***** Gas properties ***** * * ************************************************************ ************************************************************ PGAZ = 'TABLE' ; PGAZB = 'TABLE' ; * Polynomial degree of specific heats PGAZ . 'NORD' = 4 ; PGAZB . 'NORD' = 0 ; * Species explicitly treated in the Euler Equations PGAZ . 'ESPEULE' = 'MOTS' 'H2 ' 'O2 ' 'H2O ' ; PGAZB . 'ESPEULE' = PGAZ . 'ESPEULE' ; * Species non explicitly treated PGAZ . 'ESPNEULE' = 'N2 '; PGAZB . 'ESPNEULE' = PGAZ . 'ESPNEULE' ; * Gas species properties PGAZ . 'H2 ' = 'TABLE' ; PGAZ . 'H2O ' = 'TABLE' ; PGAZ . 'N2 ' = 'TABLE' ; PGAZ . 'O2 ' = 'TABLE' ; *---------------- PGAZB . 'H2 ' = 'TABLE' ; PGAZB . 'H2O ' = 'TABLE' ; PGAZB . 'N2 ' = 'TABLE' ; PGAZB . 'O2 ' = 'TABLE' ; * R (J/Kg/K) and molar masses mH2 = 2. '*' 1.00797E-3 ; mO2 = 2. '*' 15.9994E-3 ; mH2O = mh2 '+' (0.5 '*' mo2) ; mN2 = 2 '*' 14.0067E-3 ; RGAS = 8.31441 ; *------------------------------------- PGAZ . 'H2 ' . 'R' = RGAS '/' mh2 ; PGAZ . 'H2O ' . 'R' = RGAS '/' mh2o ; PGAZ . 'N2 ' . 'R' = RGAS '/' mn2 ; PGAZ . 'O2 ' . 'R' = RGAS '/' mo2 ; *------------------------------------- PGAZB . 'H2 ' . 'R' = RGAS '/' mh2 ; PGAZB . 'H2O ' . 'R' = RGAS '/' mh2o ; PGAZB . 'N2 ' . 'R' = RGAS '/' mn2 ; PGAZB . 'O2 ' . 'R' = RGAS '/' mo2 ; *------------------------------------- * Polynomials regressions coefficients PGAZ . 'H2 ' . 'A' = 'PROG' 9834.91866 0.54273926 0.000862203836 -2.37281455E-07 1.84701105E-11 ; PGAZ . 'H2O ' . 'A' = 'PROG' 1155.95625 0.768331151 -5.73129958E-05 -1.82753232E-08 2.44485692E-12 ; PGAZ . 'N2 ' . 'A' = 'PROG' 652.940766 0.288239099 -7.80442298E-05 8.78233606E-09 -3.05514485E-13 ; PGAZ . 'O2 ' . 'A' = 'PROG' 575.012333 0.350522002 -0.000128294865 2.33636971E-08 -1.53304905E-12; *------------------------------------------------------- * Formation enthalpies at 0K (J/Kg)) * * h_i(0K) = h_i(T0) '-' \int_0^{T0} cp_i(x) dx * = h_i(T0) '-' (\int_0^{T0} cv_i(x) dx '+' R_i * T0) *------------------------------------------------------- PGAZ . 'H2 ' . 'H0K' = -4.195D6 ; PGAZ . 'H2O ' . 'H0K' = -1.395D7 ; PGAZ . 'N2 ' . 'H0K' = -2.953D5 ; PGAZ . 'O2 ' . 'H0K' = -2.634D5 ; *------------------------------------- PGAZB . 'H2 ' . 'H0K' = -4.195D6 ; PGAZB . 'H2O ' . 'H0K' = -1.395D7 ; PGAZB . 'N2 ' . 'H0K' = -2.953D5 ; PGAZB . 'O2 ' . 'H0K' = -2.634D5 ; *----------------------------------- *** Names of the passive scalars *----------------------------------- PGAZ . 'SCALPASS' = 'MOTS' 'H2IN' 'H2FI' 'K0 ' ; *----------------------------------------------- PGAZB . 'SCALPASS' = 'MOTS' 'H2IN' 'H2FI' 'K0 ' ; **************************************************** **************************************************** * Initial conditions * **************************************************** **************************************************** *** inside the domain tg = 337.0D0 ; pg = 100000.0D0 ; xh2o = 0.26D0 ; xh2 = 0.1D0 ; xo2 = 0.1344D0 ; xn2 = 1.0D0 - (xh2 + xo2 + xh2O) ; uxg = 0.0 ; uyg = 0.0 ; *------------------------- td = 337.0D0 ; uxd = 0.0D0 ; uyd = 0.0D0 ; pd = 100000.0D0 ; *---- values of K_0 K0IN = 0.41D0 ; K0IG = 0.41D0 ; K03 = 3.7D0 ; *-------------------------------------------- * After burning: AIBC (same pressure) *-------------------------------------------- tc = 1086.38864D0 ; xh2c = zero ; xo2c = 0.0888421D0 ; xh2oc = 0.3789469D0 ; xn2c = 1.0D0 '-' (xh2c '+' xo2c '+' xh2oc) ; pc = pg ; *--------------------------------------------- CP1 CV1 = CPCV tg PGAZ ; *'-'*'-'*'-'*'-'*'-'*- PGAZB . 'CP' = 'TABLE' ; PGAZB . 'CP' . 'H2 ' = CP1 . 'H2' ; PGAZB . 'CP' . 'H2O ' = CP1 . 'H2O' ; PGAZB . 'CP' . 'O2 ' = CP1 . 'O2' ; PGAZB . 'CP' . 'N2 ' = CP1 . 'N2' ; *---------------------------------------- PGAZB . 'CV' = 'TABLE' ; PGAZB . 'CV' . 'H2 ' = CV1 . 'H2' ; PGAZB . 'CV' . 'H2O ' = CV1 . 'H2O' ; PGAZB . 'CV' . 'O2 ' = CV1 . 'O2' ; PGAZB . 'CV' . 'N2 ' = CV1 . 'N2' ; *--------------------------------------------- * Creating CHAMPOINTS *--------------------------------------------- *** Pressure MAILC = 'DOMA' MDOMTOT 'CENTRE' ; MAILCI = 'DOMA' MDOMI 'CENTRE' ; MAILCR = 'DOMA' MDREST 'CENTRE' ; MAILC2 = 'DOMA' MDSEC 'CENTRE' ; MAILC3 = 'DOMA' MDTHIRD 'CENTRE' ; ********************************* PINI = ('MANUEL' 'CHPO' MAILCI 1 'SCAL' pc 'NATU' 'DISCRET') 'ET' ('MANUEL' 'CHPO' MAILCR 1 'SCAL' pg 'NATU' 'DISCRET') ; *** Temperature TINI = ('MANUEL' 'CHPO' MAILCI 1 'SCAL' tc 'NATU' 'DISCRET') 'ET' ('MANUEL' 'CHPO' MAILCR 1 'SCAL' tg 'NATU' 'DISCRET') ; **** Velocity WINI = 'MANUEL' 'CHPO' MAILC 2 'UX' uxg 'UY' uyg ; *------------------------------------------------------- * Combustion is complete, i.e. \ksi = 1 after combustion *------------------------------------------------------- CSIMAX = 'MANUEL' 'CHPO' MAILC 1 'SCAL' 1.0D0 ; **** Molar fractions before combustion XH21 = ('MANUEL' 'CHPO' MAILCI 1 'H2' xh2c 'NATU' 'DISCRET') 'ET' ('MANUEL' 'CHPO' MAILCR 1 'H2' xh2 'NATU' 'DISCRET') ; XH2O1 = ('MANUEL' 'CHPO' MAILCI 1 'H2O' xh2oc 'NATU' 'DISCRET') 'ET' ('MANUEL' 'CHPO' MAILCR 1 'H2O' xh2o 'NATU' 'DISCRET') ; XO21 = ('MANUEL' 'CHPO' MAILCI 1 'O2' xo2c 'NATU' 'DISCRET') 'ET' ('MANUEL' 'CHPO' MAILCR 1 'O2' xo2 'NATU' 'DISCRET') ; XN21 = ('MANUEL' 'CHPO' MAILCI 1 'N2' xn2c 'NATU' 'DISCRET') 'ET' ('MANUEL' 'CHPO' MAILCR 1 'N2' xn2 'NATU' 'DISCRET') ; XN = XH21 'ET' XH2O1 'ET' XO21 'ET' XN21 ; **** Mass Fractions MMOL = ('MANUEL' 'CHPO' ('DOMA' MDOMTOT 'CENTRE') 1 'H2' mH2 'NATURE' 'DISCRET') 'ET' ('MANUEL' 'CHPO' ('DOMA' MDOMTOT 'CENTRE') 1 'O2' mO2 'NATURE' 'DISCRET') 'ET' ('MANUEL' 'CHPO' ('DOMA' MDOMTOT 'CENTRE') 1 'H2O' mH2O 'NATURE' 'DISCRET') 'ET' ('MANUEL' 'CHPO' ('DOMA' MDOMTOT 'CENTRE') 1 'N2' mN2 'NATURE' 'DISCRET') ; MTOT = 'PSCAL' MMOL XN ('MOTS' 'H2' 'O2' 'H2O' 'N2') ('MOTS' 'H2' 'O2' 'H2O' 'N2') ; YH21 = (XH21 '*' mH2) '/' MTOT ; YO21 = (XO21 '*' mO2) '/' MTOT ; YN21 = (XN21 '*' mN2) '/' MTOT ; YH2O1 = (XH2O1 '*' mH2O) '/' MTOT ; aaa = 'MAXIMUM' (YH21) ; *---------------------------------------------- * Burned gas * First we find how much of the burned gas * produces using 1 mole of unburned gas *---------------------------------------------- CHUN = 'MANUEL' 'CHPO' MAILC 1 'SCAL' 1.0D0 'NATU' 'DISCRET' ; H2MAG = (0.5 '*' ('EXCO' 'H2' XN)) 'MASQUE' 'SUPERIEUR' ('EXCO' 'O2' XN) ; O2MAG = CHUN '-' H2MAG ; XO22a = 'MANUEL' 'CHPO' MAILC 1 'O2' zero ; XH22a = XH21 '-' ('NOMC' 'H2' (2.0 '*' (XO21 '-' XO22a))) ; XO22a = XO22a ; XH22a = XH22a ; XH22b = 'MANUEL' 'CHPO' MAILC 1 'H2' zero ; XO22b = XO21 '-' ('NOMC' 'O2' (0.5 '*' (XH21 '-' XH22b))) ; XH22 = (XH22a '*' H2MAG) '+' (XH22b '*' O2MAG) ; XO22 = (XO22a '*' H2MAG) '+' (XO22b '*' O2MAG) ; XH2O2 = XH2O1 '+' ('NOMC' 'H2O' (2 '*' (XO21 '-' XO22))) ; XN22 = XN21 ; *--------------------------------------------------- * Verifying that mtot2 = mtot *--------------------------------------------------- MTOT2 = 'PSCAL' MMOL (XH22 '+' XO22 '+' XH2O2 '+' XN22) ('MOTS' 'H2' 'O2' 'H2O' 'N2') ('MOTS' 'H2' 'O2' 'H2O' 'N2') ; 'SI' (('MAXIMUM' (mtot2 '-' mtot) 'ABS') '>' zero) ; 'MESSAGE' 'problem1' ; 'FINSI' ; ***** xtot2 = mole of burned gas per mole of unburned gas xtot2 = 'PSCAL' (XH22 '+' XO22 '+' XH2O2 '+' XN22) CHUN ('MOTS' 'H2' 'O2' 'H2O' 'N2') ('MOTS' 'SCAL' 'SCAL' 'SCAL' 'SCAL') ; ***** X now is the real molar fractions XH22 = XH22 '/' xtot2 ; XO22 = XO22 '/' xtot2 ; XH2O2 = XH2O2 '/' xtot2 ; XN22 = XN22 '/' xtot2 ; xtot3 = 'PSCAL' (XH22 '+' XO22 '+' XH2O2 '+' XN22) CHUN ('MOTS' 'H2' 'O2' 'H2O' 'N2') ('MOTS' 'SCAL' 'SCAL' 'SCAL' 'SCAL') ; 'SI' (('MAXIMUM' (xtot3 '-' 1) 'ABS') '>' zero) ; 'MESSAGE' 'problem2' ; 'FINSI' ; **** mtot2 = weight of 1 mole of burned gas MTOT2 = 'PSCAL' MMOL (XH22 '+' XO22 '+' XH2O2 '+' XN22) ('MOTS' 'H2' 'O2' 'H2O' 'N2') ('MOTS' 'H2' 'O2' 'H2O' 'N2') ; **** Mass fractions of burned gas after complete combustion YH22 = (XH22 '*' MH2) '/' MTOT2 ; YO22 = (XO22 '*' MO2) '/' MTOT2 ; YN22 = (XN22 '*' MN2) '/' MTOT2 ; YH2O2 = (XH2O2 '*' MH2O) '/' MTOT2 ; *--------------------------------------------- * Conservative variables *--------------------------------------------- RINI = (('NOMC' 'SCAL' yh21) '*' (PGAZ . 'H2' . 'R')) '+' (('NOMC' 'SCAL' yo21) '*' (PGAZ . 'O2' . 'R')) '+' (('NOMC' 'SCAL' yh2o1) '*' (PGAZ . 'H2O' . 'R')) '+' (('NOMC' 'SCAL' yn21) '*' (PGAZ . 'N2' . 'R')) ; RN = PINI '/' (RINI '*' TINI) ; YN = YH21 '+' YO21 '+' YH2O1 ; GN RETN RYN = CONS RN WINI PINI YN MDOMTOT PGAZ ; *-------------------------------------- * Passive scalars *-------------------------------------- YININ = 'MANUEL' 'CHPO' MAILC 1 'H2' aaa 'NATU' 'DISCRET' ; YFINN = 'MANUEL' 'CHPO' MAILC 1 'H2' zero 'NATU' 'DISCRET' ; K0N = ('MANUEL' 'CHPO' MAILCI 1 'SCAL' K0IG 'NATU' 'DISCRET') 'ET' ('MANUEL' 'CHPO' MAILC2 1 'SCAL' K0IN 'NATU' 'DISCRET') 'ET' ('MANUEL' 'CHPO' MAILC3 1 'SCAL' K03 'NATU' 'DISCRET') ; *--------------------------------------- RSN = RN '*' (('NOMC' 'H2IN' YININ 'NATU' 'DISCRET') 'ET' ('NOMC' 'H2FI' YFINN 'NATU' 'DISCRET') 'ET' ('NOMC' 'K0' K0N 'NATU' 'DISCRET')) ; *------------------- ***** Verification *------------------- VN P T Y S GAMN = 'PRIM' 'PERFTEMP' PGAZ RN GN RETN RYN RSN ; *----------------------- PCEL = 'MAXIMUM' PINI ; TCEL = 'MAXIMUM' TINI ; 'SI'(('MAXIMUM' ((P '-' PINI) '/' PCEL) 'ABS') '>' (zero '*' 10.0)) ; 'MESSAGE' 'problem3' ; 'LISTE' ('MAXIMUM' ((P '-' PINI) '/' PCEL) 'ABS') ; 'FINSI' ; 'SI'(('MAXIMUM' ((T '-' TINI) '/' TCEL) 'ABS') '>' (zero '*' 10.0)) ; 'MESSAGE' 'problem4' ; 'LISTE' ('MAXIMUM' ((T '-' TINI) '/' TCEL) 'ABS') ; 'FINSI' ; 'SI' (('MAXIMUM' VN 'ABS') '>' zero) ; 'MESSAGE' 'problem5' ; 'LISTE' ('MAXIMUM' VN 'ABS') ; 'FINSI' ; 'SI'(('MAXIMUM' (Y '-' YN) 'ABS') '>' zero) ; 'MESSAGE' 'problem6' ; 'LISTE' ('MAXIMUM' (Y '-' YN) 'ABS') ; 'FINSI' ; *-------------------------------------- * Graphics *-------------------------------------- 'SI' GRAPH ; CHM_TN = 'KCHA' MDOMTOT 'CHAM' TINI ; CHM_RN = 'KCHA' MDOMTOT 'CHAM' RN ; CHM_PN = 'KCHA' MDOMTOT 'CHAM' PINI ; CHM_VN = 'KCHA' MDOMTOT 'CHAM' WINI ; CHM_H2 = 'KCHA' MDOMTOT 'CHAM' YN ; CHM_YIN = 'KCHA' MDOMTOT 'CHAM' YININ ; CHM_YFI = 'KCHA' MDOMTOT 'CHAM' YFINN ; CHM_K0 = 'KCHA' MDOMTOT 'CHAM' K0N ; 'TRACER' CHM_TN MDOMTOT 'TITR' ('CHAINE' 'Temperature at t=' 0.0) ; 'TRACER' CHM_RN MDOMTOT 'TITR' ('CHAINE' 'RN at t=' 0.0) ; 'TRACER' CHM_PN MDOMTOT 'TITR' ('CHAINE' 'PN at t=' 0.0) ; 'TRACER' CHM_VN MDOMTOT 'TITR' ('CHAINE' 'VN at t=' 0.0) ; 'TRACER' CHM_H2 MDOMTOT 'TITR' ('CHAINE' 'H2 at t=' 0.0) ; 'TRACER' CHM_YIN MDOMTOT 'TITR' ('CHAINE' 'YIN at t=' 0.0) ; 'TRACER' CHM_YFI MDOMTOT 'TITR' ('CHAINE' 'YFI at t=' 0.0) ; 'TRACER' CHM_K0 MDOMTOT 'TITR' ('CHAINE' 'K0 at t=' 0.0) ; 'FINSI' ; *----------------------------------------------- * Parameters for the computation *----------------------------------------------- *** Upwind scheme METO = 'SS' ; **** Iterations NITER = 100 ; **** Final time (no restriction) TFINAL = 10000.1 ; **** Safety factor for the time step SAFFAC = 0.7 ; *** Second order ? LOGSO = VRAI ; LOGSO = FAUX ; *** EPS_CC '-' epsilon in CREBCOM criterion EPS_CC = 0.5 ; **** "explosive" initial conditions? LOGEXP = VRAI ; LOGEXP = FAUX ; *** Are there any Boundary conditions LOGBC = VRAI ; *------------------------------- * Time intervals for the output *------------------------------- LISTT = 'PROG' 4.5 ; *** Mesh size DELTAXN = ('DOMA' MDOMTOT 'VOLUME') '**' 0.5 ; DELX = 'MAXIMUM' DELTAXN ; *------------------------------------- TABLIM = 'TABLE' ; pwin = (('MAXIMUM' RN) '*' 9.8) '*' 4.668 ; TABLIM . 'CHPOUT' = 'MANUEL' 'CHPO' ('DOMA' MFRONTH 'CENTRE') 1 'PN' (pg '-' pwin) ; TABLIM . 'H' = 50.0D0 ; TABLIM . 'TNUL' = 'MANUEL' 'CHPO' ('DOMA' MDOMTOT 'CENTRE') 1 'SCAL' tg ; TABLIM . 'FREQ' = 1 ; *----- gravity ---------- TABLIM . 'GRAV' = 'MANUEL' 'CHPO' ('DOMA' MDOMTOT 'CENTRE') 'UX' 0.0D0 'UY' (-9.8D0) ; *---------------------------------------------- * Table for evolution of the variables with time * at the fixed points *---------------------------------------------- TABC = 'TABLE' 'CAPTEUR' ; TABC . 'LPOINTS' = PCEL1 'ET' PCEL2 'ET' PCEL3 'ET' PCEL4 'ET' PCEL5 'ET' PCEL6 'ET' PCEL7 'ET' PCEL8 'ET' PCEL9 ; TABC . 'LTPS' = 'PROG'; TABC . 'NOMC' = 'MOTS' 'C1' 'C2' 'C3' 'C4' 'C5' 'C6' 'C7' 'C8' 'C9' ; TABC . 'COMP' = 'MOTS' 'PN' 'TN' ; TABC . 'FREQ' = 10 ; MSLIN2 = 'DOMA' MDLIN2 'SOMMET' ; bnb = 'NBEL' MSLIN2 ; VOLGOR = DX1 '*' DX1 '*' bnb ; *----------------------------------------------- * In MDOM1 we 'burn' the gas; * In order to do so we shall apply the * operator 'FLAM' in MDOM1 with eps=0 * and DELTAT >> DELTATC *---------------------------------------------- YN = RYN '/' RN ; SN = RSN '/' RN ; K0N = 'EXCO' 'K0' SN ; YININ = 'EXCO' 'H2IN' SN 'H2' ; YFINN = 'EXCO' 'H2FI' SN 'H2' ; DELTATC = 0.25 '*' ('MINIMUM' (DELTAXN '/' K0N)) ; LMOT1 = 'MOTS' 'H2' 'O2' 'H2O' ; **** Stocheom. coefficients LCOEF = 'PROG' 1.0 0.5 -1.0 ; *------------------------------------------- * This to activate in case of constant volume * combustion *------------------------------------------- 'SI' LOGEXP ; DELTARE DELTARY = 'FLAM' 'CREBCOM2' MDOMI PGAZ LMOT1 LCOEF ('REDU' RN MAILCI) ('REDU' YN MAILCI) ('REDU' YININ MAILCI) ('REDU' YFINN MAILCI) ('REDU' K0N MAILCI) ('REDU' DELTAXN MAILCI) 0.0 (1D10 '*' DELTATC) zero ; RYN = RYN '+' DELTARY ; YN = RYN '/' RN ; RETN = RETN '+' DELTARE ; 'FINSI' ; *----------------------------------------------- * Names of the components *----------------------------------------------- LISTINCO = 'MOTS' 'RN' 'RUX' 'RUY' 'RETN' 'H2' 'O2' 'H2O' 'H2IN' 'H2FI' 'K0' ; ****** names of the primitive variables LISTINCP = 'MOTS' 'RN' 'UX' 'UY' 'PN' 'H2' 'O2' 'H2O' 'H2IN' 'H2FI' 'K0' ; *-------------------------------------------- LMOT111 = 'MOTS' 'SCAL' ; GRADR CCC COEFSCAL = 'PENT' MDOMTOT 'CENTRE' 'EULESCAL' 'NOLIMITE' LMOT111 RN ; LMOT2 = 'MOTS' 'UX' 'UY' ; GRADV CCC COEFVECT = 'PENT' MDOMTOT 'CENTRE' 'EULEVECT' 'NOLIMITE' LMOT2 GN ; *-------------------------------------------- * Volumes of the rooms ; *-------------------------------------------- VV2 = 'DOMA' MDSEC 'VOLUME' ; CHUN = ('MANUEL' 'CHPO' MAILC2 1 'SCAL' 1.0 'NATU' 'DISCRET') ; VV2S = ('XTY' CHUN VV2 ('MOTS' 'SCAL') ('MOTS' 'SCAL')) ; *------------------------------- VV3 = 'DOMA' MDTHIRD 'VOLUME' ; CHUN = ('MANUEL' 'CHPO' MAILC3 1 'SCAL' 1.0 'NATU' 'DISCRET') ; VV3S = ('XTY' CHUN VV3 ('MOTS' 'SCAL') ('MOTS' 'SCAL')) ; *------------------------------------------- LISTPS2 = 'PROG' ; LISTPS3 = 'PROG' ; LISTTS2 = 'PROG' ; LISTTS2 = 'PROG' ; *---- 'LISTE' for the discharge at the coridor LISCOR1 = 'PROG' ; *---- 'LISTE' for the averaged velocity at the coridor LISCOR2 = 'PROG' ; *---- 'LISTE' for the discharge at the vent LISVEN1 = 'PROG' ; *---- 'LISTE' for the averaged velocity at the vent LISVEN2 = 'PROG' ; *-------------------------------------------- TPS = 0.0 ; *------------------------------------------ * Chemical time step *------------------------------------------ DT_CHEM = SAFFAC '*' DELTATC ; *------------------------------------------- **** Temporal loop *------------------------------------------- 'MESSAGE' ; 'MESSAGE' ('CHAINE' 'Methode = ' METO) ; 'MESSAGE' ('CHAINE' 'SAFFAC = ' SAFFAC) ; 'MESSAGE' ; *----------------------------------- 'REPETER' BLG ('DIME' LISTT) ; *---------------------------------- TFINAL = 'EXTRAIRE' &BLG LISTT ; *---------------------------------- VN PN TN YN SN GAMMAN = 'PRIM' 'PERFTEMP' PGAZ RN GN RETN RYN RSN ; TNM1 = 'COPIER' TN ; 'TEMPS' 'ZERO' ; 'REPETER' BL1 NITER ; * **** Primitive variables * FREQC = TABC.'FREQ' ; 'SI' (((&BL1 '/' FREQC) '*' FREQC) 'EGA' &BL1) ; ******* The "coridor" ********************** RL2 = 'REDU' RN MSLIN2 ; GL2 = 'REDU' GN MSLIN2 ; GLX2 = 'EXCO' 'UX' GL2 ; CHUN = ('MANUEL' 'CHPO' MSLIN2 1 'SCAL' (DX1) 'NATU' 'DISCRET') ; *------- rashod discharge ---------------- GLXAV = ('XTY' CHUN GLX2 ('MOTS' 'SCAL') ('MOTS' 'SCAL')) ; GLXAV = GLXAV '/' H1 ; LISCOR1 = LISCOR1 'ET' ('PROG' GLXAV) ; *---- redifining CHUN ------------------- CHUN = ('MANUEL' 'CHPO' MSLIN2 1 'SCAL' (DX1 '*' DX1) 'NATU' 'DISCRET') ; *------- averaged density ---------------- RLAV = ('XTY' CHUN RL2 ('MOTS' 'SCAL') ('MOTS' 'SCAL')) ; RLAV = RLAV '/' VOLGOR ; *------- averaged momentum --------------- GLXAV = ('XTY' CHUN GLX2 ('MOTS' 'SCAL') ('MOTS' 'SCAL')) ; GLXAV = GLXAV '/' VOLGOR ; *----------------------------------------- * UXL = GLXAV '/' RLAV ; CHVIT = GLX2 '*' ('INVERSE' RL2) ; UXL = 'MAXIMUM' CHVIT ; LISCOR2 = LISCOR2 'ET' ('PROG' UXL) ; ******* The second room ********************** *---------------------------------------------------- * Computation of averaged pressure and temperature *---------------------------------------------------- RN2 = 'REDU' RN MAILC2 ; GN2 = 'REDU' GN MAILC2 ; RETN2 = 'REDU' RETN MAILC2 ; RYN2 = 'REDU' RYN MAILC2 ; RSN2 = 'REDU' RSN MAILC2 ; *--- we average the conservative data RN2S = ('XTY' RN2 VV2 ('MOTS' 'SCAL') ('MOTS' 'SCAL')) ; RN2S = RN2S '/' VV2S ; *---------------------- RETN2S = ('XTY' RETN2 VV2 ('MOTS' 'SCAL') ('MOTS' 'SCAL')) ; RETN2S = RETN2S '/' VV2S ; *----------------------- GNX2 = 'EXCO' 'UX' GN2 'UX' ; GNY2 = 'EXCO' 'UY' GN2 'UY' ; *------------------- GNX2S = ('XTY' GNX2 VV2 ('MOTS' 'SCAL') ('MOTS' 'SCAL')) ; GNX2S = GNX2S '/' VV2S ; *------------------- GNY2S = ('XTY' GNY2 VV2 ('MOTS' 'SCAL') ('MOTS' 'SCAL')) ; GNY2S = GNY2S '/' VV2S ; *-------------------- RYNH2 = 'EXCO' 'H2' RYN2 'H2' ; RYNO2 = 'EXCO' 'O2' RYN2 'O2' ; RYNH2O = 'EXCO' 'H2O' RYN2 'H2O' ; *------------------- RYH2S = ('XTY' RYNH2 VV2 ('MOTS' 'H2') ('MOTS' 'SCAL')) ; RYH2S = RYH2S '/' VV2S ; *------------------- RYO2S = ('XTY' RYNO2 VV2 ('MOTS' 'O2') ('MOTS' 'SCAL')) ; RYO2S = RYO2S '/' VV2S ; *------------------- RYHO2S = ('XTY' RYNH2O VV2 ('MOTS' 'H2O') ('MOTS' 'SCAL')) ; RYHO2S = RYHO2S '/' VV2S ; *-------------------------------------- *-------------------------------------- R2S = ('MANUEL' 'CHPO' MAILCI 1 'SCAL' RN2S 'NATU' 'DISCRET') ; G2S = ('MANUEL' 'CHPO' MAILCI 1 'UX' GNX2S 'NATU' 'DISCRET') 'ET' ('MANUEL' 'CHPO' MAILCI 1 'UY' GNY2S 'NATU' 'DISCRET') ; RET2S = ('MANUEL' 'CHPO' MAILCI 1 'SCAL' RETN2S 'NATU' 'DISCRET') ; RY2S = ('MANUEL' 'CHPO' MAILCI 1 'H2' RYH2S 'NATU' 'DISCRET') 'ET' ('MANUEL' 'CHPO' MAILCI 1 'O2' RYO2S 'NATU' 'DISCRET') 'ET' ('MANUEL' 'CHPO' MAILCI 1 'H2O' RYHO2S 'NATU' 'DISCRET') ; RS2S = ('MANUEL' 'CHPO' MAILCI 1 'H2IN' (aaa) 'NATU' 'DISCRET') 'ET' ('MANUEL' 'CHPO' MAILCI 1 'H2FI' zero 'NATU' 'DISCRET') 'ET' ('MANUEL' 'CHPO' MAILCI 1 'K0' K0IN 'NATU' 'DISCRET') ; *-------------------------------------------- VS PS TS YS SS GAMS = 'PRIM' 'PERFTEMP' PGAZ R2S G2S RET2S RY2S RS2S ; PN2S = 'MAXIMUM' PS ; TN2S = 'MAXIMUM' TS ; *------------------------------- LISTPS2 = LISTPS2 'ET' ('PROG' PN2S) ; LISTTS2 = LISTTS2 'ET' ('PROG' TN2S) ; *------------------------------------------------- ******** IN THE BIGGEST ROOM ******************** RN2 = 'REDU' RN MAILC3 ; GN2 = 'REDU' GN MAILC3 ; RETN2 = 'REDU' RETN MAILC3 ; RYN2 = 'REDU' RYN MAILC3 ; RSN2 = 'REDU' RSN MAILC3 ; *--- we average the conservative data RN2S = ('XTY' RN2 VV3 ('MOTS' 'SCAL') ('MOTS' 'SCAL')) ; RN2S = RN2S '/' VV3S ; *---------------------- RETN2S = ('XTY' RETN2 VV3 ('MOTS' 'SCAL') ('MOTS' 'SCAL')) ; RETN2S = RETN2S '/' VV3S ; *----------------------- GNX2 = 'EXCO' 'UX' GN2 'UX' ; GNY2 = 'EXCO' 'UY' GN2 'UY' ; *------------------- GNX2S = ('XTY' GNX2 VV3 ('MOTS' 'SCAL') ('MOTS' 'SCAL')) ; GNX2S = GNX2S '/' VV3S ; *------------------- GNY2S = ('XTY' GNY2 VV3 ('MOTS' 'SCAL') ('MOTS' 'SCAL')) ; GNY2S = GNY2S '/' VV3S ; *-------------------- RYNH2 = 'EXCO' 'H2' RYN2 'H2' ; RYNO2 = 'EXCO' 'O2' RYN2 'O2' ; RYNH2O = 'EXCO' 'H2O' RYN2 'H2O' ; *------------------- RYH2S = ('XTY' RYNH2 VV3 ('MOTS' 'H2') ('MOTS' 'SCAL')) ; RYH2S = RYH2S '/' VV3S ; *------------------- RYO2S = ('XTY' RYNO2 VV3 ('MOTS' 'O2') ('MOTS' 'SCAL')) ; RYO2S = RYO2S '/' VV3S ; *------------------- RYHO2S = ('XTY' RYNH2O VV3 ('MOTS' 'H2O') ('MOTS' 'SCAL')) ; RYHO2S = RYHO2S '/' VV3S ; *-------------------------------------- *-------------------------------------- R2S = ('MANUEL' 'CHPO' MAILCI 1 'SCAL' RN2S 'NATU' 'DISCRET') ; G2S = ('MANUEL' 'CHPO' MAILCI 1 'UX' GNX2S 'NATU' 'DISCRET') 'ET' ('MANUEL' 'CHPO' MAILCI 1 'UY' GNY2S 'NATU' 'DISCRET') ; RET2S = ('MANUEL' 'CHPO' MAILCI 1 'SCAL' RETN2S 'NATU' 'DISCRET') ; RY2S = ('MANUEL' 'CHPO' MAILCI 1 'H2' RYH2S 'NATU' 'DISCRET') 'ET' ('MANUEL' 'CHPO' MAILCI 1 'O2' RYO2S 'NATU' 'DISCRET') 'ET' ('MANUEL' 'CHPO' MAILCI 1 'H2O' RYHO2S 'NATU' 'DISCRET') ; RS2S = ('MANUEL' 'CHPO' MAILCI 1 'H2IN' (aaa) 'NATU' 'DISCRET') 'ET' ('MANUEL' 'CHPO' MAILCI 1 'H2FI' zero 'NATU' 'DISCRET') 'ET' ('MANUEL' 'CHPO' MAILCI 1 'K0' K0IN 'NATU' 'DISCRET') ; *-------------------------------------------- VS PS TS YS SS GAMS = 'PRIM' 'PERFTEMP' PGAZ R2S G2S RET2S RY2S RS2S ; PN2S = 'MAXIMUM' PS ; TN2S = 'MAXIMUM' TS ; *------------------------------------------------- LISTPS3 = LISTPS3 'ET' ('PROG' PN2S) ; LISTTS3 = LISTTS3 'ET' ('PROG' TN2S) ; 'FINSI' ; *-------------------------------------------------- *------------------------------------------- VN PN TN YN SN GAMMAN = 'PRIM' 'PERFTEMP' PGAZ RN GN RETN RYN RSN TNM1 ; *-------------------------------------------- * Boundary conditions *-------------------------------------------- *** find new CP and CV ******* * CP1 CV1 = CPCV tg PGAZ ; *------------------------------------------ 'SI' LOGBC ; 'SI' (&BLG 'EGA' 1) ; SHIN = 'EXCO' SN 'H2IN' ; SHFI = 'EXCO' SN 'H2FI' ; SHK0 = 'EXCO' SN 'K0' ; RCHLIM RCHRES = 'KONV' 'VF' 'PERFMULT' 'CLIM' 'RESI' MDOMTOT MFRONTH PGAZB LISTINCO LISTINCP RN VN PN YN SHIN SHFI SHK0 (TABLIM . 'CHPOUT') 'OUTP' ; 'FINSI' ; *------------------------------------ MAILLIM = 'EXTRAIRE' RCHLIM 'MAILLAGE' ; MAIL111 = 'EXTRAIRE' RCHRES 'MAILLAGE' ; *------------------------------------ NBB = 'NBEL' MAIL111 ; TTT = 'REDU' TN MAIL111 ; UNIF = 'MANUEL' CHPO ('DOMA' MDOMTOT 'CENTRE') 1 'SCAL' 1.0 ; UNI1 = 'REDU' UNIF MAIL111 ; NUNI = UNI1 '*' (1.0 '/' NBB) ; **** computing average temperature at the boundary TAV = 'XTY' TTT NUNI ('MOTS' 'SCAL') ('MOTS' 'SCAL') ; *** find new CP and CV ******* CPN CVN = CPCV TAV PGAZ ; *** Updating CP and CV for the boundary conditions PGAZB . 'CP' . 'H2 ' = CPN . 'H2' ; PGAZB . 'CP' . 'H2O ' = CPN . 'H2O' ; PGAZB . 'CP' . 'O2 ' = CPN . 'O2' ; PGAZB . 'CP' . 'N2 ' = CPN . 'N2' ; *---------------------------------------- PGAZB . 'CV' . 'H2 ' = CVN . 'H2' ; PGAZB . 'CV' . 'H2O ' = CVN . 'H2O' ; PGAZB . 'CV' . 'O2 ' = CVN . 'O2' ; PGAZB . 'CV' . 'N2 ' = CVN . 'N2' ; * 'FINSI' ; *-------------------------------- SHIN = 'EXCO' SN 'H2IN' ; SHFI = 'EXCO' SN 'H2FI' ; SHK0 = 'EXCO' SN 'K0' ; RCHLIM RCHRES = 'KONV' 'VF' 'PERFMULT' 'CLIM' 'RESI' MDOMTOT MFRONTH PGAZB LISTINCO LISTINCP RN VN PN YN SHIN SHFI SHK0 (TABLIM . 'CHPOUT') 'OUTP' ; MAILLIM = 'EXTRAIRE' RCHLIM 'MAILLAGE' ; MAIL111 = 'EXTRAIRE' RCHRES 'MAILLAGE' ; *-------------------------------------------- * We compute the debit et vitesse *-------------------------------------------- 'SI' (((&BL1 '/' FREQC) '*' FREQC) 'EGA' &BL1) ; VOLUP = 'REDU' ('DOMA' MDOMTOT 'VOLUME') MAIL111 ; RLIM1 = 'EXCO' 'RN' RCHRES ; DEBUP = 'XTY' VOLUP RLIM1 ('MOTS' 'SCAL') ('MOTS' 'SCAL') ; DEBUP = (DEBUP '/' LW) '*' (-1.0) ; LISVEN1 = LISVEN1 'ET' ('PROG' DEBUP) ; *-------------------------------------------- RNUP = 'REDU' RN MAIL111 ; GNUP = 'REDU' GN MAIL111 ; GXUP = 'EXCO' 'UY' GNUP ; CUN = ('MANUEL' 'CHPO' MAIL111 1 'SCAL' 1.0D0 'NATU' 'DISCRET') ; *------------------------ GXUPS = 'XTY' GXUP VOLUP ('MOTS' 'SCAL') ('MOTS' 'SCAL') ; RNUPS = 'XTY' RNUP VOLUP ('MOTS' 'SCAL') ('MOTS' 'SCAL') ; VITUP = GXUPS '/' RNUPS ; LISVEN2 = LISVEN2 'ET' ('PROG' VITUP) ; 'FINSI' ; *-------------------------------------------- 'FINSI' ; *------------------------------------ *-------------------------------- TNM1 = 'COPIER' TN ; 'SI' LOGSO ; GRADR ALR = 'PENT' MDOMTOT 'CENTRE' 'EULESCAL' 'LIMITEUR' ('MOTS' 'SCAL') RN 'GRADGEO' COEFSCAL ; GRADP ALP = 'PENT' MDOMTOT 'CENTRE' 'EULESCAL' 'LIMITEUR' ('MOTS' 'SCAL') PN 'GRADGEO' COEFSCAL ; GRADV ALV = 'PENT' MDOMTOT 'CENTRE' 'EULEVECT' 'LIMITEUR' ('MOTS' 'UX' 'UY') VN 'GRADGEO' COEFVECT ; GRADY ALY = 'PENT' MDOMTOT 'CENTRE' 'EULESCAL' 'LIMITEUR' ('MOTS' 'H2' 'O2' 'H2O') YN 'GRADGEO' COEFSCAL ; GRADS ALS = 'PENT' MDOMTOT 'CENTRE' 'EULESCAL' 'LIMITEUR' ('MOTS' 'H2IN' 'H2FI' 'K0') SN 'GRADGEO' COEFSCAL ; ROF VITF PF YF SF = 'PRET' 'PERFTEMP' 2 1 MDOMTOT PGAZ RN GRADR ALR VN GRADV ALV PN GRADP ALP YN GRADY ALY SN GRADS ALS ; 'SINON' ; ROF VITF PF YF SF = 'PRET' 'PERFTEMP' 1 1 MDOMTOT PGAZ RN VN PN YN SN ; 'FINSI' ; *----------------------------------------------- 'SI' LOGBC ; RESIDU DELTAT = 'KONV' 'VF' 'PERFTEMP' 'RESI' METO MDOMTOT PGAZ LISTINCO ROF VITF PF YF SF MAILLIM ; RESIDU = RESIDU 'ET' RCHRES ; 'SINON' ; RESIDU DELTAT = 'KONV' 'VF' 'PERFTEMP' 'RESI' METO MDOMTOT PGAZ LISTINCO ROF VITF PF YF SF ; 'FINSI' ; *----------------------------------------------- DT_CON = SAFFAC '*' DELTAT ; * **** The time step linked to tps * DTTPS = (TFINAL '-' TPS) * (1. '+' ZERO) ; * **** Total time step * DTMIN = 'MINIMUM' ('PROG' DT_CON DTTPS DT_CHEM ) ; * **** Increment of the variables (convection) * RESIDU = DTMIN '*' RESIDU ; DRN = 'EXCO' 'RN' RESIDU 'SCAL' ; DGN = 'EXCO' ('MOTS' 'RUX' 'RUY') RESIDU ('MOTS' 'UX' 'UY') ; DRETN = 'EXCO' 'RETN' RESIDU 'SCAL' ; DRYN = 'EXCO' ('MOTS' 'H2' 'O2' 'H2O') RESIDU ('MOTS' 'H2' 'O2' 'H2O') ; DRSN = 'EXCO' ('MOTS' 'H2IN' 'H2FI' 'K0') RESIDU ('MOTS' 'H2IN' 'H2FI' 'K0') ; TPS = TPS '+' DTMIN ; RN = RN '+' DRN ; GN = GN '+' DGN ; RETN = RETN '+' DRETN ; RYN = RYN '+' DRYN ; RSN = RSN '+' DRSN ; YN = RYN '/' RN ; SN = RSN '/' RN ; * **** Increment of the variables (chemical source) * K0N = 'EXCO' 'K0' SN 'SCAL' ; YININ = 'EXCO' 'H2IN' SN 'H2' ; YFINN = 'EXCO' 'H2FI' SN 'H2' ; DRETN DRYN = 'FLAM' 'CREBCOM2' MDOMTOT PGAZ LMOT1 LCOEF RN YN YININ YFINN K0N DELTAXN EPS_CC DTMIN (zero '*' 1.0D3) ; RYN = RYN '+' DRYN ; RETN = RETN '+' DRETN ; *--------------------------------------------------- * Rate of loosing energy due to boundary conditions *--------------------------------------------------- * 'SI' (TPS '>' 3.0) ; VN PN TN YN SN GAMMAN = 'PRIM' 'PERFTEMP' PGAZ RN GN RETN RYN RSN ; HH = TABLIM . 'H' ; PERT = HH '*' (TN '-' (TABLIM . 'TNUL')) ; RETN = RETN '-' (DTMIN '*' PERT) ; * 'FINSI' ; *-------------------------------------------------- * Gravity forces *-------------------------------------------------- LIMOT1 = 'MOTS' 'SCAL' 'SCAL' ; LIMOT2 = 'MOTS' 'UX' 'UY' ; LIMOT3 = 'MOTS' 'UX' 'UY' ; ROG = RN '*' (TABLIM . 'GRAV') LIMOT1 LIMOT2 LIMOT3 ; GN = GN '+' (DTMIN '*' ROG) ; ROGU = 'PSCAL' ROG VN LIMOT3 LIMOT2 ; RETN = RETN '+' (DTMIN '*' ROGU) ; *--------------------------------------------------- *--------------------------------------------- * Capturing variables at the capteurs *--------------------------------------------- 'SI' (((&BL1 '/' FREQC) '*' FREQC) 'EGA' &BL1) ; CENTOT = 'DOMA' MDOMTOT 'CENTRE'; PTS = TABC . 'LPOINTS' ; ***** M1 is the number of capteurs M1 = 'NBEL' PTS ; LTPS = TABC . 'LTPS' ; NTPS = 'DIME' LTPS ; **** N1 is the number of variables captured L1 = TABC . 'COMP' ; N1 = 'DIME' L1 ; * On etudie uniquement la pression * peut evoluer par la suite !!! LOG1 = 'EXISTE' TABC 'NOMC' ; * Initialisation des tables 'SI' (NTPS EGA 0 ); 'REPETER' BOU1 M1; **** extract the names of the capteurs PT1 = PTS 'POIN' &BOU1; 'SI' (LOG1 'ET' (('DIME' (TABC.'NOMC')) 'EGA' M1 )) ; CAPI = 'EXTR' (TABC.'NOMC') &BOU1; 'SINON' ; CAPI = 'CHAINE' 'CAPT' &BOU1; 'FINSI'; PT2 = CENTOT 'POIN' 'PROCHE' PT1; 'SI' (&BOU1 'EGA' 1); TABC.'RPOINTS' = PT2 ; 'SINON'; TABC.'RPOINTS' = TABC .'RPOINTS' 'ET' PT2; 'FINSI'; TABC . CAPI = 'TABLE'; 'REPETER' BOU2 N1; MOI = EXTR L1 &BOU2 ; TABC . CAPI . MOI = 'PROG'; 'FIN' BOU2; 'FIN' BOU1; 'FINSI'; *-------------------------------------------- 'REPETER' BOU1 M1 ; PT1 = TABC.'RPOINTS' 'POIN' &BOU1; 'SI' (LOG1 'ET' (('DIME' (TABC.'NOMC')) 'EGA' M1 )) ; CAPI = 'EXTR' (TABC.'NOMC') &BOU1; 'SINON' ; CAPI = 'CHAINE' 'CAPT' &BOU1; 'FINSI'; 'REPETER' BOU2 N1; MOI = EXTR L1 &BOU2 ; *** captage de la pression 'SI' ('EGA' MOI 'PN'); X1 = 'EXTR' PN 'SCAL' PT1 ; * 'LISTE' X1 ; TABC . CAPI . MOI = (TABC . CAPI . MOI) 'ET' ('PROG' X1); 'FINSI'; *** captage de temperature 'SI' ('EGA' MOI 'TN'); X1 = 'EXTR' TN 'SCAL' PT1 ; TABC . CAPI . MOI = TABC . CAPI . MOI 'ET' ('PROG' X1); 'FINSI'; 'FIN' BOU2; 'FIN' BOU1; 'FINSI' ; *****---- creating listreel of time steps 'SI' (((&BL1 '/' FREQC) '*' FREQC) 'EGA' &BL1) ; TABC . 'LTPS' = (TABC . 'LTPS') 'ET' ('PROG' TPS); 'FINSI' ; *--------------------------------------------- 'SI' (((&BL1 '/' 20) '*' 20) 'EGA' &BL1) ; 'MESSAGE' ('CHAINE' 'ITER =' &BL1 ' TPS =' TPS) ; 'FINSI' ; *------------------------------------------------------- 'SI' (TPS '>EG' TFINAL) ; 'QUITTER' BL1 ; 'FINSI' ; 'FIN' BL1 ; *------------------------------------ *----------------------------------- 'FIN' BLG ; 'TEMPS' ; ************************************************* * Posttreatment * ************************************************* LISTP0 = 'PROG' 10 '*' 100000.0D0 ; OVPRES = (TABC . 'C1' . 'PN') '-' LISTP0 ; EVERP1 = 'EVOL' 'MANU' 'Niter' (TABC . 'LTPS' ) 'PC1' OVPRES ; *------------------------------------------------ LPGOOD = 'PROG' 1.3086 21.753 43.070 41.450 30.604 26.187 32.615 43.088 48.922 47.665 ; DIFFER = OVPRES '-' LPGOOD ; 'SI' (('MAXIMUM' DIFFER)