* fichier : ns_ouvert.dgibi
'OPTI' 'DIME' 3 'ELEM' 'QUA4' ;

*** CAS-test permettant de vérifier la conservation des
*** débits volumiques dans le cas d'un ouvert
*** Date de création : Juillet 2007 (A. Bleyer)
*** correction du bug dans normno.eso

'DEBP' EXECPDA PDA*'TABLE' NBIT*'ENTIER' ;

'SI' ('NON' ('EXIST' PDA 'TIC')) ;
  PDA.'TIC' = 'TABLE' 'INCO' ;
  PDA.'TIC'.'LTPS'  = 'PROG' 0.0 ;
  PDA.'TIC'.'NUPAT' = 0 ;
  PDA.'TIC'.'NU' = 1.8E-5   ;
  PDA.'TIC'.'Uin'   = 'PROG' 0.0 ;
'FINSI' ;

*** decentrement des inconnues
KSUPG = 'SUPG' ;
PDA.'KSUPG' = KSUPG ;
*** discretisation de la pression
KPRES = 'MSOMMET' ;
'SI' ('EXIST' PDA 'KPRES') ;
KPRES = PDA.'KPRES' ;
'FINSI' ;
PDA.'KPRES' = KPRES ;
*** discretisation de la vitesse ou des champs transportés
DISCR = 'LINE' ;
'SI' ('EXIST' PDA 'DISCR') ;
DISCR = PDA.'DISCR' ;
'FINSI' ;
PDA.'DISCR' = DISCR ;
*** freqence de re-preconditionnement
'SI' ('NON' ('EXIST' PDA 'FRPREC')) ;
PDA.'FRPREC' = 10 ;
'FINSI' ;
*** fonction de paroi + conditions de glissement
FPAROI = FAUX ;
'SI' ('EXIST' PDA 'FPAROI') ;
FPAROI = PDA.'FPAROI' ;
'FINSI' ;
PDA.'YP' = 1.E-2 ;


*----------------------------------------------------------------------*
*          Remplissage de la table GEO                                 *
*----------------------------------------------------------------------*
'SI' ('NON' ('EXIST' PDA 'GEO')) ;
  GEO = 'TABLE' ;

  'SI' ('NON' ('EXIST' PDA 'EPSI')) ;
    PDA.'EPSI' = 1.E-5 ;
  'FINSI' ;

  M_VTF    = 'CHAN' PDA.'VTF'    'QUAF' ;
  M_INLET  = 'CHAN' PDA.'INLET'  'QUAF' ;
  M_OUTLET = 'CHAN' PDA.'OUTLET' 'QUAF' ;
  'ELIM' (M_VTF 'ET' M_INLET 'ET' M_OUTLET) PDA.'EPSI' ;

  $VTF     = 'MODE' M_VTF    'NAVIER_STOKES' DISCR ;
  VTotal      = 'SOMT' ('DOMA' $VTF    'VOLUME')       ;
  GEO.'M_ENVF'  = 'DOMA' $VTF 'ENVELOPP' ;

  'ELIM' (GEO.'M_ENVF' 'ET' M_INLET 'ET' M_OUTLET) PDA.'EPSI' ;
  GEO.'$ENVF'   = 'MODE' (GEO.'M_ENVF') 'NAVIER_STOKES' DISCR ;
  GEO.'ENVF'    = 'DOMA' (GEO.'$ENVF')  'MAILLAGE' ;

  M_INLET  = 'INTE' GEO.'M_ENVF' M_INLET  ;
  M_OUTLET = 'INTE' GEO.'M_ENVF' M_OUTLET ;

  $INLET   = 'MODE' M_INLET  'NAVIER_STOKES' DISCR ;
  $OUTLET  = 'MODE' M_OUTLET 'NAVIER_STOKES' DISCR ;

  M_PEAUI  = 'DIFF' GEO.'M_ENVF' (M_INLET 'ET' M_OUTLET) ;
  $PEAUI   = 'MODE' M_PEAUI  'NAVIER_STOKES' DISCR ;

  GEO.'VTF'     = 'DOMA' $VTF    'MAILLAGE' ;
  GEO.'INLET'   = 'DOMA' $INLET  'MAILLAGE' ;
  GEO.'OUTLET'  = 'DOMA' $OUTLET 'MAILLAGE' ;

  GEO.'$VTF'    = $VTF    ;
  GEO.'$INLET'  = $INLET  ;
  GEO.'$OUTLET' = $OUTLET ;

  GEO.'M_VTF'    = M_VTF    ;
  GEO.'M_INLET'  = M_INLET  ;
  GEO.'M_OUTLET' = M_OUTLET ;

  GEO.'SINLET'  = 'MAXI' ('RESU' ('DOMA' $INLET  'XXDIAGSI')) ;
  GEO.'SOUTLET' = 'MAXI' ('RESU' ('DOMA' $OUTLET 'XXDIAGSI')) ;

  'SI' ('EGA' ('VALE' 'DIME') 2) ;
  GEO.'P0' = 0.0 0.0 ;
  IN1 = 'CHAN' GEO.'INLET' 'POI1' ;
  IN2 = 'POINT' IN1 'INITIAL' ;
  IN3 = 'POINT' IN1 'FINAL' ;
  GEO.'IN_C' = IN2 'ET' IN3 ;
  IO1 = 'CHAN' GEO.'OUTLET' 'POI1' ;
  IO2 = 'POINT' IO1 'INITIAL' ;
  IO3 = 'POINT' IO1 'FINAL' ;
  GEO.'OUT_C' = IO2 'ET' IO3 ;
  'SINON' ;
  GEO.'P0' = 0.0 0.0 0.0 ;
