Télécharger stationary_shock.dgibi

Retour à la liste

* fichier :  stationary_shock.dgibi
***********************************************************************
***********************************************************************
***********************************************************************
*** CALCUL DU TUBE A CHOC; CAS MULTIESPECE                         ****
*** GAZ mono-especes "THERMALLY PERFECT"                           ****
***                                                                ****
*** FORMULATION VF COMPRESSIBLE EXPLICITE                          ****
*** DIFFERENTS SOLVEURS                                            ****
***                                                                ****
*** Choc stationnaire                                              ****
*** A. BECCANTINI DRN/DMT/SEMT/LTMF    FEVRIER 2000                ****
***********************************************************************
***********************************************************************
***********************************************************************
 
 'MESSAGE' 'A mettre a jours' ;
 'FIN' ;
 
'OPTION'  'DIME' 2 ;
'OPTION'  'ELEM' 'QUA4' ;
'OPTION'  'ISOV' 'SULI' ;
'OPTION'  'ECHO' 0 ;
'OPTION' 'TRAC' 'X';
 
GRAPH   = FAUX ;
* GRAPH   = VRAI ;
 
*
*** Methodes possibles :
*
*   'VLH'
*   'CG'
 
METO = 'CG' ;
 
 
***********************
**** LA TABLE PGAZ ****
***********************
 
PGAZ = 'TABLE' ;
 
*
**** Ordre des polynoms cv_i
*
 
PGAZ . 'NORD' = 4 ;
 
*
**** La seul espece 
*
 
PGAZ . 'ESPNEULE' = 'H2  ';
 
*
 
PGAZ .  'H2  ' = 'TABLE'  ;
 
*
**** R (J/Kg/K)
*
 
PGAZ .  'H2  ' . 'R' = 4130.0 ;
 
*
**** Regressions polynomials
*
 
   PGAZ .  'H2  ' . 'A' = 'PROG'  9834.91866 0.54273926 0.000862203836
                               -2.37281455E-07 1.84701105E-11 ;
 
 
*
**** "Enthalpies" (ou energies) de formations a OK (J/Kg)
*     Note: ce sont des entites numeriques
*
 
PGAZ .  'H2  ' . 'H0K' = -4.195D6 ;
 
*
*** Fin PGAZ
*
 
*
**** Choc stationnaire
*
 
pg = 1.0D5 ;
rog = 1.0 ;
Tg = (pg '*' rog) '/' (PGAZ .  'H2  ' . 'R') ;
 
* sstren > 1.0
 
sstren = 4.01 ;
Td = sstren '*' Tg ;
 
* Det. ug, rod, ud
 
A0H2 = 'EXTRAIRE' (PGAZ . 'H2  ' . 'A') 1 ;
A1H2 = 'EXTRAIRE' (PGAZ . 'H2  ' . 'A') 2 ;
A2H2 = 'EXTRAIRE' (PGAZ . 'H2  ' . 'A') 3 ;
A3H2 = 'EXTRAIRE' (PGAZ . 'H2  ' . 'A') 4 ;
A4H2 = 'EXTRAIRE' (PGAZ . 'H2  ' . 'A') 5 ;
 
R = PGAZ . 'H2  ' . 'R' ;
 
T2 = Tg '*' Tg ;
T3 = T2 '*' Tg ;
T4 = T3 '*' Tg ;
T5 = T4 '*' Tg ;
 
ETHERg =  ((A0H2 * Tg) '+' ((A1H2 '/' 2.0) * T2) '+'
         ((A2H2 '/' 3.0) * T3) '+' ((A3H2 '/' 4.0) * T4) '+'
         ((A4H2 '/' 5.0) * T5)) ;
 
 
T2 = Td '*' Td ;
T3 = T2 '*' Td ;
T4 = T3 '*' Td ;
T5 = T4 '*' Td ;
 
ETHERd =  ((A0H2 * Td) '+' ((A1H2 '/' 2.0) * T2) '+'
         ((A2H2 '/' 3.0) * T3) '+' ((A3H2 '/' 4.0) * T4) '+'
         ((A4H2 '/' 5.0) * T5)) ;
 
 
RTg = R '*' Tg ;
RTd = R '*' Td ;
 
