* fichier :  aerosol3.dgibi
************************************************************************
* Section : Fluides Transitoire
************************************************************************

GRAPH = FAUX ;

************************************************************************
*                          AEROSOL3.DGIBI                              *
*                                                                      *
*    Exemple d'utilisation des FONCTIONS DE PAROI AEROSOL LAMINAIRES.  *
*    Ce jeu de données teste les operateurs TSCA, ECHI, KUET et FPAL.  *
*    On résoud une équation de transport de particules dans un         *
*    ecoulement 2D plan de Poiseuille donne, avec dépot en paroi.      *
*                                                                      *
*    resolution semi-implicite d'une equation de concentration avec    *
*    depot d'aérosols en paroi                                         *
*                                                                      *
*P. CORNET SEMT/TTMF DECEMBRE 1997                                     *
************************************************************************

*------------------------- maillage ------------------------------------

TITRE 'MAILLAGE POUR DEPOT PAR DIFFUSION' ;
OPTIO DIME 2 ELEM QUA4    ;                                                     

L = 1. ;
H = 0.1  ;

nx = 30  ;
nz = 10  ;
E = H/nx ;
                                                                                
* POINTS :                                                                      
                                                                                
PA1 = 0. 0.    ;
PA2 = 0. (E-h) ;
PA3 = 0. (0.-H);
                                                                                
PB1 = L   0.   ;
PB2 = L  (e-h) ;
PB3 = L  (0.-H);
                                                                                
PC1 = (L/2.)   0.   ;
PC2 = (L/2.)  (e-h) ;
PC3 = (L/2.)  (0.-H);
                                                                                
* SEGMENTS DE BASE :                                                            
                                                                                
DA1 = PA1 DROIT nz PA3 ;
DB1 = PB1 DROIT nz PB3 ;
DC1 = PC1 DROIT nz PC3 ;
                                                                                
DAC3 = PA3 DROIT (nx/2) PC3 ;
DAC1 = PA1 DROIT (nx/2) PC1 ;
DAC2 = PA2 DROIT (nx/2) PC2 ;

DCB3 = PC3 DROIT (nx/2) PB3 ;
DCB1 = PC1 DROIT (nx/2) PB1 ;
DCB2 = PC2 DROIT (nx/2) PB2 ;

DSOL = DAC3 ET DCB3 ;
DAXE = DAC1 ET DCB1 ;
DFLU = DAC2 ET DCB2 ;
SOL2 = DSOL ELEM (LECT 2 PAS 1 NX) ;

DENT = DA1   ;
DSOR = DB1   ;
DMIL = DC1   ;
ENT2 = DENT ELEM (LECT 1 PAS 1 (NZ-1)) ;

* DOMAINES
                                                                                
DOM1 = DA1 DAC3 (INVE DC1) (INVE DAC1) DALLER PLAN ;
DOM2 = DC1 DCB3 (INVE DB1) (INVE DCB1) DALLER PLAN ;

DOMTOT = DOM1 ET DOM2 ;

BORTOT = CONTOUR DOMTOT ;                                                       
                                                                                
*---------------------- MODELES et normales aux faces -----------

MDOMTOT = CHAN DOMTOT QUAF ;
Mdsor   = CHAN dsor   QUAF ;
Mdmil   = CHAN dmil   QUAF ;
Mdent   = CHAN dent   QUAF ;
Mdsol   = CHAN dsol   QUAF ;
Msol2   = CHAN sol2   QUAF ;
ELIM (MDOMTOT et Mdsor et Mdmil et Mdent et Mdsol et Msol2) 1.e-5 ;

$domtot =mode MDOMTOT 'NAVIER_STOKES' LINE ;
DOMA $domtot 'IMPR' ;
$dsor   =mode Mdsor 'NAVIER_STOKES' LINE ;
$dmil   =mode Mdmil 'NAVIER_STOKES' LINE ;
$dent   =mode Mdent 'NAVIER_STOKES' LINE ;
$dsol   =mode Mdsol 'NAVIER_STOKES' LINE ;
$sol2   =mode Msol2 'NAVIER_STOKES' LINE ;

NOR = DOMA  $DOMTOT 'NORMALE' ;
                                                                                
