Accueil | Cast3M

*****************************************************
************************************************************************
************************************************************************
* fichier :  vahldavis3D.dgibi                      *
** modifie le 15/06/2014 passage EQPR -> EQEX       *
*****************************************************
 
GRAPH = FAUX ;
NITER = 250 ;
CFL = 0.9 ;
DISCR='MACRO';
KPRESS='CENTRE';
BETA=1.;
BETA = 100. ;
 
************************
* CAVITE CARREE 3D     *
* Paillere TTMF Dec 97 *
************************
 
******************************
*ESTIMATION DE LA CONVERGENCE*
******************************
 
DEBPROC CALCUL;
ARGU RVX*'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);
unz = kcht $MT scal sommet (exco 'UZ' un);
unm1z = kcht $MT scal sommet (exco 'UZ' unm1);
 
ERRX = KOPS unx '-' unm1x;
ERRY = KOPS uny '-' unm1y;
ERRZ = KOPS unz '-' unm1z;
ELIX = MAXI ERRX 'ABS';
ELIY = MAXI ERRY 'ABS';
ELIZ = MAXI ERRZ 'ABS';
ELIX = (LOG (ELIX + 1.0E-10))/(LOG 10.);
ELIY = (LOG (ELIY + 1.0E-10))/(LOG 10.);
ELIZ = (LOG (ELIZ + 1.0E-10))/(LOG 10.);
MESSAGE 'ITER' RV.PASDETPS.'NUPASDT'   'ERREUR LINF' ELIX ELIY ELIZ ;
RV.INCO.'UNM1' = KCHT $MT 'VECT' 'SOMMET' (RV.INCO.'UN');
IT = PROG RV.PASDETPS.'NUPASDT';
ERX = PROG ELIX ;
ERY = PROG ELIY ;
ERZ = PROG ELIZ ;
RV.INCO.'IT' = (RV.INCO.'IT') ET IT ;
RV.INCO.'ERX' = (RV.INCO.'ERX') ET ERX ;
RV.INCO.'ERY' = (RV.INCO.'ERY') ET ERY ;
RV.INCO.'ERZ' = (RV.INCO.'ERZ') ET ERZ ;
 
FINSI;
as2 ama1  = 'KOPS' 'MATRIK' ;
FINPROC as2 ama1 ;
 
************
* MAILLAGE *
************
 
OPTI DIME 3;
OPTI ELEM CU20 ;
 
NX = 4 ;
NY = 2 ;
NZ = 4 ;
 
A1 = 0.0 0.0 0.0 ;
A2 = 1.0 0.0 0.0 ;
A3 = 1.0 0.4 0.0 ;
A4 = 0.0 0.4 0.0 ;
B1 = 0.0 0.0 1.0 ;
B2 = 1.0 0.0 1.0 ;
B3 = 1.0 0.4 1.0 ;
B4 = 0.0 0.4 1.0 ;
 
BAS = 'DALLER' (A1 D NX A2) (A2 D NY A3) (A3 D NX A4) 
      (A4 D NY A1) 'PLAN' ;
HAU = 'DALLER' (B1 D NX B2) (B2 D NY B3) (B3 D NX B4) 
      (B4 D NY B1) 'PLAN' ;
DRO = 'DALLER' (A2 D NY A3) (A3 D NZ B3) (B3 D NY B2) 
      (B2 D NZ A2) 'PLAN' ;
GAU = 'DALLER' (A1 D NY A4) (A4 D NZ B4) (B4 D NY B1) 
      (B1 D NZ A1) 'PLAN' ;
DEV = 'DALLER' (A1 D NX A2) (A2 D NZ B2) (B2 D NX B1)
      (B1 D NZ A1) 'PLAN' ;
DER = 'DALLER' (A4 D NX A3) (A3 D NZ B3) (B3 D NX B4)
      (B4 D NZ A4) 'PLAN' ;
 