* entree
  GEO.'M_IN_C'   = 'DOMA' (GEO.'$INLET') 'ENVELOPPE' ;
  GEO.'$IN_C'    = 'MODE' (GEO.'M_IN_C') 'NAVIER_STOKES' DISCR ;
  GEO.'IN_C'     = 'DOMA' (GEO.'$IN_C')  'MAILLAGE' ;
* sortie
  GEO.'M_OUT_C'  = 'DOMA' (GEO.'$OUTLET') 'ENVELOPPE' ;
  GEO.'$OUT_C'   = 'MODE' (GEO.'M_OUT_C') 'NAVIER_STOKES' DISCR ;
  GEO.'OUT_C'    = 'DOMA' (GEO.'$OUT_C')  'MAILLAGE' ;
  'FINSI' ;

  GEO.'IN_I'     = 'DIFF' ('CHAN' (GEO.'INLET') 'POI1')
                   ('CHAN' (GEO.'IN_C') 'POI1') ;

  GEO.'OUT_I'    = 'DIFF' ('CHAN' (GEO.'OUTLET') 'POI1')
                   ('CHAN' (GEO.'OUT_C') 'POI1') ;

  GEO.'M_PEAUI'  = M_PEAUI  ;
  GEO.'$PEAUI'  = $PEAUI  ;
  GEO.'PEAUI'   = 'DOMA' $PEAUI  'MAILLAGE' ;

  'DOMA' (GEO.'$VTF') 'IMPR' ;

