* fichier :  cav_ray_proj.dgibi
************************************************************************
* Section : Fluides Thermique
************************************************************************

*----------------------------------------------
*     Couplage Rayonnement/Convection/Conduction
*
*     Test d'un cavité rayonnante en 2D plan
*       Le mélange gazeux est transparent
*               Rayleigh = 1e6
*     Résolution par la méthode de Projection
*     On introduit une paroi pour le rayonnement
*     afin de récupérer une constante de temps car
*     l'algorithme de projection est un algorithme
*     transitoire (DFDT)
*     Auteur : E. Studer - Mars 2002
*     Référence : 'Numerical Analysis of Natural convection
*                  with surface radiation in a square enclosure'
*                 M. Akiyama and Q.P. Chong
*                 Numerical Heat Transfer 31:419-433,1997
*----------------------------------------------



*----------------------------------------------
* Declarations générales
*----------------------------------------------
'OPTI' 'DIME' 2 'ELEM' 'QUA8' ;
'OPTI' 'TRAC' 'PSC' ;

COMPLET = FAUX ;
GRAPH = FAUX ;

****************
* PARAMETRES : *
****************

* nombre d'itérations internes
  iter = 1 ;
  omeg = 1.0 ;

* nombre d'itérations en temps

'SI' COMPLET ;
  itma = 400 ;
  d1 = 0.01 ;
  d2 = 0.05 ;
* elements en vitesse
   DISCR = 'QUAF' ;
'SINON' ;
  itma = 20 ;
  d1 = 0.01 ;
  d2 = 0.1 ;
* elements en vitesse
   DISCR = 'LINE' ;
'FINSI' ;

* pas de temps
  DT = 1.0 ;

* option de décentrement
*  KSUPG = 'CENTREE' ;
  KSUPG = 'SUPG' ;

* support de la pression
  KPRES = 'MSOMMET' ;
*  KPRES = 'CENTRE' ;
*  KPRES = 'CENTREP1' ;


* choix de la méthode de résolution IMPLICITE ou PROJECTION
*  RESOL = MOT IMPLICITE ;
  RESOL = MOT PROJECTION ;
  
* prise en compte du rayonnement
  PRAYO = VRAI ;
*  PRAYO = FAUX ;

* Ajout d'une paroi pour conduction et stabilisation
  PPAROI = VRAI ;
*  PPAROI = FAUX ;

* Parametres du maillage:

* Dimension de la cavité
  L = 1 ;

* Parametres physiques:
* EPSI = écart de température pour paroi chaude/froide
  EPSI = 0.017 ;

* Température de référence
  T0 = 293.5 ;

* T1 = Température chaude
* T2 = Température froide
  T1 = T0 '*' (1.0 '+' EPSI) ;
  T2 = T0 '*' (1.0 '-' EPSI) ;

*
  beta = 1.0 '/' T0 ;
  gbeta = -9.81 '*' beta ;
  gb = 0. gbeta ; 

  V0 = ( 'ABS' ( (gbeta '*' (T1 '-' T2)) '*' L ) ) '**' 0.5 ;
  Pr = 0.71 ;
  Ra = 1.E6 ;
  GR = Ra '/' Pr ;
  Nu = V0 '*' L '/'(Gr '**' 0.5) ;
  alpha = Nu '/' Pr ;
* conductivité de la paroi
  alphap = 1000.0 ;

  uref = 1 ;

* Rayonnement: emissivité des surfaces

'SI' (PRAYO) ;
  e_sud = 0.5 ;
  e_nord = 0.5 ;
  e_est = 0.5 ;
  e_ouest = 0.5 ;
'FINSI' ;



**************
* MAILLAGE : *
**************

p0 = 0. 0.   ;
p1 = 0.5 0.  ;
p2 = 1. 0.   ;
p3 = 1. 0.5  ;
p4 = 1. 1.   ;
p5 = 0.5 1.  ;
p6 = 0. 1.  ;
p7 = 0. 0.5  ;

d01 = 'DROI' p0 p1 'DINI' d1 'DFIN' d2 ;
d12 = 'DROI' p1 p2 'DINI' d2 'DFIN' d1 ;
d23 = 'DROI' p2 p3 'DINI' d1 'DFIN' d2 ;
d34 = 'DROI' p3 p4 'DINI' d2 'DFIN' d1 ;
d45 = 'DROI' p4 p5 'DINI' d1 'DFIN' d2 ;
d56 = 'DROI' p5 p6 'DINI' d2 'DFIN' d1 ;
d67 = 'DROI' p6 p7 'DINI' d1 'DFIN' d2 ;
d70 = 'DROI' p7 p0 'DINI' d2 'DFIN' d1 ;

