Download blasius.dgibi

Back to the list

* fichier :  blasius.dgibi
COMPLET = FAUX ;
SI ( COMPLET ) ;
      d1= 0.5 ; d2= 0.02 ;
      nbit = 5000 ;
      err1 = 7.5E-3 ;
SINON ;
      d1 = 1.0 ; d2 = 0.1 ;
      nbit = 150 ;
      err1 = 7.E-2 ;
FINSI ;
GRAPH = 'N' ;

************************************************************************
* COUCHE LIMITE LAMINAIRE SUR PLAQUE PLANE - COMPARAISON AVEC SOLUTION *
* THEORIQUE DE BLASIUS                                                 *
************************************************************************

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) ;

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 *
************

option dime  2 elem qua8 ;

p1 = -1.5 0.0 ;
p2 =  0.0 0.0 ;
p3 =  1.0 0.0 ;
p4 =  1.0 0.2 ;
p5 =  1.0 1.0 ;
p6 =  0.0 1.0 ;
p7 = -1.5 1.0 ;
p8 = -1.5 0.2 ;

bas1 = p1 d dini d1 dfin d2 p2 ;
bas2 = p2 d dini d2 dfin d2 p3 ;
bas = bas1 et bas2 ;
dro = p3 d dini d2 dfin d2 p4 d dini d2 dfin d1 p5 ;
hau = p5 d dini d2 dfin d2 p6 d dini d2 dfin d1 p7 ;
gau = p7 d dini d1 dfin d2 p8 d dini d2 dfin d2 p1 ;

mt= bas dro hau gau daller ;

tass MT ;
orienter mt ;

mt = chan mt quaf ;
dro= chan dro quaf ;
bas = chan bas quaf ;
gau = chan gau quaf ;
hau = chan hau quaf ;

elim (mt et dro et bas et gau et hau) 1.e-3 ;
$mt = mode mt  'NAVIER_STOKES' MACRO ; doma $mt 'IMPR' ;
$dro = mode dro  'NAVIER_STOKES' MACRO ;

mt2 = mt ;
mt = doma $mt maillage ;

Nu = 1./2000. ;
tos = 0. 0. ;

RV = EQEX $MT  'ITMA' nbit 'ALFA' 0.9
 'OPTI' 'SUPG'
 'ZONE' $MT 'OPER' CALCUL
 'ZONE' $MT 'OPER' 'NS' NU 'INCO' 'UN' 
 'OPTI' 'CENTREE'  'EFM1'
 'ZONE' $MT  'OPER' 'DFDT'  1. 'UN' 'DELTAT  ' INCO 'UN'
 'ZONE' $DRO 'OPER' 'TOIM' tos 'INCO' 'UN' ;
  
RV = EQEX RV
 'CLIM' 'UN' 'UIMP' bas2 0. 
        'UN' 'UIMP' gau 1.0
        'UN' 'VIMP' bas 0.0 ;

RVP = EQPR $MT KTYPI 1
         ZONE $MT  'OPER' 'PRESSION' 0.  ;

RV.INCO=TABLE INCO ;
RV.INCO.'UN'=kcht $MT VECT SOMMET (1. 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.'NBIT' = 0 ;
RV.PRESSION =RVP ;

EXEC RV ;

* SOLUTION THEORIQUE BLASIUS *

LIYT = PROG
   0.000000 0.003162 0.006325 0.009487 0.012649 0.015811 0.018974
   0.022136 0.025298 0.028460 0.031623 0.034785 0.037947 0.041110
   0.044272 0.047434 0.050596 0.053759 0.056921 0.060083 0.063246
   0.069570 0.075895 0.082219 0.088544 0.094868 0.101193 0.107517
   0.113842 0.120167 0.126491 0.132816 0.139140 0.145465 0.151789
   0.158114 0.164438 0.170763 0.177088 0.183412 0.189737 1.000000 ;
LIUT = PROG
   0.000000 0.046960 0.093910 0.140810 0.187610 0.234230 0.280580
   0.326530 0.371960 0.416720 0.460630 0.503540 0.545250 0.585590
   0.624390 0.661470 0.696700 0.729930 0.761060 0.790000 0.816690
   0.863300 0.901070 0.930600 0.952880 0.969050 0.980370 0.987970
   0.992890 0.995940 0.997770 0.998820 0.999400 0.999700 0.999860
   0.999940 0.999971 0.999988 0.999995 0.999998 0.999999 1.000000 ;
LIVT = PROG
   0.000000 0.000037 0.000149 0.000334 0.000593 0.000925 0.001327
   0.001800 0.002338 0.002940 0.003599 0.004312 0.005070 0.005867
   0.006695 0.007545 0.008408 0.009275 0.010136 0.010982 0.011804
   0.013350 0.014721 0.015886 0.016836 0.017576 0.018129 0.018525
   0.018797 0.018975 0.019088 0.019155 0.019195 0.019216 0.019228
   0.019234 0.019237 0.019238 0.019239 0.019239 0.019239 0.019239 ;

EVOLU = EVOL 'MANU' 'Y' LIYT 'U(1.,y)' LIUT ;
EVOLV = EVOL 'MANU' 'Y' LIYT 'V(1.,y)' LIVT ;

YY = COOR 2 (DOMA $MT SOMMET) ;
EVOL1 = EVOL 'CHPO' RV.INCO.'UN' UX DRO ;
EVOL2 = EVOL 'CHPO' RV.INCO.'UN' UY DRO ;
EVOL3 = EVOL 'CHPO' YY SCAL DRO ;
LIY = 'EXTR' EVOL3 ORDO ;
LIU = 'EXTR' EVOL1 ORDO ;
LIV = 'EXTR' EVOL2 ORDO ;

LIUC = IPOL LIYT LIY LIU ;
NP = DIME LIYT ;
ERR0 = 0. ;
REPETER BLOC1 NP ;
    UCAL = EXTRAIRE LIUC &BLOC1 ;
    UREF = EXTRAIRE LIUT &BLOC1 ;
    ERR0 = ERR0 + ((UCAL-UREF)*(UCAL-UREF)) ;
FIN BLOC1 ;
ERR0 = ERR0/NP ;
ERR0 = ERR0 '**' 0.5 ;
MESSAGE 'ERREUR ' ERR0 ;
SI ( ERR0 > ERR1 ) ;
    ERREUR 5 ;
FINSI ;

SI ('EGA' graph 'O') ;
opti isov lign ;
pn = elno $MT RVP.PRESSIOn ;
TRACE Pn mt2 ;
unch = vect RV.INCO.'UN' 0.1 UX UY JAUNE ;
TRACE UNCH MT2 ;
TAB1 = TABLE ;
* TAB1.1 = 'MOT' 'TITR SOLUTION_CASTEM2000' ;
TAB1.2 = 'MARQ LOSA NOLI ' ;
TAB1.'TITRE' = TABLE ;
TAB1.'TITRE' . 1  = 'MOT' 'SOLUTION_CASTEM2000' ;
TAB1.'TITRE' . 2  = MOT 'SOLUTION_EXACTE' ;
DESS (EVOL1 ET EVOLU) 'LEGE' 'TITY' 'UX' 'TITR' 'U EN SORTIE' TAB1 ;
DESS (EVOL2 ET EVOLV) 'LEGE' 'TITY' 'UY' 'TITR' 'V EN SORTIE' TAB1 ;
FINSI ;

fin ;

 

 

 

 

© Cast3M 2003 - All rights reserved.
Disclaimer