*  CALCUL DES FACTEURS GEOMETRIQUES POUR L'ENTREE ET LA SORTIE
  diru1=(PDA.'Dir_in')/('NORM' PDA.'Dir_in') ;
  'SI' ('EGA' ('VALE' 'DIME') 2) ;
  gjc = 'MANU' 'CHPO' (GEO.'IN_C') 2 'UX' 0.0 'UY' 0.0 ;
  'SINON' ;
  gjc = 'KCHT' (GEO.'$IN_C')   'VECT' 'SOMMET' GEO.'P0' ;
  'FINSI' ;
  gjx = 'KCHT' (GEO.'$INLET')  'VECT' 'SOMMET' diru1 gjc ;
  unvf0 = 'KCHT' (GEO.'$ENVF') 'VECT' 'SOMMET' gjx ;
  src = 'DBIT' unvf0 GEO.'$ENVF' ;
  GEO.'IN_FG' = 'ABS' (GEO.'SINLET' '*' (src**(-1.0))) ;

  diru1=(PDA.'Dir_out')/('NORM' PDA.'Dir_out') ;
  'SI' ('EGA' ('VALE' 'DIME') 2) ;
  gjc = 'MANU' 'CHPO' (GEO.'OUT_C') 2 'UX' 0.0 'UY' 0.0 ;
  'SINON' ;
  gjc = 'KCHT' (GEO.'$OUT_C')  'VECT' 'SOMMET' GEO.'P0' ;
  'FINSI' ;
  gjx = 'KCHT' (GEO.'$OUTLET') 'VECT' 'SOMMET' diru1 gjc ;
  unvf0 = 'KCHT' (GEO.'$ENVF') 'VECT' 'SOMMET' gjx ;
  src = 'DBIT' unvf0 GEO.'$ENVF' ;
  GEO.'OUT_FG' = 'ABS' (GEO.'SOUTLET' '*' (src**(-1.0))) ;

  PDA.'GEO' = GEO ;
'FINSI' ;

  GEO = PDA.'GEO' ;
  $PEAUI = GEO.'$PEAUI' ;

*----------------------------------------------------------------------*
*                    Definition des équations                          *
*----------------------------------------------------------------------*
'SI' ('NON' ('EXIST' PDA 'RV')) ;

*--------->            definition QDM             <--------------------*
* gravité dirigée suivant -z.

 RV = 'EQEX'  'OMEGA' 1. 'NITER' 1  'ITMA' 1
  'OPTI' 'EF' 'IMPL'   KSUPG            KPRES 'NODIV'
  'ZONE' (GEO.'$VTF') 'OPER' 'KONV' 1.0 'UN' 0.0 'INCO' 'UN'
  'OPTI' 'EF' 'IMPL'   KSUPG            KPRES 'NODIV'
  'ZONE' (GEO.'$VTF') 'OPER' 'LAPN' 'NUEFF' 'INCO' 'UN' ;

 RV = 'EQEX' RV
  'OPTI' 'EFM1' 'IMPL' 'CENTREE'
  'ZONE' (GEO.'$VTF') 'OPER' 'DFDT' 1.0 'UNM' 'DT' 'INCO' 'UN' ;

'SI' FPAROI ;
   RV = 'EQEX' RV
     'OPTI' 'EF' 'IMPL'   KSUPG  KPRES 'NODIV'
     'ZONE' $PEAUI 'OPER' 'FPU' 1.0 'UN' 'NU' 'UET' 'YP' 'INCO' 'UN' ;
'SINON' ;
  RV = 'EQEX' RV 'CLIM'
      'UN' 'UIMP' (GEO.'PEAUI') 0.
      'UN' 'VIMP' (GEO.'PEAUI') 0. ;
  'SI' ('EGA' ('VALE' 'DIME') 3) ;
  RV = 'EQEX' RV 'CLIM'
      'UN' 'WIMP' (GEO.'PEAUI') 0. ;
  'FINSI' ;
'FINSI' ;

  RV = 'EQEX' RV 'CLIM'
      'UN' 'UIMP' (GEO.'INLET') 0.
      'UN' 'VIMP' (GEO.'INLET') 0. ;

'SI' ('EGA' ('VALE' 'DIME') 3) ;
  RV = 'EQEX' RV 'CLIM'
      'UN' 'WIMP' (GEO.'INLET') 0. ;