sud   = d01 'ET' d12 ;
est   = d23 'ET' d34 ;
nord  = d45 'ET' d56 ;
ouest = d67 'ET' d70 ;

parois = (est 'ET' ouest 'ET' nord 'ET' sud) ;

'SI' (PPAROI) ;
* ajout d'une paroi pour mettre de la conduction
     mest = est 'TRAN' 1 (0.1 0.0) ;
     mest = mest 'COUL' 'BLEU' ;
     mouest = ouest 'TRAN' 1 (-0.1 0.0) ;
     mouest = mouest 'COUL' 'ROUGE' ;
     pest = 'COTE' 3 mest ;
     pouest = 'COTE' 3 mouest ;
     msud = sud 'TRAN' 1 (0.0 -0.1) ;
     msud = msud 'COUL' 'BLEU' ;
     mnord = nord 'TRAN' 1 (0. 0.1) ;
     mnord = mnord 'COUL' 'ROUGE' ;
'SI' (PRAYO) ;
     mparoi = mest 'ET' mouest 'ET' mnord 'ET' msud ;
'SINON' ;
     mparoi = mest 'ET' mouest ;
'FINSI' ;
'FINSI' ;

* Maillage du carre fluide
carre = sud est nord ouest 'DALLER' ;

* Maillage total
'SI' (PPAROI) ;
     mt = carre 'ET' mparoi ;
'SINON' ;
     mt = carre ;
'FINSI' ;

* Tracé des maillages
'SI' GRAPH ;
    'TRAC' carre 'TITR' 'MAILLAGE du Fluide' ;
    'TRAC' mt 'TITR' 'MAILLAGE total Fluide+paroi' ;
'FINSI' ;

**************************
* MODELE NAVIER-STOKES : *
**************************

carre = 'CHAN' carre 'QUAF' ;
mt = 'CHAN' mt 'QUAF' ;

'SI' (PPAROI) ;
     mparoi = 'CHAN' mparoi 'QUAF' ;
    'ELIM' (mt 'ET' carre 'ET' mparoi) 0.0001 ;
'SINON' ;
    'ELIM' (mt 'ET' carre) 0.0001 ;
'FINSI' ;


$mt = 'MODE' mt 'NAVIER_STOKES' DISCR ;
$carre = 'MODE' carre 'NAVIER_STOKES' DISCR ;

'SI' (PPAROI) ;
$mparoi = 'MODE' mparoi 'NAVIER_STOKES' DISCR ;
'FINSI' ;

$parois = 'MODE' parois 'NAVIER_STOKES' DISCR ;

'SI' (PPAROI) ;
$pest = 'MODE' pest 'NAVIER_STOKES' DISCR ;
$pouest = 'MODE' pouest 'NAVIER_STOKES' DISCR ;
'FINSI' ;

$sud = 'MODE' sud 'NAVIER_STOKES' DISCR ;
$est = 'MODE' est 'NAVIER_STOKES' DISCR ;
$nord = 'MODE' nord 'NAVIER_STOKES' DISCR ;
$ouest = 'MODE' ouest 'NAVIER_STOKES' DISCR ;

mtf = 'DOMA' $mt 'FACE' ;
carref = 'DOMA' $mt 'FACE' ;
parois = 'DOMA' $parois 'CENTRE' ;

'SI' (PPAROI) ;
pest = 'DOMA' $pest 'CENTRE' ;
pouest = 'DOMA' $pouest 'CENTRE' ;
'FINSI' ;

est = 'DOMA' $est 'CENTRE' ;
sud = 'DOMA' $sud 'CENTRE' ;
ouest = 'DOMA' $ouest 'CENTRE' ;
nord = 'DOMA' $nord 'CENTRE' ;

'SI' (PPAROI) ;
'ELIM' carref (parois 'ET' pest 'ET' pouest
'ET' est 'ET' sud 'ET' ouest 'ET' nord) 1.D-4 ;
'ELIM' mtf (parois 'ET' pest 'ET' pouest
'ET' est 'ET' sud 'ET' ouest 'ET' nord) 1.D-4 ;
'SINON' ;
'ELIM' carref (parois 'ET' est 'ET' sud 'ET' ouest 'ET' nord) 1.D-4 ;
'ELIM' mtf (parois 'ET' est 'ET' sud 'ET' ouest 'ET' nord) 1.D-4 ;
'FINSI' ;

