Numérotation des lignes :

* fichier :  dvispassi.dgibi************************************************************************************************************************************************* CAVITE CARREE ­ VAHL­DAVIS Methode de projection implicite ** Malvina Renesson Aout 1999 * *************************************************************** FIN ;GRAPH=VRAI ; GRAPH=FAUX ;COMPLET=VRAI; COMPLET=FAUX ;TCPT   = FAUX;TKPR   = VRAI;TRESOU = FAUX;IMPARA = FAUX; Si (NON Complet);itmax= 80 ; DISCR='LINE' ;*DISCR='MACRO'; DISCR='QUAF';d1 = 0.04; d2 = 0.12; d1 = 0.02; d2 = 0.06;DT=.25;Ra = 1.e6 ;Sinon ;*Ra = 1.e5 ;*itmax= 600 ;Ra = 1.e6 ;itmax= 600 ;DISCR='LINE' ;d1 = 0.01; d2 = 0.07;DT=.25; Finsi ; KSUPG ='CENTREE' ; KPRESS='MSOMMET' ; KPRESS='CENTRE'; KPRESS='CENTREP1'; ****************************** *ESTIMATION DE LA CONVERGENCE* ****************************** DEBPROC CALCUL; ARGU RVX*'TABLE' MDZ/'MMODEL' MDZ/'TABLE';RV = RVX.'EQEX'; DD = RV.PASDETPS.'NUPASDT'; NN = DD/5; LO = (DD-(5*NN)) EGA 0;SI (LO); UN = RV.INCO.'UN'; UNM1 = RV.INCO.'UNM1'; unx = kcht $mt scal sommet (exco 'UX' un);unm1x = kcht$mt scal sommet (exco 'UX' unm1);uny = kcht $mt scal sommet (exco 'UY' un);unm1y = kcht$mt scal sommet (exco 'UY' unm1);ERRX = KOPS unx '-' unm1x;ERRY = KOPS uny '-' unm1y;ELIX = MAXI ERRX 'ABS'; ELIY = MAXI ERRY 'ABS'; ELIX = (LOG (ELIX + 1.0E-10))/(LOG 10.);ELIY = (LOG (ELIY + 1.0E-10))/(LOG 10.);MESSAGE 'ITER' RV.PASDETPS.'NUPASDT' 'ERREUR LINF' ELIX ELIY; RV.INCO.'UNM1' = KCHT $MT 'VECT' 'SOMMET' (RV.INCO.'UN'); IT = PROG RV.PASDETPS.'NUPASDT'; ER = PROG ELIY; RV.INCO.'IT' = (RV.INCO.'IT') ET IT; RV.INCO.'ER' = (RV.INCO.'ER') ET ER; FINSI; as2 ama1 = 'KOPS' 'MATRIK' ; 'FINPROC' as2 ama1 ; ************ * MAILLAGE * ************ OPTI DIME 2; OPTI ELEM QUA8;crit = 0.0000001; xxll = 1.;p1 = 0. 0.; p15 = (xxll/2.) 0.; p2 = xxll 0.; p25 = xxll 0.5; p3 = xxll 1.; p35 = (xxll/2.) 1.; p4 = 0. 1.; p45 = 0. 0.5; bas1 = p1 d dini d1 dfin d2 p15; bas2 = p15 d dini d2 dfin d1 p2;bas = (bas1 et bas2);elim crit bas; cdro = p2 d dini d1 dfin d2 p25 d dini d2 dfin d1 p3; haut = p3 d dini d1 dfin d2 p35 d dini d2 dfin d1 p4; cgau = p4 d dini d1 dfin d2 p45 d dini d2 dfin d1 p1; coupx = bas plus (0. 0.5) ;coupy=cgau plus (0.5 0.) ;cnt = bas et cdro et haut et cgau; mt = bas cdro haut cgau daller; bord = cont mt; *mt = surf bord; elim (mt et coupx) 1.e-5;Si GRAPH ; trace mt ; Finsi ;*option donn 5;orienter mt; ******** * DOMA * ******** Mmt = chan QUAF mt ;$mt = MODE Mmt 'NAVIER_STOKES' DISCR ;mt = doma $mt maillage; Mcgau=chan cgau QUAF ;$cgau=mode Mcgau 'NAVIER_STOKES' DISCR ; Mcdro=chan cdro QUAF ; $cdro=mode Mcdro 'NAVIER_STOKES' DISCR ; Mcnt=chan cnt QUAF ;$cnt=mode Mcnt 'NAVIER_STOKES' DISCR ; Mcoupx = chan coupx 'QUAF' ;$coupx=mode Mcoupx 'NAVIER_STOKES' DISCR ; Mcoupy = chan coupy 'QUAF' ;$coupy=mode Mcoupy 'NAVIER_STOKES' DISCR ; Elim (Mmt et Mcgau et Mcdro et Mcnt et Mcoupxet Mcoupy) 1.e-5 ;cgau=doma $cgau maillage ;cdro=doma$cdro maillage ;cnt=doma $cnt maillage ;coupx=doma$coupx maillage ;coupy=doma $coupy maillage ; *option donn 5; ************** * PARAMETRES * ************** Pr = 0.71; Gr = Ra/Pr; NU = 1/(Gr**0.5); ALF= NU/Pr; gb = kcht$mt vect centre (0. -1.);uref = 0.2 ; *********************** * CREATION DES TABLES ************************ RV = EQEX  'OMEGA' 1.    'NITER' 1  ITMA itmax   'FIDT' 1  'OPTI' 'EF' 'IMPL' KSUPG'ZONE' $MT 'OPER' CALCUL'ZONE'$MT 'OPER' 'NS' 1. 'UN' NU 'GB' 'TN' 0.5 'INCO' 'UN''ZONE' $MT 'OPER' 'TSCAL' 1. 'UN' ALF 0. 'INCO' 'TN''OPTI' 'CENTREE' ZONE$mt    OPER  DFDT 1. 'UN' DT  INCO 'UN' ZONE  $mt OPER DFDT 1. 'TN' DT INCO 'TN'; **** RV = EQEX RV 'CLIM' 'UN' 'UIMP' cnt 0. 'UN' 'VIMP' cnt 0. 'TN' 'TIMP' cgau 1.'TN' 'TIMP' cdro 0.; rv.'METHINV'.TYPINV=3 ; rv.'METHINV'.IMPINV=0 ; rv.'METHINV'.NITMAX=400; rv.'METHINV'.PRECOND=3 ; rv.'METHINV'.RESID =1.e-10; rv. 'METHINV' . 'FCPRECT'=1 ; rv. 'METHINV' . 'FCPRECI'=1 ; ***** La cavité est fermée il faut imposer la pression en 1 point !prep1=doma$mt kpress;bcp=elem  prep1 POI1 (lect 1 ) ; RVP = EQEX'OPTI' 'EF' KPRESSZONE  $mt OPER KBBT -1. INCO 'UN' 'PRES'CLIM PRES TIMP bcp 0. ; 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 ; RV.INCO = TABLE INCO;RV.INCO.'UN' = kcht$MT VECT SOMMET (0. 0.); *RV.INCO.'TN' = kcht $MT SCAL SOMMET (1.-corsx);RV.INCO.'TN' = kcht$MT SCAL SOMMET 0.5 ;RV.INCO.'PRES' = kcht $MT SCAL KPRESS 0. ;RV.INCO.'UNM1' = kcht$MT VECT SOMMET (1.E-5 1.E-5);RV.INCO.'IT' = PROG 1; RV.INCO.'ER' = PROG 0.; RV.INCO.'GB' =  gb ;RV.'PROJ' = RVP;  RV.'TCPT'  = TCPT  ; RV.'TKPR'  = TKPR  ; RV.'TRESOU'= TRESOU; RV.'IMPARA'= IMPARA;temps;EXEC RV; temps ; ************* * RESULTATS * ************* *option sauv 'SAV';*sauv Pr Ra Gr NU ALF gb uref RV;*************OPTI ISOV SULI;   un=(RV.INCO.'UN') * (Pr*(Gr**0.5)) ;unch = vect un (0.3/(Pr*(Gr**0.5))) UX UY JAUNE;  u1y = exco un 'UX' ;u2x = exco un 'UY' ;evu=evol 'CHPO'  u2x  coupx ;evv=evol 'CHPO'  u1y  coupy ; EVOL4 = EVOL 'MANU' 'ITERATIONS' (RV.INCO.'IT') 'LOG---E---inf' (RV.INCO.'ER'); GRADT = KOPS RV.INCO.'TN' 'GRAD' $MT; DTDX = KCHT$MT 'SCAL' 'CENTRE' (EXCO 'UX' GRADT); DTDX = (-1)*(ELNO $MT DTDX); EVOLd = EVOL 'CHPO' DTDX SCAL cdro ;EVOLg = EVOL 'CHPO' DTDX SCAL cgau ;evol1=evold et evolg ;LISTE EVOL1; TAB1 = TABLE; TAB1.1 = 'MARQ LOSA';TAB1.2 = 'MARQ CROI'; Si GRAPH ; DESS EVOL1 'XBOR' 0. 1. 'GRIL' TAB1 'TITRE' ' Nusselt' ; trac RV.INCO.'TN' mt (cont mt) 14 'TITRE' ' Isothermes' ; trace unch mt (cont mt) 'TITRE' ' Vitesses' ; dess evu 'TITRE' 'Coupe ox à y=1/2 (uy)' ;dess evv 'TITRE' 'Coupe oy à x=1/2 (ux)' ;dess EVOL4 'XBOR' 0. itmax 'YBOR' -10.0 0.0;Finsi; Nusd = (redu DTDX cdro);Nusg = (redu DTDX cgau); Numg = somt (Nusg*(doma$cgau 'XXDIAGSI')) ;Numd = somt (Nusd*(doma $cdro 'XXDIAGSI')) ; un2=KOPS un '*' (Pr*(Gr**0.5)); sw = KOPS un 'ROT'$mt ; rk = EQEX $mt 'OPTI' 'EF' 'IMPL' ZONE$mt OPER LAPN -1. INCO 'PSI' ZONE $mt OPER FIMP sw INCO 'PSI' 'CLIM' 'PSI' 'TIMP' cnt 0.; rk.'INCO'=table 'INCO' ; rk.'INCO'.'PSI'=kcht$mt scal sommet 0. ; EXec rk ;psi=rk.'INCO'.'PSI'; psi2=kops psi '*' (Pr*(Gr**0.5));  Si GRAPH ; trac psi mt CNT 14 'TITRE' 'Fonction de courant '; Finsi; opti 'ECHO' 0 ;MESS '----------- Ra = ' Ra '------------------------------------';MESS '                     NUSSELT MOYEN   NUSSELT MAX   NUSSELT MIN' ;MESS 'CAST3M : Gauche     ' NUMG '  ' (MAXI Nusg) '  ' (MINI Nusg)    ;MESS 'CAST3M : Droite     ' NUMD '  ' (MAXI Nusd) '  ' (MINI Nusd)    ;Si (EGA Ra 1.e5);MESS 'Vahl Davis : Ra 1.e6    4.519              7.717         0.729';Finsi ;Si (EGA Ra 1.e6);MESS 'Vahl Davis : Ra 1.e6    8.8               17.925         0.989';Finsi ;MESS '---------------------------------------------------------------';MESS '                         Umax            Vmax   ' ;MESS 'CAST3M      :          ' (MAXI (exco un 'UX')) (MAXI (exco un 'UY'));Si (EGA Ra 1.e5);MESS 'Vahl Davis  : Ra 1.e6    34.73            68.59  ' ;Finsi ;Si (EGA Ra 1.e6);MESS 'Vahl Davis  : Ra 1.e6    64.63           219.36  ' ;Finsi ;MESS '---------------------------------------------------------------';MESS '                          Psi max         ' ;MESS 'CAST3M      :           '  (MAXI PSI)       ;Si (EGA Ra 1.e5);MESS 'Vahl Davis  : Ra 1.e6      9.612          ' ;Finsi;Si (EGA Ra 1.e6);MESS 'Vahl Davis  : Ra 1.e6     16.75           ' ;Finsi;MESS '---------------------------------------------------------------';opti 'ECHO' 1 ;  Si (NON COMPLET) ; ern=abs( NUMD - 7.835)/7.835 ; ernm=abs((MAXI Nusg) - 15.6)/15.6 ;eru=abs((MAXI (exco un 'UY')) - 238.5)/238.5;erp=abs((MAXI PSI) - 31.7)/31.7;mess 'ern = ' ern ' ernm = ' ernm ' eru = ' eru  ' erp = ' erp ;Si (ern  > 0.005); erreur 5 ; Finsi ;Si (ernm > 0.005); erreur 5 ; Finsi ;Si (eru  > 0.005); erreur 5 ; Finsi ;Si (erp  > 0.005); erreur 5 ; Finsi ; Finsi ; FIN ;

© Cast3M 2003 - Tous droits réservés.
Mentions légales