'FINSI' ;

 RV.'METHINV'.TYPINV=3 ;
 RV.'METHINV'.IMPINV=0 ;
 RV.'METHINV'.NITMAX=300;
 RV.'METHINV'.PRECOND=3 ;
 RV.'METHINV'.ILUTLFIL=4 ;
 RV.'METHINV'.ILUTDTOL=0.;
 RV.'METHINV'.RESID  =1.e-10;
 RV.'METHINV'.'FCPRECT'=1   ;
 RV.'METHINV'.'FCPRECI'=1   ;


 RVP = 'EQEX'
    'OPTI' 'EF' KPRES
    'ZONE' (GEO.'$VTF') 'OPER' 'KBBT' (-1.) 20. 'INCO' 'UN' 'PRES' ;

 RVP = 'EQEX' RVP
    'OPTI' 'EFM1' 'INCOD' KPRES
    'ZONE' (GEO.'$VTF')
         'OPER' 'FIMP' 0.0 'INCO' 'PRES' ;

 'SI' FPAROI ;
    RVP = 'EQEX' RVP
    'OPTI' 'EF' 'IMPL' KPRES
    'ZONE' $PEAUI
       'OPER' VNIMP (GEO.'$VTF') 0.0 'INCO' 'UN' 'PRES' ;
 'FINSI' ;

 'SI' ('EGA' KPRES 'MSOMMET') ;
    RVP = 'EQEX' RVP
    'CLIM' 'PRES' 'TIMP' ('DOMA' GEO.'$OUTLET' 'MSOMMET') 0. ;
 'FINSI' ;

* TYPINV : type ENTIER
*          1 : résolution directe (Crout)
*          2 : Gradient Conjugué
*          3 : Bi-Gradient Conjugué Stabilisé (BiCGSTAB)
*          4 : BiCGSTAB(2) (expérimental)
*          5 : GMRES(m) : restarted Generalized Minimal Residual
'SI' ('EGA' DISCR 'MACRO') ;
    RVP.'METHINV'.TYPINV=3 ;
'SINON' ;
    RVP.'METHINV'.TYPINV=1 ;
'FINSI' ;
    RVP.'METHINV'.IMPINV=0 ;
    RVP.'METHINV'.NITMAX=300;
* PRECOND : type ENTIER
*          0 : pas de préconditionnement
*          1 : préconditionnement par la diagonale
*          2 : préconditionnement D-ILU
*          3 : préconditionnement ILU(0) (Crout)
*          4 : préconditionnement MILU(0) (Crout modifié)
*          5 : préconditionnement ILUT (dual truncation)
*          6 : préconditionnement ILUT2 (une variante du précédent)
    RVP.'METHINV'.PRECOND=3 ;
    RVP.'METHINV'.ILUTLFIL=10;
    RVP.'METHINV'.ILUTDTOL=0.;
    RVP.'METHINV'.RESID  =1.e-10 ;
    RVP.'METHINV'.'FCPRECT'=50000 ;
    RVP.'METHINV'.'FCPRECI'=50000 ;


    RV.'PROJ'= RVP ;


RV.'INCO' = PDA.'TIC' ;
PDA.'RV' = RV ;
** fin de la definition des equations
'FINSI' ;

*
'MESS' 'Discrétisation vitesse|pression ' DISCR ' | ' KPRES ;
'MESS' 'Décentrement utilisé : ' KSUPG ;
*

  PDA.'TIC'.'DT' = PDA.'DT0' ;
  RV = PDA.'RV' ;

** boucle temporelle
  'REPETER' BCLTPS NBIT ;

    'SI' (NBIT 'EGA' 0) ;
      'QUITTER' BCLTPS ;
    'FINSI' ;

    TPS0 = 'EXTR' (PDA.'TIC'.'LTPS') ('DIME' (PDA.'TIC'.'LTPS')) ;
    TPS0 = TPS0 + (PDA.'DT0') ;
    PDA.'TIC'.'LTPS' = (PDA.'TIC'.'LTPS') 'ET' ('PROG' TPS0) ;
    PDA.'TIC'.'NUPAT' = (PDA.'TIC'.'NUPAT') '+' 1 ;
    'MESS' ' Pas de temps ' (PDA.'TIC'.'NUPAT') ' Temps ' TPS0 ;

    Qint = 'IPOL' TPS0 PDA.'scenario'.'TIME' PDA.'scenario'.'Qinj' ;
    Uinj  = Qint / PDA.'GEO'.'SINLET' ;
    PDA.'TIC'.'Uin' = PDA.'TIC'.'Uin' 'ET' ('PROG' Uinj) ;