*------------------------ Donnees --------------------------------------*
*                                                                       *
* ECOULEMENT : UM  = vitesse debitante (m/s)                            *
*              NU  = viscosite cinematique (m2/s)                       *
*              ROF = masse vomumique du gaz (kg/m3)                     *
*                                                                       *
* PARTICULES : C0  = concentration initiale (part/m3)                   *
*              CP  = concentration a la paroi (part/m3)                 *
*              ROG = gravite x masse volumique des particules (kg/m2s2) *
*              DP  = diametre des particules (m)                        *
*              DIF = coefficient de diffusion brownienne (m2/s)         *
*                                                                       *
*-----------------------------------------------------------------------*

R0   = H      ;
L0   = L      ;

UM  = 7.5E-3  ;
NU  = 1.5E-5  ;
ROF = 1.2     ;

ROG = 0. (-9810.)  ;
DP  = 1.E-7   ;
DIF = 6.8E-10 ;
RAP = DP/2.   ;

C0   = 1.  ;
CP   = 0.  ;

*----------------- INITIALISATION champ de vitesse ---------------------

cy=coor 2 domtot ;
cy=nomc 'UX' cy ;
cy=abs cy ;
VX= kcht $domtot scal sommet comp 'UX' (1.5*UM*( 1.-((cy/R0)**2.)) ) ;
VY= kcht $domtot scal sommet comp 'UY' 0. ;

VN= kcht $domtot vect sommet comp (mots 'UX' 'UY') (vx et vy) ;


*-------------------- calcul vitesse de depot des aerosols -------------

UET = KUET  NU VN NOR  $DOMTOT $SOL2  ;

AK = FPAL  NU ROF UET NOR ROG RAP $SOL2  ;


*------------------------- Equation de CONCENTRATION -------------------

CALCUL = EQEX $domtot 'ITMA' 100    'ALFA' 0.9  'TFINAL'  50.

  OPTI 'EFM1' 'EXPL' 'SUPGDC'
  ZONE $DOMTOT  'OPER' 'TSCA'  DIF 'VN'  0.    'INCO' 'CN'
  ZONE $SOL2    'OPER' 'ECHI'  AK   CP         'INCO' 'CN'
  OPTI EFM1 'CENTREE'
  ZONE $DOMTOT  'OPER' 'DFDT' 1. 'CN' 'DELTAT' 'INCO' 'CN'

'CLIM' 'CN' TIMP DENT C0 ;

*--------------------------- Initialisations et historiques ------------

calcul.'INCO' = table 'INCO' ;
calcul.'INCO'.'VN'= vn ;
calcul.'INCO'.'CN' = kcht $domtot scal sommet c0 ;

lh = PB1 et PB2 et PB3 et PC1 et PC2 et PC3 ;
his = khis 'CN'   lh  ;
calcul.'HIST'=his ;


*--------------------------- EXECUTION ---------------------------------

EXEC CALCUL  ;

*------------------------- DESSINS -------------------------------------
SI GRAPH ;

  VNCH = VECTEUR VN 10.  UX UY VERT ;
  TITRE 'depot laminaire : VITESSES ';
  TRACE VNCH DOMTOT (BORTOT) ;

  CN = calcul.'INCO'.'CN' ;
  opti isov ligne ;
  titre 'CONCENTRATION' ;
  trace cn domtot (bortot) ;

  dessin his.'TABD' his.'CN' ;

  CMIL = KCHT $DMIL scal sommet CN ;
  CSOR = KCHT $DSOR scal sommet CN ;
  CSOL = KCHT $DSOL scal sommet CN ;

  TITRE ' CONCENTRATION MILIEU ';
  evoc1 = evol chpo CMIL      DMIL ; dess evoc1 ;

  TITRE ' CONCENTRATION SORTIE ';
  evoc2 = evol chpo CSOR      DSOR ; dess evoc2 ;

  TITRE ' CONCENTRATION PAROI  ';
  evoc3 = evol chpo CSOL      DSOL ; dess evoc3 ;

FINSI ;

*-------------------- test sur le flux deposé total --------------------

SUR = DOMA $SOL2 'VOLUME' ;
KV  = SUR*AK ;
CN1 = KCHT $SOL2 scal sommet CALCUL.INCO.'CN' ;
CE2 = NOEL $SOL2 CN1 ;
CKS = CE2*KV ;
CKS = SOMT CKS ;
DEPREL = ( CKS - 1.26289E-06 ) / 1.26289E-06 ;
mess 'FLUX DEPOSE =' CKS ;

SI ( (ABS DEPREL) > 0.05 ) ;
     ERREUR 5 ;
FINSI ;

FIN ;