parois = 'DOMA' $parois 'MAILLAGE' ;

'SI' (PPAROI) ;
pest = 'DOMA' $pest 'MAILLAGE' ;
pouest = 'DOMA' $pouest 'MAILLAGE' ;
pest = pest 'COUL' 'BLEU' ;
pouest = pouest 'COUL' 'ROSE' ;
'FINSI' ;

est = 'DOMA' $est 'MAILLAGE' ;
ouest = 'DOMA' $ouest 'MAILLAGE' ;
sud = 'DOMA' $sud 'MAILLAGE' ;
nord = 'DOMA' $nord 'MAILLAGE' ;
mt = 'DOMA' $mt 'MAILLAGE' ;
carre = 'DOMA' $carre 'MAILLAGE' ;

mt = mt 'COUL' 'VERT' ;
parois = parois 'COUL' 'ROUGE' ;

*
*****************
* MODELISATION : *
******************
'SI' (PRAYO) ;
* Création des modèles de rayonnement
mrsud = MODE sud THERMIQUE 'RAYONNEMENT'  'CAVITE' CONS 'CAV1';
mrest = MODE est  THERMIQUE  'RAYONNEMENT'  'CAVITE'  CONS 'CAV1';
mrnord = MODE nord  THERMIQUE 'RAYONNEMENT' 'CAVITE'  CONS 'CAV1';
mrouest = MODE ouest  THERMIQUE 'RAYONNEMENT' 'CAVITE'  CONS 'CAV1';

* assemblage des modèles de rayonnement
mray = mrsud 'ET' mrest 'ET' mrnord 'ET' mrouest ;

        
* création des matériaux émissivité
esud = 'MATE' mrsud 'EMIS' e_sud ;
eest = 'MATE' mrest 'EMIS' e_est ;
enord = 'MATE' mrnord 'EMIS' e_nord ;
eouest = 'MATE' mrouest 'EMIS' e_ouest ;

* assemblage des matériaux
e = esud 'ET' eest 'ET' enord 'ET' eouest ;

*************************************************
* FACTEURS DE FORME et MATRICE DE RAYONNEMENT : *
*************************************************

fft = 'FFOR' mray  e;
chamr = 'RAYE' mray fft e ;

**********************************************************
* CREATION DE LA TABLE POUR L'OPERATEUR DE RAYONNEMENT : *
**********************************************************

TAB_RAY = 'TABLE' ;
TAB_RAY.'MAILLAGE' = (sud 'ET' est 'ET' nord 'ET' ouest) ;
TAB_RAY.'MODELE' = mray ;
TAB_RAY.'MATERIAU' = e ;
TAB_RAY.'FFORM' = fft ;
TAB_RAY.'MATRICE' = chamr ;
* Modèles Navier-Stokes
TAB_RAY.'MODELENS' = ($sud 'ET' $est 'ET' $nord 'ET' $ouest) ;
'FINSI' ;


******************************
* OPERATEUR DE RAYONNEMENT : *
******************************

'DEBP' FRAY ;

'ARGU' RVX*'TABLE' ; 

RV = RVX.'EQEX' ;
* récupération du champ de température
T = RV.INCO.'TN' ;
* récupération de la table TABR associée à l'opérateur
TABR = RVX.'ARG1' ;
MAILR = TABR.'MAILLAGE' ;
* TSU : température de surface
TSU = 'REDU' T MAILR ;
TSU = 'EXCO' TSU 'SCAL' 'T' ;
TR = 'CHAN' 'CHAM' TSU TABR.'MODELE' 'GRAVITE' ;
* MRS : matrice de rayonnement
MRS = 'RAYN' TABR.'MODELE' TABR.'MATRICE' TR ;
MAILR = 'DOMA' TABR.'MODELENS' 'SOMMET' ;
* MRS : transformation Rigidité/Matrik
MRS = 'KOPS' 'RIMA' MRS MAILR ('MOTS' 'TN') ;

AS AF = 'KOPS' 'MATRIK' ;

'FINP' AS MRS ;

*************************
* CREATION DES TABLES : *
*************************