*    chercher Uinj en fonction de Qinj
    diru1=(PDA.'Dir_in')/('NORM' PDA.'Dir_in') ;
    gclim  matxx= 'KOPS' 'MATRIK' ;
    gj   ='KCHT' (GEO.'$INLET') 'VECT' 'SOMMET' 'COMP'
          ('MOTS' '1UN' '2UN' '3UN') (diru1*Uinj) ;
    'SI' ('NON' FPAROI) ;
      gj   = 'REDU' gj (GEO.'IN_I') ;
    'FINSI' ;
    gj   ='KCHT' (GEO.'$INLET') 'VECT' 'SOMMET' 'COMP'
          ('MOTS' '1UN' '2UN' '3UN') gj    ;

     RV.'CLIM'    = (RV.'CLIM' '-' RV.'CLIM') '+' gj ;

'SI' ('NON' ('EXIST' PDA.'TIC' 'UN')) ;
NCOMP = 'MOTS' 'UX' 'UY' ;
'SI' ('EGA' ('VALE' 'DIME') 3) ;
NCOMP = 'MOTS' 'UX' 'UY' 'UZ' ;
'FINSI' ;

gj2 = 'NOMC' gj ('EXTR' gj 'COMP') NCOMP
      'NATU' 'DISCRET' ;

PDA.'TIC'.'UN'  = 'KCHT' (GEO.'$VTF')
                  'VECT' 'SOMMET' (GEO.'P0') gj2 ;
PDA.'TIC'.'UNM' = 'COPIER' PDA.'TIC'.'UN' ;
PDA.'TIC'.'PRES'= 'KCHT' (GEO.'$VTF')
                  'SCAL' (PDA.'KPRES') 0.0;
PDA.'TIC'.'YP' = PDA.'YP' ;
'FINSI' ;

  UN  = PDA.'TIC'.'UN' ;
  P   = PRODT UN (PDA.'GEO'.'$VTF') ;
  LM  =PDA.'PARAM_TURB' ;
  PDA.'TIC'.'PRODT'=P;
  NUT = (LM*LM*(P**0.5)) ;
*Réalisabilité sur Nut
  NU  = PDA.'NU' ;
  a   = 5000.* NU ;
  al  =0.8 ; ala=al*a ; b=ala*((2.*(1.-al))**(-1.)) ;
  ik  =masq Nut 'INFERIEUR' ala ;
  NUT =(ik*NUT)+((1.-ik)*a*(NUT + b)*(inve (a + NUT + b))) ;
  Nutsnu=NUT / NU ;
  MESS ' Mini Maxi Nut/Nu ' ('MINI' Nutsnu) ('MAXI' Nutsnu) ;
  PDA.'TIC'.'NUEFF' ='KCHT' (PDA.'GEO'.'$VTF') 'SCAL' 'SOMMET'
                   (NU + NUT) ;

*---> re-preconditionnement des systemes toutes les 10 iterations
    'SI' ('MULT' PDA.'TIC'.'NUPAT' PDA.'FRPREC') ;
      RV.'calprec'  = VRAI ;
      'MESS' '==> On repreconditionne le systeme (QDM)' ;
    'FINSI' ;

      EXEC RV ;

  'FIN' BCLTPS ;

'FINP' ;





EPSI = 1.E-5 ;

Rinj = 0.1 ;
Linj = 0.2 ;
Rpip = 0.3 ;
Lpip = 0.5 ;
Hpip = 4.0 ;
Rcuv = 2.0 ;
Hcuv = 8.0 ;