BAS = BAS COUL TURQ ;
HAU = HAU COUL TURQ ;
DRO = DRO COUL BLEU ;
GAU = GAU COUL ROUG ;
DEV = DEV COUL VERT ;
DER = DER COUL VERT ;
 
V1 = 0.0 0.0 1.0 ;
MT = BAS 'VOLU' NZ 'TRAN' V1 ;
MMT= CHAN MT QUAF ;
 
Mbas = chan bas quaf ;
Mhau = chan hau quaf ;
Mdro = chan dro quaf ;
Mgau = chan gau quaf ;
Mdev = chan dev quaf ;
Mder = chan der quaf ;
 
elim (MMT et Mbas et Mhau et Mdro et Mgau et Mdev et Mder) 0.001;
 
$MT = 'MODE' MMT 'NAVIER_STOKES' DISCR ;
MT = DOMA $MT MAILLAGE ;
$HAU = 'MODE' MHAU 'NAVIER_STOKES' DISCR ;
HAU = DOMA $HAU MAILLAGE ;
$BAS = 'MODE' MBAS 'NAVIER_STOKES' DISCR ;
BAS = DOMA $BAS MAILLAGE ;
$GAU = 'MODE' MGAU 'NAVIER_STOKES' DISCR ;
GAU = DOMA $GAU MAILLAGE ;
$DRO = 'MODE' MDRO 'NAVIER_STOKES' DISCR ;
£DRO = 'DOMA' $DRO 'TABLE' ;
DRO = DOMA $DRO MAILLAGE ;
$DEV = 'MODE' MDEV 'NAVIER_STOKES' DISCR ;
DEV = DOMA $DEV MAILLAGE ;
$DER = 'MODE' MDER 'NAVIER_STOKES' DISCR ;
DER = DOMA $DER MAILLAGE ;
 
CNT = ENVE MT ;
 
**************
* PARAMETRES *
**************
 
Pr = 0.71D0 ;
Ra = 1.D5 ;
Gr = Ra/Pr ;
NU = 1.D0/(Gr**0.5) ;
ALF= NU/Pr ;
gb = 0.D0 0.D0 -1.0D0 ;
 
***********************
* CREATION DES TABLES *
***********************
 
RV = EQEX $MT 'ITMA' NITER 'ALFA' CFL
 'ZONE' $MT 'OPER' CALCUL
 'ZONE' $MT 'OPER' 'NS' NU GB 'TN' 0.5 'INCO' 'UN'
 'ZONE' $MT 'OPER' 'TSCAL' ALF 'UN' 0. 'INCO' 'TN'
 'OPTI' 'CENTREE'
  ZONE  $MT 'OPER' 'DFDT' 1. 'UN' 'DELTAT' INCO 'UN'
  ZONE  $MT 'OPER' 'DFDT' 1. 'TN' 'DELTAT' INCO 'TN'
;
 
RV = EQEX RV
 'CLIM' 'UN' 'UIMP' (HAU ET BAS ET DRO ET GAU) 0.
        'UN' 'VIMP' cnt 0.
        'UN' 'WIMP' (HAU ET BAS ET DRO ET GAU) 0.
        'TN' 'TIMP' gau 1.
        'TN' 'TIMP' dro 0. ;
 
RVP = EQEX 'OPTI' 'EF' KPRESS
 'ZONE' $MT OPER  KBBT  -1. beta  INCO 'UN' 'PRES'
;
 
    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.'PROJ'     =RVP ;
 
 RV.INCO = TABLE INCO;
 RV.INCO.'PRES' = kcht $MT SCAL KPRESS 0.;
 RV.INCO.'UN' = kcht $MT VECT SOMMET (0. 0. 0.);
 corx = coor 1 mt;
 RV.INCO.'TN' = kcht $MT SCAL SOMMET (1.-corx);
 
 RV.INCO.'UNM1' = kcht $MT VECT SOMMET (1.D-5 1.D-5 1.D-5);
 RV.INCO.'IT' = PROG 1;
 RV.INCO.'ERX' = PROG 0.;
 RV.INCO.'ERY' = PROG 0.;
 RV.INCO.'ERZ' = PROG 0.;
 
 