corx = 'COOR' 1 carre ;

RV = 'EQEX' $carre 'NITER' iter 'OMEGA' omeg 'ITMA' 1 'ALFA' 1.0
   'FIDT' 100
   'OPTI' 'EF' KSUPG 'IMPL' KPRES
   'ZONE' $carre 'OPER' 'NS' 1. 'UN' Nu gb 'TF' T0 'INCO' 'UN'
;

RT = 'EQEX' $mt 'NITER' iter 'OMEGA' 1.0 'ITMA' 1 'ALFA' 1.0
   'OPTI' 'EF' KSUPG 'IMPL' KPRES
   'ZONE' $carre 'OPER' 'TSCAL' 1. 'UN' alpha 0. 'INCO' 'TN'
;
'SI' (PPAROI) ;
RT = 'EQEX' RT
   'OPTI' 'EF' KSUPG 'IMPL' KPRES
   'ZONE' $mparoi   'OPER' 'LAPN' alphap 'INCO' 'TN'
;
'FINSI' ;

'SI' ('EGA' RESOL 'PROJECTION') ;
RV = 'EQEX' RV
   'OPTI' 'EF' 'CENTREE'
   'ZONE' $carre 'OPER' 'DFDT' 1.0 'UNM' DT 'INCO' 'UN'
;
RT = 'EQEX' RT
   'OPTI' 'EF' 'CENTREE'
   'ZONE' $mt    'OPER' 'DFDT' 1.0 'TNM' DT 'INCO' 'TN'
;
'FINSI' ;

SI ('EGA' RESOL 'IMPLICITE') ;
RV = 'EQEX' RV
   'OPTI' 'EFM1' 'CENTREE'
   'ZONE' $carre 'OPER' 'DFDT' 1.0 'UNM' DT 'INCO' 'UN'
;
RT = 'EQEX' RT
   'OPTI' 'EFM1' 'CENTREE'
   'ZONE' $mt    'OPER' 'DFDT' 1.0 'TNM' DT 'INCO' 'TN'
;
'FINSI' ;

'SI' (PRAYO) ;
RT = 'EQEX' RT
   'OPTI' 'EF' 'CENTREE'
   'ZONE' ($sud 'ET' $est 'ET' $nord 'ET' $ouest)
   'OPER' 'FRAY' TAB_RAY 'INCO' 'TN'
;
'FINSI' ;

RV = 'EQEX' RV
   'CLIM' 'UN' 'UIMP' parois 0.
         'UN' 'VIMP' parois 0.
;
'SI' (PPAROI) ;
RT = 'EQEX' RT
   'CLIM' 'TN' 'TIMP' pouest T1
         'TN' 'TIMP' pest   T2 ;
'SINON' ;
RT = 'EQEX' RT
   'CLIM' 'TN' 'TIMP' ouest T1
         'TN' 'TIMP' est   T2 ;
'FINSI' ;

RV.'INCO' = 'TABLE' 'INCO' ;
RV.'INCO'.'UN' = 'KCHT' $carre 'VECT' 'SOMMET' (0. 0.) ;
RV.'INCO'.'TN' = 'KCHT' $mt 'SCAL' 'SOMMET' T0 ;
RV.'INCO'.'TF' = 'KCHT' $carre 'SCAL' 'SOMMET' T0 ;
RV.'INCO'.'UNM' = 'KCHT' $carre 'VECT' 'SOMMET' (0. 0.) ;
RV.'INCO'.'TNM' = 'KCHT' $mt 'SCAL' 'SOMMET' T0 ;
RV.'INCO'.'PRES' = 'KCHT' $carre 'SCAL' KPRES (0.0);
RV.'INCO'.'UNM1' = 'KCHT' $carre 'VECT' 'SOMMET' (1.e-5 1.e-5) ;
RV.'INCO'.'TNM1' = 'KCHT' $mt 'SCAL' 'SOMMET' T0 ;
RV.'INCO'.'IT' = 'PROG' 1 ;
RV.'INCO'.'NITER' = 0 ;
RV.'INCO'.'ERX' = 'PROG' 0. ;
RV.'INCO'.'ERY' = 'PROG' 0. ;
RV.'INCO'.'ERT' = 'PROG' 0. ;

RT.'INCO' = RV.'INCO' ;

point = ('DOMA' $carre KPRES) 'ELEM' 1 ;