OO = 0.0 0.0 0.0 ;
OX = 1.0 0.0 0.0 ;
OY = 0.0 1.0 0.0 ;
OZ = 0.0 0.0 1.0 ;

PX0 = 0.0  0.0 0.0 ;
PX1 = Rinj 0.0 0.0 ;
PX2 = Rcuv 0.0 0.0 ;

DX1 = 'DROIT' PX0 PX1 1 ;
DX2 = 'DROIT' PX1 PX2 5 ;
DZ0 = DX1 'MOINS' (Linj * OZ) ;
DZ1 = DX1 ET DX2 ;
DZ2P = 'DROIT' PX0 PX2 ('NBEL' DZ1) ;
DZ2 = DZ2P 'PLUS' ((Hpip - Rpip) * OZ) ;
DZ3 = DZ2P 'PLUS' ((Hpip + Rpip) * OZ) ;
DZ4 = DZ2P 'PLUS' (Hcuv * OZ) ;

SZ1 = 'REGL' DZ0 DX1 2 ;
SZ2 = 'REGL' DZ1 DZ2 6 ;
SZ3 = 'REGL' DZ2 DZ3 4 ;
SZ4 = 'REGL' DZ3 DZ4 6 ;
SZ5 = 'DROIT' (Rcuv 0.0 (Hpip - Rpip)) (Rcuv 0.0 (Hpip + Rpip)) 4 ;
Lsort = SZ5 'PLUS' (Lpip * OX) ;
SZ5 = 'REGL' SZ5 Lsort 3 ;

SZT = SZ1 ET SZ2 ET SZ3 ET SZ4 ET SZ5 ;
'ELIM' SZT EPSI ;

SS1 SS2 SS3 SS4 = SZ1 SZ2 SZ3 SZ4 'SYME' 'PLAN' OO OZ OY ;
SY0 = SZ1 ET SZ2 ET SZ3 ET SZ4 ET SZ5 ET SS1 ET SS2 ET SS3 ET SS4 ;
'ELIM' SY0 EPSI ;

SY1  = SY0 'PLUS' (Rinj * OY) ;
SY1P = ('DIFF' SY0 (SZ1 ET SS1)) 'PLUS' (Rinj * OY) ;
SY2  = ('DIFF' SY0 (SZ1 ET SS1)) 'PLUS' (Rpip * OY) ;
SY2P = ('DIFF' SY0 (SZ1 ET SS1 ET SZ5)) 'PLUS' (Rpip * OY) ;
SY3  = ('DIFF' SY0 (SZ1 ET SS1 ET SZ5)) 'PLUS' (Rcuv * OY) ;

'OPTI' 'ELEM' 'CUB8' ;

VY1 = 'VOLU' SY0  SY1 1 ;
VY2 = 'VOLU' SY1P SY2 1 ;
VY3 = 'VOLU' SY2P SY3 4 ;

VV1 VV2 VV3 = VY1 VY2 VY3 'SYME' 'PLAN' OO OX OZ ;

VVT = VY1 ET VY2 ET VY3 ET VV1 ET VV2 ET VV3 ;
'ELIM' VVT EPSI ;

S_INLET = ('ENVE' VVT) 'POINT' 'PROC' (0.0 0.0 (-1.0 * Linj)) ;
S_INLET = 'ELEM' ('ENVE' VVT) 'APPUYE' 'LARGEMENT' S_INLET ;
S_INLET = S_INLET 'COUL' 'ROUG' ;

S_OUTLET = ('ENVE' VVT) 'POINT' 'PROC' ((Rcuv + Lpip) 0.0 Hpip) ;
S_OUTLET = 'ELEM' ('ENVE' VVT) 'APPUYE' 'LARGEMENT' S_OUTLET ;
S_OUTLET = 'ELEM' ('ENVE' VVT) 'APPUYE' 'LARGEMENT' S_OUTLET ;
S_OUTLET = S_OUTLET 'COUL' 'BLEU' ;

PL_TOT = 'INTE' ('ENVE' (VY1 ET VY2 ET VY3))
                ('ENVE' (VV1 ET VV2 ET VV3)) ;