EXEC RV ;
 
DYDZ = DOMA $DRO 'VOLUME' ;
DYDZT = SOMT DYDZ ;
 
GRADT = KOPS RV.INCO.'TN' 'GRAD' $MT;
DTDX = KCHT $MT 'SCAL' 'CENTRE' (EXCO 'UX' GRADT);
DTDX = ELNO $MT DTDX;
 
DTDXp = KCHT $DRO 'SCAL' 'SOMMET' DTDX;
DTDXp = KOPS DTDXp '*' (-1.);
 
DTDXpc = NOEL $DRO DTDXp;
num = KOPS DTDXpc '*' DYDZ ;
num = SOMT NUM;
num = NUM/DYDZT;
 
 
MESSAGE 'NUSSELT MOYEN' NUM ' 2D ' 2.7124 ;
MESSAGE 'NUSSELT MAX' (MAXI DTDXpc) ' 2D ' 3.7240 ;
MESSAGE 'NUSSELT MIN' (MINI DTDXpc) ' 2D ' 0.93978 ;
 
ER1 = ABS( NUM - 2.7124D0) / 2.7124 ;
ER2 = ABS( (MAXI DTDXpc) - 3.7240) / 3.7240 ;
ER3 = ABS( (MINI DTDXpc) - 0.93978) / 0.93978 ;
ER4 = MINI (RV.INCO.'ERX') ;
ER5 = MINI (RV.INCO.'ERY') ;
ER6 = MINI (RV.INCO.'ERZ') ;
 
MESSAGE ER1 ER2 ER3 ER4 ER5 ER6  ;
 
SI ( ER4 > -2.0 ) ;
        ERREUR 5 ;
FINSI ;
SI ( ER5 > -2.0 ) ;
        ERREUR 5 ;
FINSI ;
SI ( ER6 > -2.0 ) ;
        ERREUR 5 ;
FINSI ;
SI ( ER1 > 0.04 ) ;
        ERREUR 5 ;
FINSI ;
 
SI ( GRAPH ) ;
 
TRACE MT 'CACH' ;
TRACE RV.INCO.'TN' MT CNT 14 ;
UNCH = VECT RV.INCO.'UN' 1. UX UY UZ VERT ;
TRACE UNCH MT CNT ;
 
EVOLX = EVOL 'MANU' 'ITERATIONS' (RV.INCO.'IT') 'LOG|E(u)|inf'
     (RV.INCO.'ERX') ;
EVOLY = EVOL 'MANU' 'ITERATIONS' (RV.INCO.'IT') 'LOG|E(v)|inf'
     (RV.INCO.'ERY') ;
EVOLZ = EVOL 'MANU' 'ITERATIONS' (RV.INCO.'IT') 'LOG|E(w)|inf'
     (RV.INCO.'ERZ') ;
EVOLX = EVOLX COUL ROUG ;
EVOLY = EVOLY COUL JAUN ;
EVOLZ = EVOLZ COUL VERT ;
TAB1 = TABLE ;
TAB1.1 = 'MARQ TRIA REGU' ;
TAB1.2 = 'MARQ LOSA REGU' ;
TAB1.3 = 'MARQ CARR REGU' ;
TAB1.'TITRE' = TABLE ;
TAB1.'TITRE' . 1 = 'MOT' 'UX' ;
TAB1.'TITRE' . 2 = 'MOT' 'UY' ;
TAB1.'TITRE' . 3 = 'MOT' 'UZ' ;
dess (EVOLX ET EVOLY ET EVOLZ) 'XBOR' 0. (2*NITER) 'YBOR' -10.0 0.0
        'TITR' 'CONVERGENCE VERS LE STATIONNAIRE' LEGE TAB1 ;
 
FINSI ;
 
FIN ;