* x = rog '/' rod `in (0,1]
 
NDEF = 100 ;
DX = 1.0 /NDEF ;
X = 0.0 ;
AEQ2 = RTg ;
AEQ1 = (2.0 '*' (ETHERd '-' ETHERg)) '+' RTd '-' RTg ;
AEQ0 = -1.0 '*' RTd ;
 
LIX = 'PROG' X ;
LIF = 'PROG' AEQ0 ;
 
'REPETER' BL1 NDEF ;
    X = X '+' DX ;
    X2 = X '*' X ;
    FUN = (AEQ2 '*' X2) '+' (AEQ1 '*' X) '+' AEQ0 ;
    LIX = LIX 'ET' ('PROG' X) ;
    LIF = LIF 'ET' ('PROG' FUN) ;
'FIN' BL1 ;
 
EVFUN = 'EVOL' 'MANU' 'rog/rod' LIX 'fun' LIF ;
 
'SI' FAUX ;
   'DESSIN' EVFUN 'MIMA' ;
'FINSI' ;   
 
* Selection de la racine:
* XSOL \in (0,1]
 
XSOL = ((((AEQ1 '**' 2) '-' (4.0 '*' AEQ0 * AEQ2)) '**' 0.5) '-' AEQ1)
       '/' (2.0 * AEQ2) ;
 
*
 
rod = rog '/' XSOL ;
 
pd = RTd '*' rod;
 
RC1 = 2.0 '*' (ETHERd '+' (R '*' Td) '-' (ETHERg '+' (R '*' Tg))) ;
ud = RC1 '*' (XSOL '**' 2) '/' (1.0 '-' (XSOL '**' 2)) ;
ud = ud '**' 0.5 ;
 
ug = (ud '*' rod) '/' rog ;
 
 
* Control de Rankine-Hugoniot
 
err1 = 'ABS' (((rog '*' ug) '-' (rod '*' ud)) '/' (rog '*' ug)) ;
err2 = 'ABS' (((rog '*' ug '*' ug) '+' pg '-'
               (rod '*' ud '*' ud) '-' pd)
        '/' ((rog '*' ug '*' ug) '+' pg)) ;
err3 = 'ABS' (((ETHERg '+' RTg '+' (0.5 '*' ug '*' ug)) '-'
        (ETHERd '+' RTd '+' (0.5 '*' ud '*' ud))) '/' 
        (ETHERd '+' RTd '+' (0.5 '*' ud '*' ud))) ;
 
ERMAX = 'MAXIMUM' ( 'PROG' err1 err2 err3) ;        
 
'SI' (ERMAX > 1.0d-10);
    'ERREUR' 'No Rankine-Hugoniot';
'FINSI' ;
 
 
*******************************************************
***** PROCEDURE POUR TESTER CONVERGENCE           *****
*******************************************************
 
'DEBPROC'  CALCUL ;
'ARGUMENT'  RVX*'TABLE' ;
 
RV    = RVX . 'EQEX' ;
 
*KDOMA = RV . 'DOMAINE' ;
KDOMA = RV . 'MODTOT' ;
 
RNi   = RV . 'INCO' . 'RNI' ;
RN0   = RV . 'INCO' . 'RN0' ;
 
*
*** RN0 = solution à t = tn_M1;
*
 
DD = RV . 'PASDETPS' . 'NUPASDT' ;
NINTE = 5;
NN = DD/NINTE ;
LO = (DD-(NINTE*NN)) EGA 0 ;
 
'SI'  ( LO ) ;
 
   RNi = 'KCHT' KDOMA 'SCAL' 'CENTRE' RNi ;
   RN0 = 'KCHT' KDOMA 'SCAL' 'CENTRE' RN0 ;
 
   ERR = 'KOPS'  RNi '-' RN0 ;
 
   ERR = 'KOPS' ERR '/' RN0;
 
   ELI = 'MAXIMUM'  ERR 'ABS' ;
   'MESSAGE' eli ;
   ELI = ('LOG'  (ELI + 1.0E-10))/('LOG'  10.) ;
   'MESSAGE' ('CHAINE'
      'ITER  ' DD  '   ERREUR LINF ' ELI );
   IT = 'PROG' DD ;
   ER = 'PROG'  ELI ;
   RV . 'INCO' . 'IT' = (RV . 'INCO' . 'IT') 'ET'  IT ;
   RV . 'INCO' . 'ER' = (RV . 'INCO' . 'ER') 'ET'  ER ;
'FINSI'  ;
 
RV . 'INCO' . 'RN0' = 'COPIER' RNi;
 
'FINPROC'  ;
 
*****************************************************
*****************************************************
** PROCEDURE EXEX POUR FORMULATION VF COMPRESSIBLE **
** CAS MULTIESPECES "THERMALLY PERFECT"            **
*****************************************************
*****************************************************
 
'DEBPROC'  EXEX ;
'ARGUMENT' RV*TABLE ;
 
*******************************************
**** Recherche de RV . *KONV . IDEUL   ****
*******************************************
 
*
**** Nom de la table RV . *'KONV' -> NOMT
*
 
  NBOP = 'DIME'  (RV . 'LISTOPER' ) ;
 
  'REPETER'  BL1 NBOP;
      MCEL = 'EXTRAIRE'  &BL1 RV . 'LISTOPER';
      'SI' ( 'EGA' MCEL 'KONV    ');
       NOMT = 'MOT' ('TEXTE' ('CHAINE'  &BL1 MCEL));
       'QUITTER' BL1;
     'FINSI' ;
  'FIN' BL1;
 
  IEUL = RV . NOMT . 'KOPT' . 'IDEUL';
  'SI' ('NON' (IEUL 'EGA' 3));
      'MESSAGE' 'EULERMST ???';
      'ERREUR' 21;
  'FINSI' ;
 
* Mono-espece ou multi-especes  
 
  'SI' ('EXISTE' (RV . NOMT) 'ARG5') ;
      LOGMUL = VRAI ;
  'SINON' ;
      LOGMUL = FAUX ;
  'FINSI' ;
 
*
**** CL
*
 
 
  LOGLIM = RV . 'INCO' . 'CLIM' ;
 
 
******************************************
**** Ordre en espace, ordre en temps  ****
******************************************
 
ORD_ESP = RV . 'ORDREESP' ;
ORD_TPS = RV . 'ORDRETPS' ;
 
'MESSAGE'  '--------------------------';
'MESSAGE'  'Ordre en Espace :' ord_esp;
'MESSAGE'  'Ordre en Temps  :' ord_tps;
'MESSAGE'  '--------------------------';
 
'SI' ((ORD_ESP 'EGA' 1) 'ET' (ORD_TPS 'EGA' 2));
    'MESSAGE' ;
    'MESSAGE' (CHAINE 'Ordre en Espace doit etre 2');
    'MESSAGE' (CHAINE 'On impose ça.');
    'MESSAGE' ;
      RV . 'ORDREESP' = 2 ;
    'MESSAGE' ;
    'MESSAGE'  '--------------------------';
    'MESSAGE'  'Ordre en Espace :' ord_esp;
    'MESSAGE'  'Ordre en Temps  :' ord_tps;
    'MESSAGE'  '--------------------------';
'FINSI' ;
 
 
******************************
**** La table 'PASDETPS'  ****
******************************
 
 TPSI = RV . 'TPSI' ;
 TFIN = RV . 'TFINAL'; 
 RV . 'PASDETPS' . 'TPS' = TPSI;
 
*
**** DELTAT-1 est un argument de PRET (prediction)
*    Donc on doit l'initialiser.
*
 
 RV . 'PASDETPS' . 'DELTAT-1' = 0.0D0;
 CFL = rv.'ALFA' ;
 
 
*********************
**** Les TABLES *****
*********************
 
*   
**** RV . 'INCO'      
*    RV . 'DOMAINE'
*    RV . 'KIZD'
*    RV . 'KIZG'
 
 
*   
**** RV . 'INCO'      -> KINCO
*
 
KINCO   = (RV . 'INCO') ;
 
 
*
**** RV . 'DOMAINE'   -> KDOMA
*
 
*KDOMA  = (RV . 'DOMAINE') ;
KDOMA  = (RV . 'MODTOT') ;
KDOMA2  = (RV . 'DOMAINE') ;
 
*
**** RV . 'KIZD' contient les volumes des ELTs
*
 
'SI'  ('NON' ('EXISTE'   RV  'KIZD')) ;
    'KDIA'  RV ;
'FINSI'  ;
 
*
***** RV . 'KIZG' contient les flux aux interfaces. 
*
 
'SI'  ('NON' ('EXISTE'  RV  'KIZG')) ;
    RV . 'KIZG' = 'TABLE'   'KIZG' ;
'FINSI'  ;
 
'SI' LOGMUL ;
 
*********************************************************
****   Multi-especes, boucle Sur les Pas de Temps    ****
*********************************************************
 
*
**** Evaluation de coeff pour le calcule des pentes
*
 
 KINCO . 'V'  KINCO . 'P'  KINCO . 'T'  KINCO . 'Y'
      KINCO . 'GAMMA' =   'PRIM' 'PERFTEMP' (KINCO . 'IPGAZ')
      (KINCO . 'RNI') (KINCO . 'GNI') (KINCO . 'ENI')
      (KINCO. 'RYNI');
 
 GRADR ALR COEFR = 'PENT' KDOMA 'CENTRE' 'EULESCAL'   'LIMITEUR'
                     (KINCO . 'RNI');
 
 GRADP ALP COEFP = 'PENT' KDOMA 'CENTRE' 'EULESCAL'   'LIMITEUR'
                     (KINCO . 'P');
 
 GRADV ALV COEFV = 'PENT' KDOMA 'CENTRE' 'EULEVECT'   'LIMITEUR'
                     (KINCO . 'V');
 
 GRADY ALY COEFY = 'PENT' KDOMA 'CENTRE' 'EULESCAL'   'LIMITEUR'
                     (KINCO . 'Y');
 
 
I = 0 ;
 
'REPETER'  BLOC1 (RV . 'ITMA')  ;
 
    I = I + 1 ;
 
*
***** Les variables primitives
*
         KINCO . 'V'  KINCO . 'P'  KINCO . 'T'  KINCO . 'Y'
         KINCO . 'GAMMA' =   'PRIM' 'PERFTEMP' (KINCO . 'IPGAZ')
         (KINCO . 'RNI') (KINCO . 'GNI') (KINCO . 'ENI')
         (KINCO. 'RYNI');
 
         'SI'  (ORD_ESP 'EGA'  1) ;
 
             ROF VITF PF YF  =  'PRET' 'PERFTEMP'
                 ORD_ESP ORD_TPS KDOMA (KINCO . 'IPGAZ') (KINCO . 'RNI')  
                 (KINCO . 'V') (KINCO . 'P')  (KINCO . 'Y') ;
 
         'SINON';
 
*       
***** Ordre 2 en espace => calcul des pentes
*
 
            GRADR ALR  = 'PENT' KDOMA 'CENTRE' 'EULESCAL'   'LIMITEUR'
                     (KINCO . 'RNI')  'GRADGEO' COEFR ;
 
            GRADP ALP  = 'PENT' KDOMA 'CENTRE' 'EULESCAL'   'LIMITEUR'
                     (KINCO . 'P')  'GRADGEO' COEFP ;
 
            GRADV ALV  = 'PENT' KDOMA 'CENTRE' 'EULEVECT'   'LIMITEUR'
                     (KINCO . 'V')  'GRADGEO' COEFV ;
 
            GRADY ALY  = 'PENT' KDOMA 'CENTRE' 'EULESCAL'   'LIMITEUR'
                     (KINCO . 'Y')  'GRADGEO' COEFY ;
 
           'SI' (ORD_TPS 'EGA' 1);
 
               ROF VITF PF YF  = 'PRET' 'PERFTEMP'  ORD_ESP ORD_TPS
                                   KDOMA (KINCO . 'IPGAZ')
                                   (KINCO . 'RNI') GRADR ALR
                                   (KINCO . 'V')   GRADV ALV
                                   (KINCO . 'P')   GRADP ALP
                                   (KINCO . 'Y')   GRADY ALY ;
 
 
           'SINON' ;              
*       
********* Ordre 2 en temps
*
               ROF VITF PF YF  = 'PRET' 'PERFTEMP'  ORD_ESP ORD_TPS
                                  KDOMA (KINCO . 'IPGAZ')                                
                                  (KINCO . 'RNI') GRADR ALR
                                  (KINCO . 'V')   GRADV ALV
                                  (KINCO . 'P')   GRADP ALP
                                  (KINCO . 'Y')   GRADY ALY
                                   (KINCO . 'GAMMA')                                   
                                  ((RV . 'PASDETPS' . 'DELTAT-1')/2.0);
           'FINSI' ;
 
        'FINSI'  ;
 
*
*********** Creation de MCHAML de type 'FACEL' pour les
*           calcul de flux aux interfaces
 
 
 
 KINCO . 'RNF'    =  ROF ;
 KINCO . 'VITNF'  =  VITF ;
 KINCO . 'PNF'    =  PF ;
 KINCO . 'YF'     =  YF ;
 
 
*
********* Boucle sur les operateurs
*
 
   'REPETER'   BLOC2  NBOP ;
        NOMPER  = 'EXTRAIRE'  &BLOC2  (RV . 'LISTOPER');
        NOTABLE = 'MOT'  ('TEXTE'  ('CHAINE'  &BLOC2  NOMPER) ) ;
        ('TEXTE'  NOMPER) (RV . NOTABLE) ;
   'FIN'  BLOC2 ;
 
*
********* Mise a jour de la table RV . 'PASDETPS' 
*
 
   'SI' ('EXISTE'  RV  'DTI');
        DTI = 'MINIMUM'
            ('PROG' ((RV . 'DTI') '/' CFL)
            (RV . 'PASDETPS' . 'DTCONV') );
   'SINON';
        DTI = (RV . 'PASDETPS' . 'DTCONV');
   'FINSI';
 
    RV .  'PASDETPS' . 'DELTAT'  = DTI ;
 
    TMPS = RV . 'PASDETPS' . 'TPS';
    DTI0 = TFIN '-' TMPS;
    DTI0 = DTI0 '/' CFL;
 
    'SI' (DTI0 '<EG' DTI);
        DTI = DTI0;
        RV . 'PASDETPS' . 'DELTAT' = DTI ;
        LOGQUIT = VRAI;
    'SINON' ;
        LOGQUIT = FAUX ;
    'FINSI' ;
 
 
*
******** On avance au pas de temps suivant 
*
*   N.B. Tn+1 = Tn + (CFL *  RV . 'PASDETPS' . 'DELTAT');
*
 
    'AVCT'  RV CFL  'IMPR' (RV.'FIDT') ;
 
*
******** On detrui les choses qui ne servent plus
*
*
*   Les variables primitives
*
    'DETR'  ( KINCO . 'V');
    'DETR'  ( KINCO . 'P');
    'DETR'  ( KINCO . 'T');
    'DETR'  ( KINCO . 'Y');
    'DETR'  ( KINCO . 'GAMMA');
    'OUBL'  KINCO  'V';
    'OUBL'  KINCO  'P';
    'OUBL'  KINCO  'T';
    'OUBL'  KINCO  'Y';
    'OUBL'  KINCO  'GAMMA';
*
*   Les MCHAML faces
*
    'DETR'  ROF ;
    'DETR'  VITF ;
    'DETR'  PF ;
    'DETR'  YF;
    'OUBL'  KINCO  'RNF';
    'OUBL'  KINCO  'VITNF';
    'OUBL'  KINCO  'PNF';
    'OUBL'  KINCO  'YF';
 
*
*   Les pentes 'ET' les limiteurs
*
    'SI' (ORD_ESP 'EGA' 2);
       'DETR' GRADR; 
       'DETR' GRADP;
       'DETR' GRADV;
       'DETR' GRADY;
       'DETR' ALR;
       'DETR' ALP;
       'DETR' ALV;
       'DETR' ALY;
    'FINSI' ;
 
 
*
******** On impose le CL
*
*
   'SI' LOGLIM;
      PROLIM RV ;
   'FINSI' ;   
 
 
    'SI' LOGQUIT;
        'QUITTER' BLOC1;
    'FINSI';
 
 
    'MENAGE' ;
 
'FIN'  BLOC1 ;
 
********************************************
**** FIN de Boucle Sur les Pas de Temps  ***
********************************************
 
'SINON'  ;
 
 
*********************************************************
****   Mono-espece,   boucle sur les pas de temps    ****
*********************************************************
 
 
 KINCO . 'V'  KINCO . 'P'  KINCO . 'T'  
   KINCO . 'GAMMA' =   'PRIM' 'PERFTEMP' (KINCO . 'IPGAZ')
  (KINCO . 'RNI') (KINCO . 'GNI') (KINCO . 'ENI') ;
 
 
 GRADR ALR COEFR = 'PENT' KDOMA 'CENTRE' 'EULESCAL'   'LIMITEUR'
                     (KINCO . 'RNI');
 
 GRADP ALP COEFP = 'PENT' KDOMA 'CENTRE' 'EULESCAL'   'LIMITEUR'
                     (KINCO . 'P');
 
 GRADV ALV COEFV = 'PENT' KDOMA 'CENTRE' 'EULEVECT'   'LIMITEUR'
                     (KINCO . 'V');
 
 
I = 0 ;
 
'REPETER'  BLOC1 (RV . 'ITMA')  ;
 
    I = I + 1 ;
 
*
***** Les variables primitives
*
         KINCO . 'V'  KINCO . 'P'  KINCO . 'T'  
         KINCO . 'GAMMA' =   'PRIM' 'PERFTEMP' (KINCO . 'IPGAZ')
         (KINCO . 'RNI') (KINCO . 'GNI') (KINCO . 'ENI') ;
 
         'SI'  (ORD_ESP 'EGA'  1) ;
 
             ROF VITF PF   =  'PRET' 'PERFTEMP'
                 ORD_ESP ORD_TPS KDOMA (KINCO . 'IPGAZ') (KINCO . 'RNI')  
                 (KINCO . 'V') (KINCO . 'P')   ;
 
         'SINON';
 
*       
***** Ordre 2 en espace => calcul des pentes
*
 
            GRADR ALR  = 'PENT' KDOMA 'CENTRE' 'EULESCAL'   'LIMITEUR'
                     (KINCO . 'RNI')  'GRADGEO' COEFR ;
 
            GRADP ALP  = 'PENT' KDOMA 'CENTRE' 'EULESCAL'   'LIMITEUR'
                     (KINCO . 'P')   'GRADGEO' COEFP ;
 
            GRADV ALV  = 'PENT' KDOMA 'CENTRE' 'EULEVECT'   'LIMITEUR'
                     (KINCO . 'V')   'GRADGEO' COEFV ;
 
 
           'SI' (ORD_TPS 'EGA' 1);
 
               ROF VITF PF   = 'PRET' 'PERFTEMP'  ORD_ESP ORD_TPS
                                   KDOMA (KINCO . 'IPGAZ')
                                   (KINCO . 'RNI') GRADR ALR
                                   (KINCO . 'V')   GRADV ALV
                                   (KINCO . 'P')   GRADP ALP ;
 
           'SINON' ;              
*       
********* Ordre 2 en temps
*
               ROF VITF PF     = 'PRET' 'PERFTEMP'  ORD_ESP ORD_TPS
                                  KDOMA (KINCO . 'IPGAZ')                                
                                  (KINCO . 'RNI') GRADR ALR
                                  (KINCO . 'V')   GRADV ALV
                                  (KINCO . 'P')   GRADP ALP
                                  (KINCO . 'GAMMA')                                   
                                  ((RV . 'PASDETPS' . 'DELTAT-1')/2.0);
           'FINSI' ;
 
        'FINSI'  ;
 
*
*********** Creation de MCHAML de type 'FACEL' pour les
*           calcul de flux aux interfaces
 
 
 
 KINCO . 'RNF'    =  ROF ;
 KINCO . 'VITNF'  =  VITF ;
 KINCO . 'PNF'    =  PF ;
 
*
********* Boucle sur les operateurs
*
 
   'REPETER'   BLOC2  NBOP ;
        NOMPER  = 'EXTRAIRE'  &BLOC2  (RV . 'LISTOPER');
        NOTABLE = 'MOT'  ('TEXTE'  ('CHAINE'  &BLOC2  NOMPER) ) ;
        ('TEXTE'  NOMPER) (RV . NOTABLE) ;
   'FIN'  BLOC2 ;
 
*
********* Mise a jour de la table RV . 'PASDETPS' 
*
 
   'SI' ('EXISTE'  RV  'DTI');
        DTI = 'MINIMUM'
            ('PROG' ((RV . 'DTI') '/' CFL)
            (RV . 'PASDETPS' . 'DTCONV') );
   'SINON';
        DTI = (RV . 'PASDETPS' . 'DTCONV');
   'FINSI';
 
    RV .  'PASDETPS' . 'DELTAT'  = DTI ;
 
    TMPS = RV . 'PASDETPS' . 'TPS';
    DTI0 = TFIN '-' TMPS;
    DTI0 = DTI0 '/' CFL;
 
    'SI' (DTI0 '<EG' DTI);
        DTI = DTI0;
        RV . 'PASDETPS' . 'DELTAT' = DTI ;
        LOGQUIT = VRAI;
    'SINON' ;
        LOGQUIT = FAUX ;
    'FINSI' ;
 
 
*
******** On avance au pas de temps suivant 
*
*   N.B. Tn+1 = Tn + (CFL *  RV . 'PASDETPS' . 'DELTAT');
*
 
    'AVCT'  RV CFL  'IMPR' (RV.'FIDT') ;
 
*
******** On detrui les choses qui ne servent plus
*
*
*   Les variables primitives
*
    'DETR'  ( KINCO . 'V');
    'DETR'  ( KINCO . 'P');
    'DETR'  ( KINCO . 'T');
    'DETR'  ( KINCO . 'GAMMA');
    'OUBL'  KINCO  'V';
    'OUBL'  KINCO  'P';
    'OUBL'  KINCO  'T';
    'OUBL'  KINCO  'GAMMA';
*
*   Les MCHAML faces
*
    'DETR'  ROF ;
    'DETR'  VITF ;
    'DETR'  PF ;
    'OUBL'  KINCO  'RNF';
    'OUBL'  KINCO  'VITNF';
    'OUBL'  KINCO  'PNF';
*
*   Les pentes 'ET' les limiteurs
*
    'SI' (ORD_ESP 'EGA' 2);
       'DETR' GRADR; 
       'DETR' GRADP;
       'DETR' GRADV;
       'DETR' ALR;
       'DETR' ALP;
       'DETR' ALV;
    'FINSI' ;
 
 
*
******** On impose le CL
*
*
   'SI' LOGLIM;
      PROLIM RV ;
   'FINSI' ;   
 
 
    'SI' LOGQUIT;
        'QUITTER' BLOC1;
    'FINSI';
 
 
    'MENAGE' ;
 
'FIN'  BLOC1 ;
 
********************************************
**** FIN de Boucle Sur les Pas de Temps  ***
********************************************
 
'FINSI' ;
 
'FINPROC'  ;
 
*****************************************************
*****************************************************
** FIN PROCEDURE EXEX                              **
*****************************************************
*****************************************************
 
 
*****************************************************
*****************************************************
***** PROCEDURE PROLIM                          *****
*****************************************************
*****************************************************
 
*
**** Cas Euler mono-espece
*
 
 
'DEBPROC'  PROLIM ;
'ARGUMENT'  RV*TABLE ;
 
*
**** RV . 'DOMAINE'   -> KDOMA
*
 
*KDOMA  = RV . 'DOMAINE' ;
KDOMA  = RV . 'MODTOT' ;
KDOMA2 = RV . 'DOMAINE' ;
KINCO  = RV . 'INCO' ;
 
KINCO . 'RNI' =    'KCHT' KDOMA 'SCAL' 'CENTRE'
      (KINCO . 'RNI') (KINCO . 'RNLIM');
 
KINCO . 'GNI' =    'KCHT' KDOMA 'VECT' 'CENTRE'
      (KINCO . 'GNI') (KINCO . 'GNLIM');
 
KINCO . 'ENI' =    'KCHT' KDOMA 'SCAL' 'CENTRE'
      (KINCO . 'ENI') (KINCO . 'ENLIM');
 
 
'FINPROC' ;
 
*****************************************************
*****************************************************
***** FIN PROCEDURE PROLIM                      *****
*****************************************************
*****************************************************
 
 
************ 
* MAILLAGE *
************
 
RAF = 4 ;
NY = RAF ;
NX = 6 '*' RAF ;
 
L   = 1.0D0  ;
DX = L '/' NX '/' 2.0D0;
H   = NY '*' DX  ;
 
xD =  0.5D0 '*' L ;
xG = -1.0D0 '*' xD  ;
yH =  0.5D0 '*' H   ;
yB = -1.0D0 '*' yH ;
 
A1 = xG yB     ;
A2 = 0.0D0 yB  ;
A3 = xD yB     ;
A4 = xD yH     ;
A5 = 0.0D0 yH  ;
A6 = xG yH     ;
VECTG = XG 0.0D0 ;
VECTD = XD 0.0D0 ;
 
xBG = xG '-' DX;
XBD = xD '+' DX;
 
B1 = xBG yB;
B2 = xBG yH;
B3 = xBD yB;
B4 = xBD yH;
 
BAS1  = A1 'DROI' NX A2 ;
BAS2  = A2 'DROI' NX A3 ;
HAU2  = A4 'DROI' NX A5 ;
HAU1  = A5 'DROI' NX A6 ;
LATG  = B1 'DROI' NY B2 ;
LAT1  = A1 'DROI' NY A6 ;
LAT12 = A2 'DROI' NY A5 ;
LAT2  = A3 'DROI' NY A4 ;
LATD  = B3 'DROI' NY B4 ;
 
 
DOM1  = LAT12 'TRANSLATION' NX VECTG ;
 
DOM2  = LAT12 'TRANSLATION' NX VECTD ;
 
VECTFG = (-2.0D0 '*' DX) 0.0D0;
VECTFD = (2.0 '*' DX) 0.0D0;
 
FRONTG = LAT1 'TRANSLATION' 2 VECTFG;
FRONTG = FRONTG 'COUL' ROUG ;
 
FRONTD =  LAT2  'TRANSLATION' 2 VECTFD;
FRONTD = FRONTD 'COUL' VERT ;
 
*
*** Rotation
*
 
ANGLE = 30.0D0;
ORIG = 0.0D0 0.0D0;
 
'MESSAGE';
'MESSAGE' (CHAIN 'ANGLE = ' ANGLE);
'MESSAGE';
 
DOM1   = DOM1   'TOURNER' ANGLE ORIG;
DOM2   = DOM2   'TOURNER' ANGLE ORIG;
FRONTG = FRONTG 'TOURNER' ANGLE ORIG;
FRONTD = FRONTD 'TOURNER' ANGLE ORIG;
 
DOMINT = DOM1 'ET' DOM2 ;
'ELIMINATION' DOMINT 1D-6;
 
FRONT = FRONTG 'ET'  FRONTD ;
'ELIMINATION' FRONT 1D-6;
 
DOMTOT = DOMINT 'ET'  FRONT;
'ELIMINATION' DOMTOT 1D-6;
 
**********************
*** OBJETS MODELES ***
**********************
MDOMTOT = 'CHANGER' DOMTOT 'QUAF' ;
MDOMINT = 'CHANGER' DOMINT 'QUAF' ;
MDOM1   = 'CHANGER' DOM1   'QUAF' ;
MDOM2   = 'CHANGER' DOM2   'QUAF' ;
MFRONTG = 'CHANGER' FRONTG 'QUAF' ;
MFRONTD = 'CHANGER' FRONTD 'QUAF' ;
MFRONT  = 'CHANGER' FRONT  'QUAF' ;
'ELIM' (MDOMTOT 'ET' MDOMINT 'ET' MDOM1 'ET' MDOM2 'ET' 
        MFRONTG 'ET' MFRONTD 'ET' MFRONT) 1.E-5;
MDNS    = 'EULER' ;
$DOMTOT = 'MODE' MDOMTOT MDNS  ;
$DOMINT = 'MODE' MDOMINT MDNS  ;
$DOM1   = 'MODE' MDOM1   MDNS  ;
$DOM2   = 'MODE' MDOM2   MDNS  ;
$FRONTG = 'MODE' MFRONTG MDNS  ;
$FRONTD = 'MODE' MFRONTD MDNS  ;
$FRONT  = 'MODE' MFRONT  MDNS  ;
 
*
******* Creation de la ligne Utilisée pour le Post-Traitement
*       reliant les points centres
*
 
  XINIT = XG '+' (0.5D0 '*' DX) ;
  YINIT = YB '+' (0.5D0 '*' DX) ;
  XFIN  = XD '-' (0.5D0 '*' DX) ;
  YFIN = YINIT ;
  PINI = XINIT YINIT;
  PFIN = XFIN YFIN ;
  IAUX = (2 '*' NX) '-' 1 ;
  COURB = PINI 'DROIT' IAUX PFIN;
  COURB = COURB 'TOURNER' ANGLE ORIG;
  COURB = COURB 'COULEUR' 'VERT';
  'ELIMINATION'  0.001 Courb ('DOMA' $DOMTOT 'CENTRE') ;
 
'SI' GRAPH ;
   'TRACER' (('DOMA' $DOMTOT 'MAILLAGE') 'ET' COURB)
   'TITRE' ('CHAINE' 'Maillage ');
'FINSI' ;   
 
*
*** Etats exacte
*
 
ro1 = rog ;
u1 = ug ;
p1 = pg ;
ETHER1 = ETHERg ;
T1 = Tg ;
 
ro2 = rod ;
u2 = ud ;
p2 = pd ;
ETHER2 = ETHERd ;
T2 = Td ;
 
*
**** On inverse les etats, pour controller qu'on a le meme resultats
*
 
'SI' VRAI ;
   rog = ro2 ;
   ug = -1.0 '*' u2 ;
   pg = p2 ;
   ETHERg = ETHER2 ;
   Tg = T2 ;
 
   rod = ro1 ;
   ud = -1.0 '*' u1 ;
   pd = p1 ;
   ETHERd = ETHER1 ;
   Td = T1 ;
 
'FINSI' ;   
 
*
*** Etat gauche (rog, ug, ETHERg, Tg , Pg deja definis)
*
 
ung   = ug  ;
utg   = 0.0 ;
retg  = rog '*' (ETHERg '+' (0.5 '*' ung '*' ung)) ;
 
*
 
rouxg  = ((ung '*' ('COS' ANGLE)) '-'
          (utg '*' ('SIN' ANGLE))) '*' rog ;
 
rouyg  = ((ung '*' ('SIN' ANGLE)) '+'
          (utg '*' ('COS' ANGLE))) '*' rog;
 
 
*
*** Etat droite
*
 
und  = ud ;
utd  = 0.0D0 ;
retd = rod '*' (ETHERd '+' (0.5 '*' und '*' und)) ;
 
rouxd = ((und '*' ('COS' ANGLE)) '-'
         (utd '*' ('SIN' ANGLE))) '*' rod;
rouyd = ((und '*' ('SIN' ANGLE)) '+'
         (utd '*' ('COS' ANGLE))) '*' rod;
 
*
*** ro
*
 
RO_f1  = 'KCHT'  $FRONTG  'SCAL' 'CENTRE' rog ;
RO_f2  = 'KCHT'  $FRONTD  'SCAL' 'CENTRE' rod ;
RO_f   = 'KCHT'  $FRONT   'SCAL' 'CENTRE' (RO_f1 'ET' RO_f2) ; 
 
RO_i1  = 'KCHT'  $DOM1   'SCAL' 'CENTRE' rog;
RO_i2  = 'KCHT'  $DOM2   'SCAL' 'CENTRE' rod;
RO_i   = 'KCHT'  $DOMINT 'SCAL' 'CENTRE' (RO_i1 'ET' RO_i2);
 
RN  = 'KCHT'  $DOMTOT 'SCAL' 'CENTRE' (RO_i 'ET'  RO_f) ;
 
*
*** ro u, ro v
*
 
GN_f1  = 'KCHT'  $FRONTG   'VECT' 'CENTRE' (rouxg rouyg);
GN_f2  = 'KCHT'  $FRONTD   'VECT' 'CENTRE' (rouxd rouyd);
GN_f   = 'KCHT'  $FRONT    'VECT' 'CENTRE' (GN_f1 'ET' GN_f2);
 
GN_i1  = 'KCHT'  $DOM1    'VECT' 'CENTRE' (rouxg rouyg);
GN_i2  = 'KCHT'  $DOM2    'VECT' 'CENTRE' (rouxd rouyd);
GN_i   = 'KCHT'  $DOMINT  'VECT' 'CENTRE' (GN_i1 'ET' GN_i2);
 
GN  = 'KCHT'  $DOMTOT 'VECT' 'CENTRE' (GN_i 'ET'  GN_f) ;
 
*
*** ro e 
*
 
RE_f1  = 'KCHT'  $FRONTG  'SCAL' 'CENTRE'  retg ;
RE_f2  = 'KCHT'  $FRONTD  'SCAL' 'CENTRE'  retd ;
RE_f   = 'KCHT'  $FRONT   'SCAL' 'CENTRE' (RE_f1 'ET' RE_f2) ; 
 
RE_i1  = 'KCHT'  $DOM1   'SCAL' 'CENTRE'  retg ;
RE_i2  = 'KCHT'  $DOM2   'SCAL' 'CENTRE'  retd ;
RE_i   = 'KCHT'  $DOMINT 'SCAL' 'CENTRE' (RE_i1 'ET' RE_i2);
 
REN = 'KCHT'  $DOMTOT 'SCAL' 'CENTRE' (RE_i ET RE_f) ;
 
 
********************************************************
*** CREATION DE 'MODE' POUR GRAPHIQUER LE CHAMELEM   ***
********************************************************
 
MOD1     =  'MODELISER'  ('DOMA' $DOMTOT 'MAILLAGE') 'THERMIQUE' ;
 
*
**** Les variables primitives
*
 
VN PN TN GAMN =   'PRIM' 'PERFTEMP' PGAZ  RN GN REN  ;
 
*
**** Enthalpie totale
*
 
 
HTN = ((REN '+' PN) '/' RN) ;
 
*
**** Les vitesses dans le repaire tube
*
 
 
VX = 'EXCO' 'UX  ' VN ;
VY = 'EXCO' 'UY  ' VN ;
 
VNN = (VX * ('COS' ANGLE)) '+'
      (VY * ('SIN' ANGLE)) ;
 
VNT = (VY * ('COS' ANGLE)) '-'
      (VX * ('SIN' ANGLE)) ;
 
*
*** GRAPHIQUE DES C.I.
*
 
'SI' GRAPH ;
* 'SI' FAUX ;
 
*
*** CREATION DE CHAMELEM
*
  CHM_RN   =  'KCHA' $DOMTOT 'CHAM' RN ;
  CHM_VNN  =  'KCHA' $DOMTOT 'CHAM' VNN ;
  CHM_VNT  =  'KCHA' $DOMTOT 'CHAM' VNT ;
  CHM_PN   =  'KCHA' $DOMTOT 'CHAM' PN ;
  CHM_TN   =  'KCHA' $DOMTOT 'CHAM' TN ;
  CHM_HTN  =  'KCHA' $DOMTOT 'CHAM' HTN ;
 
  TRAC CHM_RN  MOD1 'TITR'  ('CHAINE'  'RO at t=' 0.0);
  TRAC CHM_PN MOD1 'TITR'  ('CHAINE'   'P at t=' 0.0);
  TRAC CHM_VNN MOD1 'TITR'  ('CHAINE'  'VNN at t=' 0.0);
  TRAC CHM_VNT MOD1 'TITR'  ('CHAINE'  'VNT at t=' 0.0);
  TRAC CHM_TN MOD1 'TITR'  ('CHAINE'   'T at t=' 0.0);
  TRAC CHM_HTN MOD1 'TITR' ('CHAINE'   'ht at t=' 0.0);
'FINSI' ;
 
 
***********************
**** La table EQEX **** 
***********************
 
 
 CFL   = 0.8D0  ;
 NITER = 50 ;
 TFIN  = 10000 ;
 
 RV = 'EQEX' ('DOMA' $DOMTOT 'TABLE') 'ITMA' NITER 'ALFA' CFL
      'TPSI' 0. 'TFINAL' TFIN
*
***** Option VF       = volumes finis
*            CONS     = conservative
*            EULER    = euler monoespece
*            VANLEER  = methode
*
      'OPTI' 'VF' 'CONS' 'EULERMST' METO
*
***** Procedure CALCUL
*      
      'ZONE' $DOMTOT 'OPER' 'CALCUL'
*
***** Operateur 'KONV'
*      
      'ZONE' $DOMTOT 'OPER' 'KONV'
*
***** Arguments de 'KONV' (maximum 8 chatacters)
*      
      'IPGAZ' 'RNF' 'VITNF' 'PNF'
*
***** LIST des inconnues
*
      'INCO' 'RNI' 'GNI' 'ENI'   ;
 
 
*
*** La table RV . INCO (de soustype INCO);
*
 
RV . 'INCO'  = TABLE  'INCO' ;
 
*
*** Stokage des inconnues des equations d'Euler.
*
 
RV . 'INCO' . 'RNI' = 'COPIER' RN ;
RV . 'INCO' . 'GNI' = 'COPIER' GN ;
RV . 'INCO' . 'ENI' = 'COPIER' REN ;
 
 
*
*** Le gaz
*
 
RV . 'INCO' . 'IPGAZ' =  PGAZ ;
 
*
*** CONDITIONS LIMITS
*
 
RV . 'INCO' . 'CLIM'   = VRAI  ;
 
RV . 'INCO' . 'RNLIM' = 'COPIER' RO_f ;
RV . 'INCO' . 'GNLIM' = 'COPIER' GN_f ;
RV . 'INCO' . 'ENLIM' = 'COPIER' RE_f ;
 
*
**** Ordre en espace
*    Ordre en temps     
*
 
IE = 2 ;
IT = 2 ;
 
RV . 'ORDREESP' = IE ;
RV . 'ORDRETPS' = IT ;
 
*
**** Test de la convergence
*
 
RN0 = 'COPIER' RN ;
RV . 'INCO' . 'RN0' = 'COPIER' RN ;
RV . 'INCO' . 'IT' = 'PROG' ;
RV . 'INCO' . 'ER' = 'PROG' ;
 
*???
RV . 'MODTOT' = $DOMTOT ;
 
*
***********************************
*** Execution EXEX              ***
***********************************
*
 
 
 'MESSAGE'  ;
 'MESSAGE'  '**************************';
 'MESSAGE'   ('CHAINE'  'METHODE : ' METO) ;
 'MESSAGE'  '**************************';
 'MESSAGE' ;
 
 
 
'TEMPS' ;
   EXEX RV ;
'TEMPS' ;
 
 
TFIN = RV. 'PASDETPS'. 'TPS'; 
 
*
**** SOLUTIONS
*
 
*
**** Les variables conservatives
*
 
RN  = 'COPIER' (RV . 'INCO' . 'RNI');
GN  = 'COPIER' (RV .  'INCO'. 'GNI'); 
REN = 'COPIER' (RV .  'INCO'. 'ENI');
 
 
*
**** Les variables primitives
*
 
VN PN TN GAMN =   'PRIM' 'PERFTEMP' PGAZ  RN GN REN  ;
 
*
**** Enthalpie totale
*
 
HTN = ((REN '+' PN) '/' RN)  ;
 
*
**** Les vitesses dans le repaire tube
*
 
VX = 'EXCO' 'UX  ' VN ;
VY = 'EXCO' 'UY  ' VN ;
 
VNN = (VX * ('COS' ANGLE)) '+'
      (VY * ('SIN' ANGLE)) ;
 
VNT = (VY * ('COS' ANGLE)) '-'
      (VX * ('SIN' ANGLE)) ;
 
*
*** GRAPHIQUE DES SOLUTIONS
 
*
 
'SI' GRAPH ;
 
* 'SI' FAUX ;
*
*** CREATION DE CHAMELEM
*
  CHM_RN   =  'KCHA' $DOMTOT 'CHAM' RN ;
  CHM_VNN  =  'KCHA' $DOMTOT 'CHAM' VNN ;
  CHM_VNT  =  'KCHA' $DOMTOT 'CHAM' VNT ;
  CHM_PN   =  'KCHA' $DOMTOT 'CHAM' PN ;
  CHM_TN   =  'KCHA' $DOMTOT 'CHAM' TN ;
  CHM_HTN  =  'KCHA' $DOMTOT 'CHAM' HTN ;
  TRAC CHM_RN  MOD1 'TITR'  ('CHAINE'  'RO at t=' TFIN);
  TRAC CHM_PN MOD1 'TITR'  ('CHAINE'   'P at t=' TFIN);
  TRAC CHM_VNN MOD1 'TITR'  ('CHAINE'  'VNN at t=' TFIN);
  TRAC CHM_VNT MOD1 'TITR'  ('CHAINE'  'VNT at t=' TFIN);
  TRAC CHM_TN MOD1 'TITR'  ('CHAINE'   'T at t=' TFIN);
  TRAC CHM_HTN MOD1 'TITR' ('CHAINE'   'ht at t=' TFIN);
'FINSI' ;
 
 
*
*** Objects evolutions
*
 
xx yy = 'COORDONNEE' Courb;
ss =  'KOPS' ('KOPS' ('COS' ANGLE) '*' xx) '+'
             ('KOPS' ('SIN' ANGLE) '*' yy);
 
evxx = 'EVOL' 'CHPO' ss Courb ;
lxx = 'EXTRAIRE'  evxx 'ORDO' ;
 
x0 = 'MINIMUM' lxx;
x1 = 'MAXIMUM' lxx ;
 
* ro
 
evro = 'EVOL'  'CHPO' RN Courb ;
lro  = 'EXTRAIRE'  evro 'ORDO' ;
evro = 'EVOL' 'MANU' 'x' lxx 'ro' lro;
tro  = CHAINE  '1D ' METO  ' : RO   IT '  IT     '  IE ' IE   
      '  tmps ' TFIN  ' elem ' 'QUA4' ;
 
* u
 
evu  = 'EVOL' 'CHPO' VNN Courb ;
lu   = 'EXTRAIRE'  evu 'ORDO' ;
evu  = 'EVOL' 'MANU' 'x' lxx 'u' lu ;
tu   = CHAINE  '1D ' METO ' : U  IT '  IT     '  IE ' IE   
      '  tmps ' TFIN  ' elem ' 'QUA4' ;
 
* v
 
evv  = 'EVOL' 'CHPO' VNT Courb ;
lv   = 'EXTRAIRE'  evv 'ORDO' ;
evv  = 'EVOL' 'MANU' 'x' lxx 'v' lv ;
tv   = CHAINE  '1D ' METO ' : V  IT '  IT     '  IE ' IE   
      '  tmps ' TFIN  ' elem ' 'QUA4' ;
 
'SI' GRAPH ;      
  'DESSIN' evv 'TITRE' tv 'XBOR' x0 x1;
'FINSI' ;
 
* p
 
evp  = 'EVOL' 'CHPO' PN Courb ;
lp   = 'EXTRAIRE'  evp 'ORDO' ;
evp  = 'EVOL' 'MANU' 'x' lxx 'p' lp ;
tp   = CHAINE  '1D ' METO ' : P  IT '  IT     '  IE ' IE   
 
      '  tmps ' TFIN  ' elem ' 'QUA4' ;
 
 
* ht
 
evht = 'EVOL' 'CHPO' HTN Courb ;
lht  = 'EXTRAIRE'  evht 'ORDO' ;
evht = 'EVOL' 'MANU' 'x' lxx 'ht' lht ;
tht   = CHAINE  '1D ' METO ' : ht  IT '  IT     '  IE ' IE   
 
      '  tmps ' TFIN  ' elem ' 'QUA4' ;
 
*
*** Solution analitique 
 
lpan = 'PROG' ;
lroan = 'PROG' ;
luan = 'PROG' ;
lhtan = 'PROG' ;
 
hg = ETHERg '+' (pg '/' rog) '+' (0.5 '*' ug '*' ug) ;
 
'REPETER' BLOC2 ('DIME' lxx);
   x  = 'EXTRAIRE' lxx &BLOC2;
   'SI' (x < 0) ;
      lroan = lroan 'ET' ('PROG' rog);
      lpan = lpan 'ET' ('PROG' pg);
      luan = luan 'ET' ('PROG' ug);
      lhtan = lhtan 'ET' ('PROG' hg);
   'SINON' ;
      lroan = lroan 'ET' ('PROG' rod);
      lpan = lpan 'ET' ('PROG' pd);
      luan = luan 'ET' ('PROG' ud);
      lhtan = lhtan 'ET' ('PROG' hg);
   'FINSI' ;      
'FIN' BLOC2;   
 
 
evroa  = 'EVOL' 'MANU' 'x' lxx 'ro' lroan;
evpa   = 'EVOL' 'MANU' 'x' lxx 'p'  lpan;
evua   = 'EVOL' 'MANU' 'x' lxx 'ua' luan;
evha   = 'EVOL' 'MANU' 'x' lxx 'ht' lhtan;
 
'SI' GRAPH ;
    TAB1=TABLE;
    TAB1.'TITRE'= TABLE ;
    TAB1.1='MARQ TRIB REGU';
    TAB1.'TITRE' . 1  = MOT 'Numerique' ;
    TAB1.2='MARQ CROI REGU';
    TAB1.'TITRE' . 2  = MOT 'Exacte' ;
    'DESSIN' (evro 'ET' evroa) 'TITRE' tro 'XBOR' x0 x1
    'LEGE' TAB1;
    'DESSIN' (evp 'ET' evpa) 'TITRE' tp 'XBOR' x0 x1
    'LEGE' TAB1;
    'DESSIN' (evu 'ET' evua) 'TITRE' tu 'XBOR' x0 x1
    'LEGE' TAB1;
    'DESSIN' (evht 'ET' evha) 'TITRE' tht 'XBOR' x0 x1
    'LEGE' TAB1;
'FINSI' ;
 
*
**** Tests
*
 
* Convergence
 
ERCON = 'EXTRAIRE' ('DIME' (RV . 'INCO' . 'ER')) (RV . 'INCO' . 'ER') ;
 
'SI' (ERCON > -3) ;
   'MESSAGE' 'Probleme de convergence' ;
   'ERREUR' 5;
'FINSI' ;
 
'SI' GRAPH ;
      EVOL4 = EVOL 'MANU' 'ITERATIONS' (RV.INCO.'IT') 'LOG|E|inf'
              (RV.INCO.'ER') ;
     'DESSIN' EVOL4 'XBOR' 0. (2. '*' niter) 'YBOR' -12.0 0.0
     'TITRE' ('CHAINE' METO ' IE =' IE '  IT =' IT) ;
'FINSI' ;
 
 
ERRO  = 'ABS' (LRO '-' LROAN) ;
ERP   = 'ABS' (LP  '-' LPAN)  ;
ERU   = 'ABS' (LU  '-' LUAN)  ;
ERH   = 'ABS' (LHT '-'  LHTAN)  ;
 
L1RO = 0.0 ;
L1P  = 0.0 ;
L1U  = 0.0 ;
L1H  = 0.0 ;
 
NDIM = 'DIME' erro ;
 
'REPETER' BL1 NDIM ;
   L1RO = L1RO '+' ('EXTRAIRE' &BL1 ERRO) ;
   L1P  = L1P  '+' ('EXTRAIRE' &BL1 ERP) ;
   L1U  = L1U  '+' ('EXTRAIRE' &BL1 ERU) ;
   L1H  = L1H  '+' ('EXTRAIRE' &BL1 ERH) ;
'FIN' BL1 ;
 
L1RO = L1RO '/' (NDIM * ('MAXIMUM' lroan)) ;
 
 
L1P  = L1P '/' (NDIM * ('MAXIMUM' lpan)) ;
 
L1U  = L1U '/' (NDIM * ('MAXIMUM' luan 'ABS' )) ;
L1H  = L1H '/' (NDIM * ('MAXIMUM' lhtan)) ;
 
 
ERRMAX = 'MAXIMUM' ('PROG' L1RO L1P L1U L1H) ;
 
 
'SI' (ERRMAX > 1.0D-4);
   'ERREUR' 5;
'FINSI' ;  
 
'FIN' ;
 
 
 
 
 
 
 
 
 

© Cast3M 2003 - Tous droits réservés.
Mentions légales