Télécharger konv_fmm_test.dgibi
* fichier : konv_fmm_test.dgibi ************************************************************************ ************************************************************************ *********************************************************************** * * * 'KONV' OPERATOR * * FREE MATRIX METHOD implicitation. * * * * VF "cell-centered" discretization of the Euler equations. * * Unknowns : U (density, momentum, total energy per volume unity * * (conservative variables)) * * In the i-th cell we have to compute * * * * (U_i^{n+1} - U_i^{n}) * AN_i(U^{n}) = * * RES_i(U^{n}) + BN_i(U^{n}) - BN_i(U_i^{n+1}) * * * * In this test case we compute AN et BN in GIBIANE et en ESOPE. * * * * A. Beccantini * * SFME/LTMF * * 03/05/01 * *********************************************************************** * TYEL = 'QUA4' ; * TYEL = 'TRI3' ; GRAPH = VRAI ; GRAPH = FAUX ; ****************** **** MAILLAGE **** ****************** RAF = 1 ; NX = 4 ; NY = 4 ; DX = 1. ; A0 = 0.11 0.0 ; A1 = 3.11 0.0 ; A2 = 3.13 4.1 ; A3 = 0.0 4.12 ; * **** LIGB * LIGB= A0 'DROIT' NX A1 ; * **** LIGH * LIGH = A3 'DROIT' NX A2 ; * **** DOMINT * DOMINT = LIGH 'REGLER' NY LIGB ; LIGCON = 'CONTOUR' DOMINT ; * *** LIGG * * **** LIGD * FROD = LIGD 'TRANSLATION' 1 (DX 0.0) ; FROG = LIGG 'TRANSLATION' 1 ((-1.0 '*' DX) 0.0) ; DOMTOT = DOMINT 'ET' FROG 'ET' FROD ; 'ELIMINATION' DOMTOT (1.0D-3 '/' RAF) ; * **** Creation of DOMAINE tables via the MODEL object * $DOMTOT = 'MODELISER' DOMTOT 'EULER' ; $DOMINT = 'MODELISER' DOMINT 'EULER' ; $FROD = 'MODELISER' FROD 'EULER' ; $FROG = 'MODELISER' FROG 'EULER' ; $LIGG = 'MODELISER' LIGG 'EULER' ; $LIGD = 'MODELISER' LIGD 'EULER' ; QDOMTOT = TDOMTOT . 'QUAF' ; QDOMINT = TDOMINT . 'QUAF' ; QFROD = TFROD . 'QUAF' ; QFROG = TFROG . 'QUAF' ; QLIGG = TLIGG . 'QUAF' ; QLIGD = TLIGD . 'QUAF' ; 'ELIMINATION' QDOMTOT (1.0D-3 '/' RAF) QDOMINT ; 'ELIMINATION' QDOMTOT (1.0D-3 '/' RAF) QFROD ; 'ELIMINATION' QDOMTOT (1.0D-3 '/' RAF) QFROG ; 'ELIMINATION' QDOMTOT (1.0D-3 '/' RAF) QLIGD ; 'ELIMINATION' QDOMTOT (1.0D-3 '/' RAF) QLIGG ; * **** SEG2 mesh for BC * ELP1 ; ELP1 ; 'SI' (NBL1 > 0) ; 'REPETER' BL1 NBL1 ; ELP1 ; ELP1 ; 'FIN' BL1 ; 'FINSI' ; ELP1 ; ELP1 ; 'SI' (NBL1 > 0) ; 'REPETER' BL1 NBL1 ; ELP1 ; ELP1 ; 'FIN' BL1 ; 'FINSI' ; 'SI' GRAPH ; 'Domaine total' ; 'FINSI' ; ***************************************************** ***************************************************** ******* Initial conditions ************************** ***************************************************** ***************************************************** *************************************************************** ***** PROCEDURE POUR CALCULER LES VARIABLES CONSERVATIVES ***** *************************************************************** 'DEBPROC' CONS ; 'ARGUMENT' RN*'CHPOINT' VN*'CHPOINT' PN*'CHPOINT' GAMN*'CHPOINT' ; RECIN = 0.5 '*' CELL ; REIN = PN '/' (GAMN '-' 1.0) ; RETN = RECIN '+' REIN ; DETR CELL ; DETR RECIN ; DETR REIN ; 'RESPRO' RVN RETN ; 'FINPROC' ; *************************************************************** *************************************************************** *************************************************************** 'NATURE' 'DISCRET') '*' 'UY' 0.0 'NATURE' 'DISCRET') '*' 'NATURE' 'DISCRET') '*' 'NATURE' 'DISCRET' ; GN0 RETN0 = CONS RN0 VN0 PN0 GAMN ; ERRO = 'MAXIMUM' (PN1 '-' PN0) 'ABS' ; 'SI' (ERRO > 1.0D-6) ; 'MESSAGE' 'Problem in the ic file!!!' 'ERREUR' 5 ; 'FINSI' ; * **** Plot of IC * 'SI' GRAPH ; 'TRACER' CHM_RN MOD1 'TRACER' CHM_PN MOD1 'TRACER' CHM_VN MOD1 'FINSI' ; MOTRN = 'RN' ; MOTRNX = 'RUX' ; MOTRNY = 'RUY' ; MOTRETN = 'RETN' ; ********************************************************************* ********************** TEST ***************************************** ********************************************************************* * **** Parameter for the computations * * * Safety Factor fot the time step SAFFAC = 0.5 ; * ********** Procedure qui calcule F \cdot n = * * COEFF1 = \rho u_n * COEFF2 = \rho u_n ux '+' p nx * COEFF2 = \rho u_n uy '+' p ny * COEFF4 = \rho u_n h_t * 'DEBPROC' FSCALN; 'ARGUMENT' UX'*'FLOTTANT UY'*'FLOTTANT NX'*'FLOTTANT NY'*'FLOTTANT RHO'*'FLOTTANT P'*'FLOTTANT GAMMA'*'FLOTTANT; UN = (UX '*' NX) '+' (UY '*' NY) ; RECIN = 0.5 '*' RHO '*' ((UX '*' UX) '+' (UY '*' UY)) ; RH = (GAMMA '/' (GAMMA '-' 1.)) '*' P ; COEF1 = RHO '*' UN ; COEF2 = (RHO '*' UN '*' UX) '+' (P '*' NX) ; COEF3 = (RHO '*' UN '*' UY) '+' (P '*' NY) ; COEF4 = UN '*' (RH '+' RECIN) ; 'FINPROC' COEF1 COEF2 COEF3 COEF4 ; * **** Beginning of the main program * RN = 'COPIER' RN0 ; GN = 'COPIER' GN0 ; RETN = 'COPIER' RETN0 ; * We define DIAMIN = CHAMPOINT 'CENTRE' (diameter minimum of the cell). * One component, 'SCAL'. * DOMVOL = CHAMPOINT 'CENTRE' (volume of the cell). * One component, 'SCAL'. * SURDOM = CHAMPOINT 'FACE' (surface of each interface). * One component, 'SCAL'. * NORDOM = CHAMPOINT 'FACE' (normales to each interface) * 2 components, 'UX', 'UY' * * Names of the residuum components * * **** Primitive variables * * **** Local time step UNSDT = 1 '/' DT * * **** Loop to solve the FMM linear system with JACOBI * * a_i (DELTA U)_i^n = RES_{expl,i} '+' b_i^n '-' b_i^{n+1} * * a_i is a CHAMPOINT centre with one componenent ('SCAL) * * b_i is a CHAMPOINT centre with four components (see LISTINCO) * * b_i^{n+1} a CHAMPOINT centre with four components (see LISTINCO) * * First of all we compute a_i et b_i^n * * We also compute SIGMAF = CHAMPOINT face with one component. It * contains (un '+' c) at each face. * MOTRNX 0.0 MOTRNY 0.0 MOTRETN 0.0 ; * **** LOOP ON 'FACEL' * 'REPETER' BLFAC NPOIFA ; * * Each elements has tree points * | * +L +Cf +R * | * L,R are cell centers * Cf is a the face center * Normal are always directed frol L to Cf * ELT = 'CHANGER' 'POI1' ELTINI ; * ********** Murs * ELTL = 'MANUEL' 'POI1' PL ; ELTF = 'MANUEL' 'POI1' PF ; * Geometric parameters VOLL = 'EXTRAIRE' DOMVOL PL 'SCAL' ; DIAL = 'EXTRAIRE' DIAMIN PL 'SCAL' ; SURC = 'EXTRAIRE' SURDOM PF 'SCAL' ; NORX = 'EXTRAIRE' NORDOM PF 'UX' ; NORY = 'EXTRAIRE' NORDOM PF 'UY' ; * Physical values UNXL = 'EXTRAIRE' VN PL 'UX' ; UNYL = 'EXTRAIRE' VN PL 'UY' ; RNL = 'EXTRAIRE' RN PL 'SCAL' ; PNL = 'EXTRAIRE' PN PL 'SCAL' ; GAML = 'EXTRAIRE' GAMN PL 'SCAL' ; UNNL = (UNXL '*' NORX) '+' (UNYL '*' NORY) ; UNTL = (UNYL '*' NORX) '-' (UNXL '*' NORY) ; * RNR = RNL ; PNR = PNL ; GAMR = GAML ; UNNR = -1.0 '*' UNNL ; UNTR = UNTL ; UNXR = (UNNR '*' NORX) '-' (UNTR '*' NORY) ; UNYR = (UNNR '*' NORY) '+' (UNTR '*' NORX) ; * Computations RETL = ((1. '/' (GAML '-' 1.)) '*' PNL) '+' (0.5 '*' ((UNXL '*' UNXL) '+' (UNYL '*' UNYL)) '*' RNL) ; RETR = ((1. '/' (GAMR '-' 1.)) '*' PNR) '+' (0.5 '*' ((UNXR '*' UNXR) '+' (UNYR '*' UNYR)) '*' RNR) ; CSONF = ((0.5 '*' (GAML '+' GAMR)) '*' (PNL '+' PNR) '/' (RNL '+' RNR)) ; CSONF = CSONF '**' 0.5 ; UNF = 0.5 '*' (('ABS' UNNL) '+' ('ABS' UNNR)) ; SIGF = CSONF '+' UNF ; * SIGMAF SIGMAF = SIGMAF '+' SIGNF ; * UNSDT (SIGF '/' DIAL) ; UNSDT = 0.5 '*' ( (UNSDT '+' UNSDTL) '+' ('ABS' (UNSDT '-' UNSDTL)) ) ; * ANI (2. '*' VOLL)) ; ANI = ANI '+' ANL ; * BNI * ********** Procedure qui calcule F \cdot n = * * COEFF1 = \rho u_n * COEFF2 = \rho u_n ux '+' p nx * COEFF2 = \rho u_n uy '+' p ny * COEFF4 = \rho u_n h_t * * COEF1L COEF2L COEF3L COEF4L = FSCALN UNXL UNYL NORX NORY RNL PNL GAML ; COEF1R COEF2R COEF3R COEF4R = FSCALN UNXR UNYR NORX NORY RNR PNR GAMR ; SSDVL = SURC '/' (2. '*' VOLL) ; MOTRN ((COEF1L '+' COEF1R '-' (SIGF '*' RNR)) '*' SSDVL) MOTRNX ((COEF2L '+' COEF2R '-' (SIGF '*' RNR '*' UNXR)) '*' SSDVL) MOTRNY ((COEF3L '+' COEF3R '-' (SIGF '*' RNR '*' UNYR)) '*' SSDVL) MOTRETN ((COEF4L '+' COEF4R '-' (SIGF '*' RETR)) '*' SSDVL) ; BNI = BNI '+' BNL ; 'SINON' ; * ********** Domaine interieur * ELTL = 'MANUEL' 'SEG2' PL PR ; ELTL1 = 'CHANGER' 'POI1' ELTL ; ELTF = 'MANUEL' 'POI1' PF ; * Geometric parameters VOLL = 'EXTRAIRE' DOMVOL PL 'SCAL' ; VOLR = 'EXTRAIRE' DOMVOL PR 'SCAL' ; DIAL = 'EXTRAIRE' DIAMIN PL 'SCAL' ; DIAR = 'EXTRAIRE' DIAMIN PR 'SCAL' ; SURC = 'EXTRAIRE' SURDOM PF 'SCAL' ; NORX = 'EXTRAIRE' NORDOM PF 'UX' ; NORY = 'EXTRAIRE' NORDOM PF 'UY' ; * Physical values UNXL = 'EXTRAIRE' VN PL 'UX' ; UNYL = 'EXTRAIRE' VN PL 'UY' ; RNL = 'EXTRAIRE' RN PL 'SCAL' ; PNL = 'EXTRAIRE' PN PL 'SCAL' ; GAML = 'EXTRAIRE' GAMN PL 'SCAL' ; * UNXR = 'EXTRAIRE' VN PR 'UX' ; UNYR = 'EXTRAIRE' VN PR 'UY' ; RNR = 'EXTRAIRE' RN PR 'SCAL' ; PNR = 'EXTRAIRE' PN PR 'SCAL' ; GAMR = 'EXTRAIRE' GAMN PR 'SCAL' ; * Computations UNNL = (UNXL '*' NORX) '+' (UNYL '*' NORY) ; UNNR = (UNXR '*' NORX) '+' (UNYR '*' NORY) ; RETL = ((1. '/' (GAML '-' 1.)) '*' PNL) '+' (0.5 '*' ((UNXL '*' UNXL) '+' (UNYL '*' UNYL)) '*' RNL) ; RETR = ((1. '/' (GAMR '-' 1.)) '*' PNR) '+' (0.5 '*' ((UNXR '*' UNXR) '+' (UNYR '*' UNYR)) '*' RNR) ; CSONF = ((0.5 '*' (GAML '+' GAMR)) '*' (PNL '+' PNR) '/' (RNL '+' RNR)) ; CSONF = CSONF '**' 0.5 ; UNF = 0.5 '*' (('ABS' UNNL) '+' ('ABS' UNNR)) ; SIGF = CSONF '+' ('ABS' UNF) ; * SIGMAF SIGMAF = SIGMAF '+' SIGNF ; * UNSDT UNSDT = 0.5 '*' ( (UNSDT '+' UNSDTL) '+' ('ABS' (UNSDT '-' UNSDTL)) ) ; * ANI ((SURC '*' SIGF) '/' (2. '*' VOLR))) ; ANI = ANI '+' ANL ; * BNI COEF1L COEF2L COEF3L COEF4L = FSCALN UNXL UNYL NORX NORY RNL PNL GAML ; COEF1R COEF2R COEF3R COEF4R = FSCALN UNXR UNYR NORX NORY RNR PNR GAMR ; SSDVL = SURC '/' (2. '*' VOLL) ; SSDVR = SURC '/' (-2. '*' VOLR) ; MOTRN ('PROG' ((COEF1L '+' COEF1R '-' (SIGF '*' RNR)) '*' SSDVL) ((COEF1L '+' COEF1R '+' (SIGF '*' RNL)) '*' SSDVR)) MOTRNX ('PROG' ((COEF2L '+' COEF2R '-' (SIGF '*' RNR '*' UNXR)) '*' SSDVL) ((COEF2L '+' COEF2R '+' (SIGF '*' RNL '*' UNXL)) '*' SSDVR)) ; MOTRNY ('PROG' ((COEF3L '+' COEF3R '-' (SIGF '*' RNR '*' UNYR)) '*' SSDVL) ((COEF3L '+' COEF3R '+' (SIGF '*' RNL '*' UNYL)) '*' SSDVR)) MOTRETN ('PROG' ((COEF4L '+' COEF4R '-' (SIGF '*' RETR)) '*' SSDVL) ((COEF4L '+' COEF4R '+' (SIGF '*' RETL)) '*' SSDVR)) ; BNI = BNI '+' BNL1 '+' BNL2 ; 'FINSI' ; 'FIN' BLFAC ; UNSDT = UNSDT '/' (0.5 * SAFFAC) ; ANI = ANI '+' UNSDT ; ANI1 = 'KONV' 'VF' 'PMONOFMM' 'AN' LISTINCO $DOMTOT RN GN RETN GAMN SAFFAC ; BNI1 = 'KONV' 'VF' 'PMONOFMM' 'BN' LISTINCO $DOMTOT RN GN RETN GAMN ; 'SI' GRAPH ; 'TRAC' CHM_ER1 MOD1 ; 'TRAC' CHM_ER2 MOD1 ; 'FINSI' ; ERR1 = 'MAXIMUM' (BNI '-' BNI1) 'ABS' ; 'MESSAGE' 'Probleme en KONV' ; 'ERREUR' 5 ; 'FINSI' ; ERR2 = 'MAXIMUM' (BNI '-' BNI1) 'ABS' ; 'MESSAGE' 'Probleme en KONV' ; 'ERREUR' 5 ; 'FINSI' ; 'FIN' ;
© Cast3M 2003 - Tous droits réservés.
Mentions légales