*'TRACE' 'CACH' (VVT ET S_INLET ET S_OUTLET) ;

*** ------------------------- ***
***  remplissage de la table  ***
*** ------------------------- ***
 PDA = 'TABLE' ;

 PDA.'FPAROI' = VRAI ;
 PDA.'DISCR'='QUAF'    ;
 PDA.'KPRES'='CENTREP1' ;
 PDA.'GRAPH'=VRAI ;
 PDA.'EPSI' =EPSI ;

*--> Definition de l'injection
 PDA.'scenario' = 'TABLE' ;
 PDA.'scenario'.'TIME' = 'PROG'   0.0   1.E4  ;
 PDA.'scenario'.'Qinj' = 'PROG'   0.01  0.01  ;

*--> Viscosité du mélange
 PDA.'NU' = 1.8e-5 ;

* utilisé dans le cas glissement + fonctions de paroi
*YP=1.E-2;

* definition d'un modele de turbulence
 PDA.'PARAM_TURB' = (Rinj / 2.0) ;

*--> definition des objets geometriques
* maillage total fluide
 PDA.'VTF' = VVT ;
* maillage de l'entrée
 PDA.'INLET'  = S_INLET   ;
* direction de l'injection à l'entrée
 PDA.'Dir_in'  = OZ ;
* maillage de la sortie
 PDA.'OUTLET' = S_OUTLET ;
* direction d'ejection de gaz à la sortie
 PDA.'Dir_out' = OX ;


 PDA.'DT0' = 0.05 ;
EXECPDA PDA 1 ;

'OPTI' 'ECHO' 0 ;

'SI' FAUX ;

MPL1 = 'CHAN' PL_TOT 'QUAF' ;
'ELIM' (MPL1 ET PDA.'GEO'.'M_VTF') EPSI ;
$PL1 = 'MODE' MPL1 'NAVIER_STOKES' PDA.'DISCR' ;
PL1 = 'DOMA' $PL1 'MAILLAGE' ;
UX1 = EXCO PDA.'TIC'.'UN' 'UX' ;
UY1 = EXCO PDA.'TIC'.'UN' 'UY' ;
UZ1 = EXCO PDA.'TIC'.'UN' 'UZ' ;
'MESS' (MAXI UX1) (MAXI UY1) (MAXI UZ1) ;
TRACE UX1 PL1 ;
TRACE UY1 PL1 ;
TRACE UZ1 PL1 ;
UN = PDA.'TIC'.'UN' ;
VN = 'VECT' ('REDU' UN PL1) 1.0 'UX' 'UY' 'UZ' 'ROUG' ;
TRACE VN PL1 ;
NUN = (UN (mots 'UX' 'UY') 'PSCA' UN (mots 'UX' 'UY')) ** 0.5 ;
TRACE NUN PL1 ;

'FINSI' ;

UN      = PDA.'TIC'.'UN'      ;
INLET   = PDA.'GEO'.'INLET'   ;
$INLET  = PDA.'GEO'.'$INLET'  ;
OUTLET  = PDA.'GEO'.'OUTLET'  ;
$OUTLET = PDA.'GEO'.'$OUTLET' ;
Qin  = 'DBIT' ('REDU' UN  INLET) $INLET  ;
Qout = 'DBIT' ('REDU' UN OUTLET) $OUTLET ;
'MESS' 'Debits E/S : ' Qin Qout ;

QStot = Qin + Qout ;

ERROR = 0 ;
'SI' (('ABS' QStot) > 1.E-10) ;
ERROR = ERROR + 1 ;
'FINSI' ;

'SI' (ERROR > 0) ;
'MESS' 'Problème de conservation du débit volumique sur l ' ;
'MESS' 'enveloppe fluide !!' ;
'ERREUR' 5 ;
'SINON' ;
'MESS' 'Tout s est bien passé ' ;
'FINSI' ;

'FIN' ;