'SI' ('EGA' RESOL 'PROJECTION') ;
 RVP = 'EQEX' 
    'OPTI' 'EF' KPRES 
    'ZONE' $carre 'OPER' KBBT (-1.0) 'INCO' 'UN' 'PRES'
    'CLIM' 'PRES' 'TIMP' point 0. ;
'FINSI' ;

'SI' ('EGA' RESOL 'IMPLICITE') ;
 RV = 'EQEX' RV 
    'OPTI' 'EF' 'IMPL' KPRES KSUPG
    'ZONE' $carre 'OPER' KBBT (1.0) 'INCO' 'UN' 'PRES'
    'CLIM' 'PRES' 'TIMP' point 0. ;
'FINSI' ;

'SI' ('EGA' RESOL 'PROJECTION') ;
  RV.'PROJ' = RVP ;
'FINSI' ;


rv.'METHINV'.TYPINV=1 ;
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   ;

rt.'METHINV'.TYPINV=1 ;
rt.'METHINV'.IMPINV=0 ;
rt.'METHINV'.NITMAX=400;
rt.'METHINV'.PRECOND=3 ;
rt.'METHINV'.RESID  =1.e-10;
rt. 'METHINV' . 'FCPRECT'=1   ;
rt. 'METHINV' . 'FCPRECI'=1   ;

'SI' ('EGA' RESOL 'PROJECTION') ;
rvp.'METHINV'.TYPINV=1 ;
rvp.'METHINV'.IMPINV=0 ;
rvp.'METHINV'.NITMAX=400;
rvp.'METHINV'.PRECOND=3 ;
rvp.'METHINV'.RESID  =1.e-8 ;
rvp.'METHINV' . 'FCPRECT'=100 ;
rvp.'METHINV' . 'FCPRECI'=100 ;
'FINSI' ;

****************
* RESOLUTION : *
****************

'REPETER' BOU1 ITMA ;
* projection de la température fluide
TN = RV.'INCO'.'TN' ;
RV.'INCO'.'TF' = 'KCHT' $carre 'SCAL' 'SOMMET' 
    ('REDU' TN ('DOMA' $CARRE 'MAILLAGE')) ;
RV.INCO.'NITER' = RV.INCO.'NITER' + 1 ;
DD = RV.INCO.'NITER' ;
NN = DD '/' 5;
LO = (5 '*' NN '-' DD) 'EGA' 0 ;

'SI' (LO) ;

 UN = RV.INCO.'UN' ;
 UNM1 = RV.INCO.'UNM1' ;
 unx = 'KCHT' $carre 'SCAL' 'SOMMET' (EXCO 'UX' UN) ;
 unm1x = 'KCHT' $carre 'SCAL' 'SOMMET' (EXCO 'UX' UNM1) ;
 uny = 'KCHT' $carre 'SCAL' 'SOMMET' (EXCO 'UY' UN) ;
 unm1y = 'KCHT' $carre 'SCAL' 'SOMMET' (EXCO 'UY' UNM1) ;
 TN = RV.INCO.'TN' ;
 TNM1 = RV.INCO.'TNM1' ;

 ERRX = 'KOPS' unx '-' unm1x ;
 ELIX = 'MAXI' ERRX 'ABS' ;
 ELIX = ('LOG' (ELIX '+' 1.e-10)) '/' ('LOG' 10.) ;
 ERRY = 'KOPS' uny '-' unm1y ;
 ELIY = 'MAXI' ERRY 'ABS' ;
 ELIY = ('LOG' (ELIY '+' 1.e-10)) '/' ('LOG' 10.) ;
 ERRT = 'KOPS' TN '-' TNM1 ;
 ELIT = 'MAXI' ERRT 'ABS' ;
 ELIT = ('LOG' (ELIT '+' 1.e-10)) '/' (LOG 10.) ;
 
 'MESS' 'ITER' DD  'ERREUR LINF TEMPERATURE' ELIT ;
 'MESS' 'ITER' DD  'ERREUR LINF VITESSE' ELIX ELIY ;

 RV.INCO.'TNM1' = 'KCHT' $mt 'SCAL' 'SOMMET' (RV.INCO.'TN') ;
 RV.INCO.'UNM1' = 'KCHT' $carre 'VECT' 'SOMMET' (RV.INCO.'UN') ;
 IT = 'PROG' DD ;
 ERX = 'PROG' ELIX ;
 ERY = 'PROG' ELIY ;
 ERT = 'PROG' ELIT ;
 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.'ERT' = (RV.INCO.'ERT') 'ET' ERT ;

