* fichier : ccar5.dgibi ************************************************************************ * Section : Fluides Transitoire ************************************************************************ ** ** --- 6 JUIN 1999 --- ** ** TEST CAVITE CARRE A PAROI DEFILANTE RE=400 ** ** Comparaison des algorithmes Implicite, projection ** et semi explicite. ** teste LAPL KONV NS (CENTREE) KBBT ** Formulation MACRO CENTRE (Stabilisation) ELEMENT QUA9 ** La comparaison fine se fait avec COMPLET=VRAI ** prevoir espace et temps cpu OPTION ISOV 'SULI' ; *OPTION TRACE PS; GRAPH = 'N' ; GRAPHI= 'N' ; COMPLET = FAUX ; *COMPLET = VRAI ; Rey = 400.; SI ( COMPLET ) ; ds1=0.02; ds2=0.1 ; ds1=0.01; ds2=0.08; ITMAX = 40 ; ITMAI = 30 ; ITMAE = 5000 ; KSUPG = 'CENTREE'; SINON ; err1=1.e-4 ; ds1=0.2; ds2=0.2 ; ITMAX = 30 ; ITMAI = 30 ; ITMAE = 2000; ITMAX = 10; ITMAI = 10; ITMAE = 100; KSUPG = 'SUPG'; FINSI ; KPRESS = 'CENTRE' ; DISCR = 'LINE' ; TYPK = 'QUA8' ; ALGO = 'SPLT' ; DEBPROC TEST KPRESS*MOT TYPK*MOT DISCR*MOT GRAPH*MOT ALGO*MOT ITMAX*ENTIER ITMAI*ENTIER ITMAE*ENTIER ; option dime 2 elem TYPK ; p1= 0 0 ; p12=0.5 0. ; p2= 1 0 ; ab=p1 d dini ds1 dfin ds2 p12 d dini ds2 dfin ds1 p2 ; mt= ab trans dini ds1 dfin ds2 (0 0.5) trans dini ds2 dfin ds1 (0 0.5) ; bc=cote 2 mt ; cd=cote 3 mt ; da=cote 4 mt ; AA=bc moin (0.5 0.) ; BB=ab plus (0. 0.5) ; elim (AA et BB et mt) 1.e-3 ; mt= chan mt QUAF ; $mt=mode mt 'NAVIER_STOKES' DISCR ; MU =1. ; RO= Rey ; EPSS=-1.e-8 ; DT=1. ; * La cavité est fermée il faut imposer la pression en 1 point ! prep1=doma $mt kpress; bcp=elem prep1 POI1 (lect 15); bcp = bcp coul rouge ; Mmt= doma $mt 'MAILLAGE' ; doma $mt 'IMPR' ; Si (EGA GRAPH 'O') ; trace (Mmt et bcp); Finsi ; * Mab = chan quaf ab ; Mcd = chan quaf cd ; Mda = chan quaf da ; Mbc = chan quaf bc ; elim (Mmt et Mab et Mcd et Mda et Mbc) 1.e-5 ; $ab = model Mab 'NAVIER_STOKES' DISCR ; $cd = model Mcd 'NAVIER_STOKES' DISCR ; $da = model Mda 'NAVIER_STOKES' DISCR ; $bc = model Mbc 'NAVIER_STOKES' DISCR ; ab = doma $ab maillage ; cd = doma $cd maillage ; da = doma $da maillage ; bc = doma $bc maillage ; Si (EGA ALGO 'EXPL'); RV= eqex 'OMEGA' 1. 'NITER' 1 ITMA itmax 'FIDT' 1 'OPTI' 'EFM1' KSUPG ZONE $mt OPER NS (MU/RO) INCO 'UN' 'OPTI' 'EFM1' 'CENTREE' ZONE $mt OPER DFDT 1. 'UN' 'DELTAT' INCO 'UN' ; Sinon ; RV= eqex 'OMEGA' 1. 'NITER' 1 ITMA itmax 'FIDT' 1 'OPTI' 'EF' 'IMPL' KSUPG * ZONE $mt OPER NS (MU/RO) INCO 'UN' ZONE $mt 'OPER' 'LAPN' MU 'INCO' 'UN' ZONE $mt 'OPER' 'KONV' RO 'UN' MU DT 'INCO' 'UN' 'OPTI' 'EFM1' 'CENTREE' ZONE $mt OPER DFDT RO 'UN' DT INCO 'UN' ; Finsi ; RV= eqex RV CLIM UN UIMP CD 1. UN VIMP CD 0. UN UIMP DA 0. UN VIMP DA 0. UN UIMP AB 0. UN VIMP AB 0. UN UIMP BC 0. UN VIMP BC 0. ; rv.'METHINV'.TYPINV=3 ; rv.'METHINV'.IMPINV=0 ; rv.'METHINV'.NITMAX=400; rv.'METHINV'.PRECOND=3 ; rv.'METHINV'.RESID =1.e-8 ; rv. 'METHINV' . 'FCPRECT'=1 ; rv. 'METHINV' . 'FCPRECI'=1 ; 'SI' (EGA ALGO 'IMPL') ; RV.'ITMA' = ITMAI ; RV.'OMEGA' = 1. ; RV= eqex RV 'OPTI' 'EF' 'IMPL' KSUPG KPRESS ZONE $mt OPER KBBT 1. INCO 'UN' 'PRES' CLIM PRES TIMP bcp 0. ; 'SINON' ; Si (EGA ALGO 'EXPL') ; betastab = 1. ; Sinon ; betastab= 1.e2 ; Finsi ; Si (EGA TYPK 'TRI6') ; betastab = 1.e5 ; Finsi ; RVP= EQEX 'OPTI' 'EFM1' KPRESS ZONE $mt OPER KBBT -1. betastab INCO 'UN' 'PRES' CLIM PRES TIMP bcp 0. ; rvp.'METHINV'.TYPINV=2 ; rvp.'METHINV'.TYPINV=1 ; rvp.'METHINV'.IMPINV=0 ; rvp.'METHINV'.NITMAX=300; rvp.'METHINV'.PRECOND=3 ; rvp.'METHINV'.RESID =1.e-8 ; rvp.'METHINV' . 'FCPRECT'=100 ; rvp.'METHINV' . 'FCPRECI'=100 ; si (EGA ALGO 'SPLT') ; RV.'PROJ'= RVP ; RV.'ITMA' = ITMAX ; Finsi ; si (EGA ALGO 'EXPL') ; RV.'PROJ' = RVP ; RV.'TYPROJ'='PSCT'; RV.'ITMA' = ITMAE ; RV.'ALFA' = 0.5 ; RV.'FIDT' = 50 ; Finsi ; 'FINSI' ; rv.inco= table inco ; rv.inco.un = kcht $mt vect sommet (1.e-5 1.e-5) ; rv.inco.'PRES'= kcht $mt scal kpress 0.; exec rv ; exec rv ; AA = chan AA quaf ; BB = chan BB quaf ; elim (mt et aa et bb ) 1.e-3 ; $AA=mode AA 'NAVIER_STOKES' DISCR ; $BB=mode BB 'NAVIER_STOKES' DISCR ; srti=doma $AA 'MAILLAGE' ; srth=doma $BB 'MAILLAGE' ; evolV = EVOL 'CHPO' (rv.'INCO'.'UN') UX (srti ) ; evolH = EVOL 'CHPO' (rv.'INCO'.'UN') UY (srth ) ; evx='EXTR' evolV 'ORDO' ; list evx ; evy='EXTR' evolV 'ABSC' ; evolV= evol 'MANU' 'Vitesse' evx 'Hauteur' evy ; rv.'EVOLV'=evolV ; rv.'EVOLH'=evolH ; si ('EGA' graph 'O' ); TAB1=TABLE; TAB1.'TITRE'=TABLE ; TAB1 . 1 ='MARQ REGU ' ; TAB1.'TITRE' . 1 = mot 'Composante_UX ' ; DESS evolV 'TITX' 'R (m)' 'TITY' 'V (m/s)' LEGE TAB1 ; TAB1.'TITRE' . 1 = mot 'Composante_UY ' ; DESS evolH 'TITX' 'R (m)' 'TITY' 'V (m/s)' LEGE TAB1 ; c1=vect (rv.inco.'UN') 0.3 ux uy jaune ; trace c1 mt ; Si (EGA ALGO 'IMPL') ; pn= elno $mt (kcht $mt scal centre ( nomc scal (rv.inco.'PRES'))); sinon ; pn= elno $mt (kcht $mt scal centre ( nomc scal (rv.inco.pression))); finsi ; trace pn mt ; mess ' MAX P ' (maxi pn) 'MIN P ' (mini pn) ; finsi ; FINPROC RV ; RV1= test 'CENTRE' TYPK 'MACRO' GRAPH 'SPLT' ITMAX ITMAI ITMAE ; evv= (rv1.'EVOLV'); evh= (rv1.'EVOLH'); SI ( NON COMPLET ) ; lr = 'EXTRAIRE' (rv1.'EVOLV') 'ABSC' ; list lr ; lrr = prog -8.26684E-38 -2.74863E-02 -5.08159E-02 -7.49178E-02 -9.42076E-02 -9.99882E-02 -8.10852E-02 -4.32674E-02 7.12863E-03 5.44891E-02 0.10073 0.30942 1.0000; ER=SOMM( abs (lr - lrr) ) ; mess ' Ecart sur CENTRE QUA8 MACRO Projection ' er ; Si ( er > err1 ) ; erreur 5 ; finsi ; FINSI ; RV2= test 'CENTRE' TYPK 'MACRO' GRAPH 'IMPL' ITMAX ITMAI ITMAE ; evv= evv et (rv2.'EVOLV'); evh= evh et (rv2.'EVOLH'); SI ( NON COMPLET ) ; lr = 'EXTRAIRE' (rv2.'EVOLV') 'ABSC' ; list lr ; lrr=prog 0.0000 -3.10947E-02 -5.38085E-02 -7.39196E-02 -8.54276E-02 -8.30593E-02 -6.45006E-02 -3.23991E-02 3.13243E-03 3.90764E-02 8.72317E-02 0.29477 1.0000; ER=SOMM( abs (lr - lrr) ) ; mess ' Ecart sur CENTRE QUA8 MACRO Implicite ' er ; Si ( er > err1 ) ; erreur 5 ; finsi ; FINSI ; RV3= test 'CENTRE' TYPK 'MACRO' GRAPH 'EXPL' ITMAX ITMAI ITMAE ; list RV3.'TYPROJ' ; evv= evv et (rv3.'EVOLV'); evh= evh et (rv3.'EVOLH'); SI ( NON COMPLET ) ; lr = 'EXTRAIRE' (rv3.'EVOLV') 'ABSC' ; list lr ; lrr=prog 3.29048E-37 -1.67015E-02 -2.71125E-02 -3.64044E-02 -4.66237E-02 -5.61821E-02 -5.96575E-02 -5.28560E-02 -3.42689E-02 -2.85057E-03 5.39080E-02 0.27875 1.0000; ER=SOMM( abs (lr - lrr) ) ; mess ' Ecart sur CENTRE TRI6 MACRO Projection ' er ; Si ( er > err1 ) ; erreur 5 ; finsi ; FINSI ; *FIN ; si ('EGA' graphi 'O' ); * Solution de Ghia et al (1982) : y = 'PROG' 0. 0.0547 0.0625 0.0703 0.1016 0.1719 0.2813 0.4531 0.5 0.6172 0.7344 0.8516 0.9531 0.9609 0.9688 0.9766 1.; 'SI' (Rey 'EGA' 100.) ; u = 'PROG' 0. -0.03717 -0.04192 -0.04775 -0.06434 -0.10150 -0.15662 -0.2109 -0.20581 -0.13641 0.00332 0.23151 0.68717 0.73722 0.78871 0.84123 1.; 'FINSI' ; 'SI' (Rey 'EGA' 400.) ; u = 'PROG' 0. -0.08186 -0.09266 -0.10338 -0.14612 -0.24299 -0.32726 -0.17119 -0.11477 0.02135 0.16256 0.29093 0.55892 0.61756 0.68439 0.75837 1.; 'FINSI' ; 'SI' (Rey 'EGA' 1000.) ; u = 'PROG' 0. -0.18109 -0.20196 -0.22220 -0.29730 -0.38289 -0.27805 -0.10648 -0.06080 0.05702 0.18719 0.33304 0.46604 0.51117 0.57492 0.65928 1.; 'FINSI' ; evolug = 'EVOL' 'MANU' u y; TAB1=TABLE; TAB1.'TITRE'=TABLE ; TAB1 . 1 ='MARQ CROI REGU ' ; TAB1 . 2 ='MARQ TRIA REGU ' ; TAB1 . 3 ='MARQ LOSA REGU ' ; TAB1 . 4 ='TIRR MARQ PLUS REGU'; TAB1.'TITRE' . 1 = mot ' Proj ' ; TAB1.'TITRE' . 2 = mot ' Impl ' ; TAB1.'TITRE' . 3 = mot ' 1/2 Expl' ; TAB1.'TITRE' . 4 = MOT 'V Ghia et al (1982)'; evv = evv et evolug ; DESS evv 'TITX' 'R (m)' 'TITY' 'V (m/s)' LEGE TAB1 DATE MIMA TITR ' Cavite carree (Qua9 Macro stab comp ux) ; Comparaison Algorithmes '; * Solution de Ghia et al (1982) : x = 'PROG' 0. 0.0625 0.0703 0.0781 0.0938 0.1563 0.2266 0.2344 0.5 0.8047 0.8594 0.9063 0.9453 0.9531 0.9609 0.9688 1.; 'SI' (Rey 'EGA' 100.); v = 'PROG' 0. 0.09233 0.10091 0.10890 0.12317 0.16077 0.17507 0.17527 0.05454 -0.24533 -0.22445 -0.16914 -0.10313 -0.08864 -0.07391 -0.05906 0.; 'FINSI' ; 'SI' (Rey 'EGA' 400.); v = 'PROG' 0. 0.18360 0.19713 0.20920 0.22965 0.28124 0.30203 0.30174 0.05186 -0.38598 -0.44993 -0.33827 -0.22847 -0.19254 -0.15663 -0.12146 0.; 'FINSI' ; 'SI' (Rey 'EGA' 1000.); v = 'PROG' 0. 0.27485 0.29012 0.30353 0.32627 0.37095 0.33075 0.32235 0.02526 -0.31966 -0.42665 -0.51550 -0.39188 -0.33714 -0.27669 -0.21388 0.; 'FINSI' ; evolvg = 'EVOL' 'MANU' x v; evh = evh et evolvg ; DESS evh 'TITX' 'R (m)' 'TITY' 'V (m/s)' LEGE TAB1 DATE MIMA TITR ' Cavite carree (Qua9 Macro stab comp uy) ; Comparaison Algorithmes '; finsi ; *Si complet ; *opti sauv 'ccar5.dgibi.sauv' ; *sauv ; *finsi ; FIN ;