'FINSI' ;
     EXEC RV ;
     EXEC RT ;
'FIN' BOU1 ;

******************************
* EXPLOITATION DES RESULTATS *
******************************

'OPTI' 'ISOV' 'SULI';

'SI' GRAPH ;
    'TRAC' carre 'TITR' 'MAILLAGE' ;
'FINSI' ;

* Adimensionnement des températures pour dessin
LTF = 'PROG' -0.5 -0.4 -0.3 -0.2 -0.1 0.0 0.1
  0.2 0.3 0.4 0.5 ;
TN = ((RV.INCO.'TN') '-' T0 ) '/' (T1 '-' T2) ;
'SI' GRAPH ;
    'TRAC' TN carre (parois) 'TITR' 'TEMPERATURE' LTF ;
'FINSI' ;

'SI' GRAPH ;
    vit = 'VECT' RV.INCO.'UN' (1. '/' V0) UX UY 'TURQ' ;
    'TRAC' vit carre (parois) 'TITR' 'VITESSE';
'FINSI' ;

* Fonction de courant
un = RV.INCO.'UN' ;
un2 = 'KOPS' un '*' (Pr '*' (Gr '**' 0.5)) ;
sw = 'KOPS' un 'ROT' $carre;
rk = 'EQEX' $carre 'OPTI' 'EF' 'IMPL'
     'ZONE' $carre 'OPER' 'LAPN' -1. INCO 'PSI'
     'ZONE' $carre 'OPER' 'FIMP' sw INCO 'PSI'
     'CLIM' 'PSI' 'TIMP' parois 0.
;

rk.'INCO' = 'TABLE' 'INCO' ;
rk.'INCO'.'PSI' = 'KCHT' $carre 'SCAL' 'SOMMET' 0. ;
EXEC rk ;
psi = rk.INCO.'PSI' '/' V0 ;
psi2 = 'KOPS' psi '*' (Pr '*' (Gr '**' 0.5)) ;
'SI' GRAPH ;
    'TRAC' psi2 carre parois 14 'TITR' 'FONCTION DE COURANT';
'FINSI' ;


EVOL3 = 'EVOL' 'MANU' 'ITERATIONS' (RV.INCO.'IT') 'LOG|E|inf'
     (RV.INCO.'ERT') 'COUL' VERT ;
EVOL4 = 'EVOL' 'MANU' 'ITERATIONS' (RV.INCO.'IT') 'LOG|E|inf'
     (RV.INCO.'ERX') ;
EVOL5 = 'EVOL' 'MANU' 'ITERATIONS' (RV.INCO.'IT') 'LOG|E|inf'
     (RV.INCO.'ERY') 'COUL' 'ROUGE' ;
'SI' GRAPH ;
'DESS' (EVOL4 ET EVOL3) 'YBOR' -10.0 0.0 
        'TITR' 'CONVERGENCE VERS LE STATIONNAIRE' ;
'FINSI' ;

'SI' ('NON' COMPLET) ;
* Test de non régression
* Erreur de référence sur la température
ELITM =  -8.30554E-02 ; 
* Valeur de référence sur le maximum de la vitesse verticale
UYM = 'MAXI' ('EXCO' RV.'INCO'.'UN' 'UY') ;
UYMAX = 0.18220; 
* Valeur de référence sur le maximum de la fonction de courant
PSIM = 'MAXI' (PSI2) ;
PSIMAX = 27.652 ; 
*
ERT = 'ABS'( ELIT '-' ELITM) '/' ELITM ;
ERU = 'ABS'( UYM '-' UYMAX) '/' UYMAX ;
ERP = 'ABS'( PSIM '-' PSIMAX) '/' PSIMAX ;
'MESS' 'Erreur sur la température = ' ERT ;
'MESS' 'Erreur sur le max de UY   = ' ERU ;
'MESS' 'Erreur sur le max de Psi  = ' ERP ;

'SI' (ERT  '>' 0.002); 'ERREUR' 5 ; 'FINSI' ;
'SI' (ERU  '>' 0.002); 'ERREUR' 5 ; 'FINSI' ;
'SI' (ERP  '>' 0.002); 'ERREUR' 5 ; 'FINSI' ;

'FINSI' ;
'FIN' ;
*
*
*