$$$$ EXEC
* EXEC      PROCEDUR  MAGN      11/09/12    21:16:08     7117           
*X EXEC (Procedure)
*   Procedure EXEC
*
*   Objet   :  Execute un algorithme décrit dans une table RV
*              de type EQEX.
*              Cette table est créée par l'opérateur EQEX
*
*   Syntaxe :  EXEC RV ;
*
***********************************************************************
* VERSION    : ????
* HISTORIQUE : 20/12/99: gounand
* Rajout de la gestion de la matrice servant à l'assemblage
* (rv . 'METHINV' . 'MATASS') 
*
* HISTORIQUE : 12/05/06 : ajout STOPITER, NUITER et STOPPDT
* HISTORIQUE : 21/12/07 : ajout projection algébrique incrémentale
*                         (cf. GRESP)
* HISTORIQUE :
************************************************************************
* DISCPRES = (d'apres KMIC)
* KPRE=3 pression P0   KPRE=4 pression P1  KPRE=2 cas macro 1ère génération
* KPRE=5 pression MSOMMET
************************************************************************
'DEBPROC'  EXEC ;
'ARGUMENT' rv*'TABLE   '  nparti/'ENTIER' lopp/'LISTMOTS';
*
*
*
* Préparation pour paralléliser le comportement
*
Si(EXIST RV 'KPME');
KPME=rv.'KPME';
Sinon;
KPME=FAUX;
Finsi;

si (exist nparti);npart=nparti;sinon;npart=valeur assis;finsi;
Si(NON (KPME));npart=1;Finsi;

* On partitionne en fonction du nombre de proc ou plus voir moins
* si précisé dans l'instruction exec rv n ;
* Si Npart <eg 1 : On ne partitionne pas
Si (Npart <eg 1); mess ' Npart <eg 1 : On ne partitionne pas ';
Sinon; mess ' On partitionne en ' Npart; Finsi;
si (exist lopp);loperpar=lopp;sinon;
*loperpar=mots 'NS' 'TSCA' 'DFDT' 'KONV' 'LAPN' 'FROT' 'FIMP' 'TOIMP';
loperpar=mots 'NS' 'TSCA' 'DFDT' ;
finsi;
* La liste des opérateurs à partitionner est donnée dans loperpar
* ou surchargée si lopp est donné.
*********
*                                                                               

KTBPART=VRAI;
Si('EXIST' RV 'TBPART');
nparta=RV.'TBPART'.'NPART';
Si (Npart EGA nparta);KTBPART=FAUX; Finsi;
Finsi ;

Si KTBPART ;
mess ' On crée ou recrée la table TBPART';
* On crée une table des opérateurs partitionnés
* Faire attention si on veut inclure la procédure CALCUL, de ne pas répéter
* les opérations
TBPART='TABLE';
TBPART.'NPART'=npart;
rv.'TBPART'=TBPART;
nop=dime (rv.'LISTOPER');
repeter BPART nop;                                                              
nomper=extr (rv.'LISTOPER') &BPART ;
notable= 'MOT' ('TEXTE' ('CHAINE' &BPART nomper)) ;
TBPART.notable=rv.notable;
KPAR='EXIST' loperpar nomper;
kform=rv.notable.'KOPT'.'KFORM';
SI (NON ( EGA kform 1)); KPAR=FAUX; Finsi;
TNOT='TYPE' (rv.notable.'DOMZ');
SI(NON(EGA TNOT 'MMODEL'));KPAR=FAUX;Finsi;
*Si Kpar ; mess ' KPAR = VRAI'; sinon ; mess ' KPAR = FAUX'; finsi ;
Si (('EGA' (rv.notable.'SOUSTYPE') 'KIZX') et KPAR);
$mti=rv.notable.'DOMZ';
Si ('EXIST' TBPART $mti);                                                       
$mtx=TBPART.$mti;
Sinon;
Si (Npart > 1);
 $mtx= part npart $mti 'NOOPT';
Sinon ;
 $mtx= $mti;
Finsi;
TBPART.$mti=$mtx;
Finsi ;

Si (Npart > 1);
$tns=table 'ESCLAVE';   

repeter bloc1 npart;
$mti=$mtx.&bloc1;
tnsi=copie (rv.notable);
tnsi.'DOMZ'=$mti;
tnsi= enlev tnsi 'TDOMZ';
$tns.&bloc1= tnsi;
fin bloc1;

TBPART.notable=$tns;
*Sinon;
*TBPART.notable=notable;
Finsi;

Finsi;
fin BPART;

Finsi;


TBPART=rv.'TBPART';
*
* Fin Préparation pour paralléliser le comportement
*
*
*

*Logique d'arret
'SI' ('NON' ('EXISTE' rv 'STOPITER')) ;
   rv . 'STOPITER' = FAUX ;
'FINSI' ;
'SI' ('NON' ('EXISTE' rv 'STOPPDT')) ;
   rv . 'STOPPDT' = FAUX ;
'FINSI' ;

* Logique pilotant le recalcul de la matrice de pression
 CALPRE=FAUX                                                   ;
 Si ('EXIST' rv 'CALPRE')                                      ;
  vertytab rv 'CALPRE' 'LOGIQUE'                               ;
  CALPRE=rv.'CALPRE'                                           ;
 Finsi                                                         ;

 Si(NON('EXIST' rv 'XEQUA'));
 rv.'XEQUA'= FAUX ;
 Finsi ;


'SI' ('EGA'  (rv . 'NAVISTOK') 0) ;
  EXAC rv ;
'QUITTER' EXEC ;
'FINSI'  ;

*******PROJ*************************************************************
*******PROJ Trois possibilités *****************************************
*******PROJ TMDM1 (vrai)_> Correction Gresho                  NON
*******PROJ TGRAD (vrai)-> Formulation Gradient               NON
*******PROJ TPNM2 (vrai)-> elimination end of step velocity   OUI
*******PROJ TMDM1 = VRAI TGRAD = FAUX TPNM2 = FAUX  ancien algo
         TYPROJ='VPI1'              ;
         ROW   = 1.                 ;
       'SI' ('EXIST' RV 'TYPROJ')   ;
         TYPROJ=RV.'TYPROJ'         ;
          'SI' (non (exist (MOTS 'PSCT' 'VPI1' 'VPI2' 'PENA') TYPROJ));
 Mess '*********************************************************' ;
 Mess ' ERREUR   ERREUR   ERREUR   ERREUR   ERREUR   ERREUR     ' ;
 Mess '                                                         ' ;
             Mess 'Le mot ' TYPROJ ' n existe pas dans la liste.' ;
             Mess 'Les méthodes de projection autorisées sont :'  ;
             Mess 'VPI1, VPI2, PSCT et PENA'                      ;
 Mess '*********************************************************' ;
             erreur 21              ;
             quitter EXEC           ;
          'FINSI'                   ;
       'FINSI'                      ;

       'SI' ('EGA' TYPROJ 'PSCT')   ;
            TGRAD = FAUX ;
            TMDM1 = FAUX ;
            TPNM2 = FAUX ;
       'FINSI'                      ;
       'SI' ('EGA' TYPROJ 'VPI1')   ;
            TGRAD = FAUX ;
            TMDM1 = VRAI ;
            TPNM2 = FAUX ;
       'FINSI'                      ;
       'SI' ('EGA' TYPROJ 'VPI2')   ;
            TGRAD = FAUX ;
            TMDM1 = FAUX ;
            TPNM2 = VRAI ;
       'FINSI'                      ;
       'SI' ('EGA' TYPROJ 'PENA')   ;
            TGRAD = FAUX ;
            TMDM1 = FAUX ;
            TPNM2 = FAUX ;
       'FINSI'                      ;

*******PROJ Quatres possibilités ***************************************
*******PROJ*************************************************************

  nomvi=rv . 'NOMVI'  ;

'SI'  ('EGA' ('VALEUR' 'DIME') 2) ;
  vnul=0.D0 0.D0 ;
  vuni=1.D0 1.D0 ;
  cnvi1= 'MOT' ('TEXTE' ('CHAINE' 1 nomvi)) ;
  cnvi2= 'MOT' ('TEXTE' ('CHAINE' 2 nomvi)) ;
  lc= mots cnvi1 cnvi2 ;
  lcu=mots 'UX' 'UY' ;
'SINON' ;
  vnul=0.D0 0.D0 0.D0 ;
  vuni=1.D0 1.D0 1.D0 ;
  cnvi1= 'MOT' ('TEXTE' ('CHAINE' 1 nomvi)) ;
  cnvi2= 'MOT' ('TEXTE' ('CHAINE' 2 nomvi)) ;
  cnvi3= 'MOT' ('TEXTE' ('CHAINE' 3 nomvi)) ;
  lc= mots cnvi1 cnvi2 cnvi3 ;
  lcu=mots 'UX' 'UY' 'UZ' ;
'FINSI' ;
'SI' ('NON' ('EXISTE' rv 'OMEGA')) ;
  omeg=1.D0 ;
'SINON' ;
  omeg=rv . 'OMEGA' ;  
'FINSI' ;

testpr ='EXISTE' rv 'PRESSION' ;
testprj ='EXISTE' rv 'PROJ' ;
testran=testpr 'ET' ('EXISTE' rv 'CO') ;

  si testpr ;  rvp  = rv.'PRESSION' ;
  achp  matpr= 'KOPS' 'MATRIK' ;
  rvp.'MATP'= matpr;
  finsi ;
  si testprj;  rvp  = rv.'PROJ' ; Finsi ;

  'SI' (rv.'IMPR' >EG 1) ;
  si testpr ;
  mess 'Algorithme semi explicite (ANCIEN)'; 
  mess '==================================';
  Finsi ;
  si testprj;
  mess 'Algorithme de Projection';
  mess '========================';
  finsi ;
  si ( non (testpr ou testprj)) ;
  mess 'Algorithme standard implicite ou explicite ';
  mess '===========================================';
  finsi ;
  'FINSI' ;

'SI' ('NON' ('EXISTE' rv 'HIST')) ;
  rv . 'HIST' = 'TABLE' ;
'FINSI' ;


ITMA=(rv . 'ITMA') ;
IMPTCRR=0 ;
IMPKRES=0 ;
'SI' ('<EG' ITMA 1) ;
  ITMA=1 ;
  IMPTCRR=1 ;
'FINSI'  ;
 Si (testprj ) ; IMPTCRR=1         ;  finsi ;
 Si ( non (testpr ou testprj)) ; IMPTCRR=1 ; finsi ;


* Gestion de la matrice de préconditionnement
*
* calprec : doit-on recalculer le préconditionneur
* fcprect : fréquence de recalcul du préc. en fn du pas de temps
* fcpreci : fréquence de recalcul du préc. dans la boucle
*           d'itérations internes pour les non-linéarités
* fcprectp : idem pour la matrice de pression
* fcprecip : idem pour la matrice de pression
* resmn : le résidu au pas de temps ou à l'itération interne précédente
*
calprec = VRAI ;
fcprect = rv. 'METHINV' . 'FCPRECT' ;
fcpreci = rv. 'METHINV' . 'FCPRECI' ;
resmn maprec = 'KOPS' 'MATRIK' ;
resmn1 maprec1 = 'KOPS' 'MATRIK' ;
resmn2 maprec2 = 'KOPS' 'MATRIK' ;
'SI'(exist rv 'resmn');resmn=rv.'resmn' ; finsi ;
'SI' ('NON' ('EXIST' rv.'METHINV' 'CALPREC')) ;
  rv.'METHINV'.'CALPREC'= VRAI ;
'FINSI' ;
'SI' testprj ;
  fcprectp = rvp. 'METHINV' . 'FCPRECT' ;
  fcprecip = rvp. 'METHINV' . 'FCPRECI' ;
  'SI' ('NON' ('EXIST' rvp.'METHINV' 'CALPREC')) ;
    rvp.'METHINV'.'CALPREC'= VRAI ;
  'FINSI' ;
'FINSI' ;

'REPETER' bloc1 ITMA ;
  'SI' ('MULT' &bloc1 fcprect) ;
    rv.'METHINV'.'CALPREC' = VRAI ;
  'FINSI' ;
  'SI' (testprj) ;
    'SI' ('MULT' &bloc1 fcprectp) ;
      rvp.'METHINV'.'CALPREC' = VRAI ;
    'FINSI' ;
  'FINSI' ;
   testp1 = FAUX  ;
  si (testpr ou testprj) ;
   Si CALPRE  ; testp1=VRAI ; Finsi ;
   Si (non (exist rvp 'MATC')) ; testp1=VRAI ; Finsi ;
  finsi ;

  'REPETER' bloci (rv . 'NITER') ;
    rv . 'NUITER' = &bloci ;
   'SI' ('MULT' &bloci fcpreci) ;
      rv.'METHINV'.'CALPREC' = VRAI ;
    'FINSI' ;
  
    'SI' (testprj) ;
      'SI' ('MULT' &bloci fcprecip) ;
        rvp.'METHINV'.'CALPREC' = VRAI ;
      'FINSI' ;
    'FINSI' ;
    st mat  = 'KOPS' 'MATRIK' ;
    sf mau  = 'KOPS' 'MATRIK' ;

      mdfdt = 0 ;
    'REPETER' bloc2 ('DIME' (rv . 'LISTOPER')) ;
      nomper = 'EXTRAIRE' &bloc2 (rv . 'LISTOPER') ;
      notable= 'MOT' ('TEXTE' ('CHAINE' &bloc2 nomper)) ;

*     mess 'Operateur Bloc2 mdfdt ? ' nomper ;

         rvn=rv . notable;
         rvna=TBPART. notable;
         TASSI=FAUX;
      Si('EXIST' rvna 'SOUSTYPE');
       Si('EGA' (rvna.'SOUSTYPE') 'ESCLAVE');
        Si(Npart > 1); TASSI=VRAI;
        finsi;
       finsi;
      finsi;
      mdfdt = mdfdt + rvn . 'KOPT' . 'KFORM' ;

      si (ega nomper 'DFDT    ');
      ISCHT= (rvn.kopt.'ISCHT') ;
      Si TASSI ;
      msi mai=ASSI TOUS ('TEXTE' nomper) (rvna) ;
      msi=et msi;
      mai=et mai;
      Sinon;
      msi mai=('TEXTE' nomper) (rvn) ;
      Finsi;
      mat = mat 'ET' mai ;
      st  = st  'ET' msi ;

      sinon ;

*     mess ' 1OPER =' nomper;
      Si TASSI ;
      msi mai=ASSI TOUS ('TEXTE' nomper) (rvna) ;
      msi=et msi;
      mai=et mai;
      Sinon;
      msi mai=('TEXTE' nomper) (rvn) ;
      Finsi;
      mau = mau 'ET' mai ;
      sf  = sf  'ET' msi ;

      finsi ;

    'FIN' bloc2 ;

      s2 = sf et st ;
      ma1 = mau 'ET' mat ;

**********  Traitement PRESSION ****************************
*mess 'Traitement PRESSION';

'SI' (exist rv 'rvpd') ;
  rvpd = rv.'rvpd' ;
  IDigv=rv.'IDigv';
   Digv=rv.'Digv';
  matpr= rvp.'MATP' ;
'FINSI' ;

*/1    CALCUL de CMCT  pour PRESSION
'SI' (testpr et testp1) ;

  rvpd = (rvp.'DOMAINE') ;
  Diago = doma rvpd 'XXDIAGSI' ;
  si ( ega ('VALEUR' 'DIME') 2) ;
  Digv= ( exco Diago 'SCAL' cnvi1 ) et ( exco Diago 'SCAL' cnvi2 ) ;
  sinon ;
  Digv= ( exco Diago 'SCAL' cnvi1 ) et ( exco Diago 'SCAL' cnvi2 ) et
        ( exco Diago 'SCAL' cnvi3 ) ;
  finsi ;
  Digv = kcht rvpd vect sommet comp lc Digv ;

    'SI' ('EXISTE' rv 'CLIM') ;
    rvp . 'CLIM' = rv . 'CLIM' ;
    Digv = 'KOPS' Digv 'CLIM' (rv . 'CLIM') -1;
    'FINSI' ;

  IDigv= 'INVERSE' digv ;

    'SI' ('NON' ('EXISTE' rvp 'DIAGV')) ;
    rvp . 'DIAGV' = digv ;
    'FINSI' ;

  rvp . 'MATC' = 'KMAB' rvp ;
  rvp . 'PRESSION' = 'KCHT' rvpd 'SCAL'
                     'CENTRE' 0.D0 ;
  rvp . 'GRADP'    = 'KCHT' rvpd 'VECT'
                     'SOMMET' vnul ;

  achp  matpr= 'KOPS' 'MATRIK' ;
  rvp . 'MATP' = matpr ;
  rv.'rvpd'=rvpd;
  rv.'IDigv'=IDigv;
  rv.'Digv'=Digv;

'FINSI' ;

*/3   CALCUL de CMCT  pour PROJ
'SI' (testprj et testp1) ;

    Diago mma = 'KOPS' 'MATRIK' ;

      idfdt = 0 ;
    'REPETER' blocj ('DIME' (rv . 'LISTOPER')) ;
      nomper = 'EXTRAIRE' &blocj (rv . 'LISTOPER') ;
      notable= 'MOT' ('TEXTE' ('CHAINE' &blocj nomper)) ;
*     mess 'Operateur-> ' nomper ;
      si ((ega nomper 'DFDT    ') et
       ('EXIST' (rv.notable) 'LISTINCO'));
       si ('EXIST' (rv. notable.'LISTINCO') nomvi);
        idfdt=idfdt + 1 ;
        Diago=Diago et (doma (rv. notable .'DOMZ') 'XXDIAGSI');
        rvpd=rv. notable . 'DOMZ';
       finsi ;
      finsi ;
    'FIN' blocj ;
    si (ega idfdt 0) ; mess ' Pas de DFDT ?? ' ; erreur 21 ; finsi ;

  si ( ega ('VALEUR' 'DIME') 2) ;
  Digv= ( exco Diago 'SCAL' cnvi1 ) et ( exco Diago 'SCAL' cnvi2 ) ;
  sinon ;
  Digv= ( exco Diago 'SCAL' cnvi1 ) et ( exco Diago 'SCAL' cnvi2 ) et
        ( exco Diago 'SCAL' cnvi3 ) ;
  finsi ;
* Digv = kcht rvpd vect sommet comp lc (Digv*Ro)  ;                             
  Digv = kcht rvpd vect sommet comp lc Digv  ;
* Matrice diagonale ne contenant pas les conditions aux limites : Diag
     Diag  = copier Digv;
     IDiag = 'INVERSE' Diag ;

  Digv = 'KOPS' Digv 'CLIM' (rv . 'CLIM') -1;

  IDigv= 'INVERSE' digv ;

    'SI' ('NON' ('EXISTE' rvp 'DIAGV')) ;
    rvp . 'DIAGV' = digv ;
    'FINSI' ;

  rvp.'INCO'=rv.'INCO' ;

    sp map  = 'KOPS' 'MATRIK' ;
    sr mar  = 'KOPS' 'MATRIK' ;
    svnpc mvnpc = 'KOPS' 'MATRIK' ;
    'REPETER' blocpj ('DIME' (rvp . 'LISTOPER')) ;
      nomper = 'EXTRAIRE' &blocpj (rvp . 'LISTOPER') ;
      notable= 'MOT' ('TEXTE' ('CHAINE' &blocpj nomper)) ;
*   mess 'Procedure PROJ Operateur : ' nomper ;

      si (EGA nomper 'KBBT');
      rvp . notable . 'KOPT' . 'IKOMP' = 1 ;
      finsi ;

      msi mai=('TEXTE' nomper) (rvp . notable) ;

*******PROJ*************************************************************
*******PROJ TGRAD  calcul maig  (Formulation en Gradient)       DEBUT*** *******
  Si TGRAD ;
    mess 'Procedure PROJCT (Gradient) Operateur : ' nomper ;
      rvp . notable . 'KOPT' . 'IKOMP' = 0 ;
      msig maig=('TEXTE' nomper) (rvp . notable) ;
  Finsi ;
*******PROJ                                             -> MATG FIN*****
*******PROJ/////////////////////////////////////////////////////////////

      TVNP=FAUX;
      TVNPC=FAUX;
     Si (('EGA' nomper  'VNIMP   ') et (non (ega (rvp.'DISCPRES') 5)));
      TVNP=VRAI ;
      mar = mar 'ET' mai ;
      sr  = sr  'ET' msi ;
     Finsi ;

*     mess 'nomper=' nomper (rvp.'DISCPRES');
     Si (('EGA' nomper  'VNIMP   ') et (ega (rvp.'DISCPRES') 5));
      TVNPC=VRAI;
      mvnpc=mvnpc et mai ;
      svnpc=svnpc et msi ;
     Finsi ;

     Si (non ('EGA' nomper  'VNIMP   '));
      map = map 'ET' mai ;
      sp  = sp  'ET' msi ;
     Finsi ;

    'FIN' blocpj ;

  matpc mac  = kops 'CMCTSPLT' map ;
*******PROJ*************************************************************
*******PROJ TGRAD  (vrai  Formulation en Gradient) CMCTSPLT     DEBUT*** *******
  Si TGRAD ;
  matpcg macg  = kops 'CMCTSPLT' maig ;
  Finsi ;
*******PROJ                                             -> MATG FIN*****
*******PROJ/////////////////////////////////////////////////////////////

    Si TVNPC; mac=mac et mvnpc; Finsi;
  rvp . 'MATC' = mac ;

*******PROJ*************************************************************
*******PROJ TGRAD  (vrai  Formulation en Gradient) -> MATG      DEBUT*** *******
    SI TGRAD ;   rvp . 'MATG' = macg;  FINSI ;
*******PROJ                                             -> MATG FIN*****
*******PROJ/////////////////////////////////////////////////////////////

  scr mcr  = 'KOPS' 'MATRIK' ;
  Si TVNP ;
  rvp . 'MBTR' = mar ;
  Dunit=Idigv;
  crt= kops 'CMCT' mac mar (Dunit)  ;
  ctr= kops 'CMCT' mar mac (Dunit)  ;
  rrt= kops 'CMCT' mar mar (Dunit)  ;
  mcr= mcr et crt et ctr et rrt ;
  Finsi ;
  rvp . 'TVNP' = TVNP ;

  Si( ega (rvp.'DISCPRES') 5 ) ;
  mess ' Cas des pressions continues' ;
  matpr = matpc et mcr ;
  Sinon ;
  mess ' Cas des pressions discontinues' ;

  matpr = kops 'CMCT' mac mac IDigv ;
  matpr = matpc et matpr et mcr ;

  Finsi ;


  rvp . 'MATP' = matpr ;
  rv.'rvpd'=rvpd;
  rv.'IDigv'=IDigv;
  rv.'Digv'=Digv;
'FINSI' ;
*************
** Fin calcul CMCT pour testp1 et testprj
******************************************
*                           t
***************  Calcul de C p   Pression
'SI' testpr ;

      dt = (rv . 'PASDETPS' . 'DELTAT') '*' (rv . 'ALFA') ;
      rvp . 'DELTAT' = dt ;
      f = 'COPIER' s2 ;
      u = rv . 'INCO' . nomvi ;

      lc = 'EXTRAIRE' digv  'COMP' ;
      fu = 'KCHT'     rvpd 'VECT' 'SOMMET'
                      'COMP' lc ('EXCO' f lc) ;

      'SI' ('EXISTE' rv 'CLIM') ;
        dti = -1.D0 '/' dt ;
        dm1f= 'KOPS' (dt '*' ('KOPS'  fu '*' IDigv))
                     'CLIM' (rv . 'CLIM') 3 ;
        dm1f= dti * dm1f ;
      'SINON'  ;
        dm1f= (-1.D0) '*' ('KOPS' fu '*' IDigv)  ;
      'FINSI'  ;

      rvp . 'PRESSION' = 'KMF' (rvp . 'MATC') dm1f ;

      'KRES' rvp (rvp . 'PRESSION')
             'BETA' (rvp . 'KBETA') (rvp . 'BETA')
             'PIMP' (rvp . 'KPIMP') (rvp . 'PIMP') ;
      rvp . 'GRADP' = 'KMTP' 1 (rvp . 'MATC')
                               (rvp . 'PRESSION') lc ;
      s2 = s2 + (rvp . 'GRADP') ;
      rv.'INCO'.'PRESSION'=rvp . 'PRESSION';

'FINSI' ;

*******PROJ*************************************************************
*******PROJ MDM1 Ici on calcule M D-1 pour ensuite calculer Ctp DEBUT*** *******
*******PROJ on ne calcule M D-1 que si (rv 'MDM1') n'existe pas DEBUT*** *******
*******PROJ ou CALPRE VRAI                                      DEBUT*** *******
'SI' testprj ;

      dt = (rv . 'PASDETPS' . 'DELTAT') '*' (rv . 'ALFA') ;
      rvp . 'DELTAT' = dt ;

** Produit M D-1
 'SI' TMDM1                                  ;
  'SI' (('EXIST' rv 'MDM1') et (NON CALPRE)) ;
     MDM1=rv.'MDM1' ;
  'SINON' ;
     mess ' On calcule MD-1 ' ;
     stn matn= 'KOPS' 'MATRIK' ;

        ROW = 0.;
    'REPETER' blocpj1 ('DIME' (rvp . 'LISTOPER')) ;
      nomper1 = 'EXTRAIRE' &blocpj1 (rvp . 'LISTOPER') ;
      notable1= 'MOT' ('TEXTE' ('CHAINE' &blocpj1 nomper1)) ;
*    mess 'Operateur ' nomper1 ;

      si (ega nomper1 'KBBT    ');
      nomiv1= extr  (rvp . notable1 . 'LISTINCO')  1 ;

       'REPETER' blocpj2 ('DIME' (rv . 'LISTOPER')) ;
         nomper2 = 'EXTRAIRE' &blocpj2 (rv . 'LISTOPER') ;
         notable2= 'MOT' ('TEXTE' ('CHAINE' &blocpj2 nomper2)) ;
*        mess 'Operateur ' nomper2 ;

         si (ega nomper2 'DFDT    ');
          nomiv2= extr  (rv . notable2 . 'LISTINCO')  1 ;
           si ('EGA' nomiv2 nomiv1) ;

          Si (EGA TYPROJ 'PENA');
            rwi=(rv . notable2) . ARG1;
            trwi = type rwi ; mess ' type du coeff rwi ' trwi;

            Si (EGA trwi 'FLOTTANT');
            ROW = maxi (prog rwi ROW);
            Finsi ;
            Si ((EGA trwi 'CHPOIN') ou  (EGA trwi 'MCHAML')) ;
            ROW = maxi (prog (maxi rwi) ROW);
            Finsi ;
            Si (EGA (type rwi) 'MOT');
              rwi = rv . 'INCO' . rwi ;
              Si (EGA (type rwi) 'FLOTTANT');
              ROW = maxi (prog rwi ROW);
              Finsi ;
              Si ((EGA trwi 'CHPOIN') ou  (EGA trwi 'MCHAML')) ;
              ROW = maxi (prog (maxi rwi) ROW);
              Finsi ;
            Finsi ;
          Sinon                 ;
              ROW = 1.          ;
          Finsi                 ;
*           mess 'ROW = ' ROW ;

            domzp=rv . notable2 . 'DOMZ' ;
            msi mai='DFDT' (rv . notable2) ;
            matn = matn et mai ;
            sti = kcht domzp vect sommet comp lc vuni;
           finsi ;
         finsi ;
       'FIN' blocpj2 ;
      finsi ;
    'FIN' blocpj1 ;

  Si (EGA TYPROJ 'PENA');
  MDM1=1.               ;
  SINON                 ;
  MDM1='KMF' matn IDiag ;
  FINSI                 ;

    Si(EGA ISCHT 1);
    MDM1=MDM1 * ((dt*2.)/3.) ;
    Sinon ;
    MDM1=MDM1 * dt ;
    Finsi ;

  rv.'MDM1' = MDM1 ;

*******PROJ                                   MDM1 -> rv.'MDM1' FIN*****
*******PROJ/////////////////////////////////////////////////////////////
  'FINSI';
 'FINSI';
*******PROJ*************************************************************
*******PROJ            t  n-1                                   DEBUT*** *******
*******PROJ Calcul de C  P                                      DEBUT*** *******

  TVNP=rvp.'TVNP' ;

  'SI' ('EXIST' (rv.'INCO') 'PRESSION') ;
       PPI = rv.'INCO'.'PRESSION'       ;
     'SI' ('EXIST' (rv.'INCO') 'PNM2')  ;
* TPNM2 elimination of end of step velocity (2 Pn - Pn-1)
       PPI = 2*(PPI) - rv.'INCO'.'PNM2' ;
     'FINSI'                            ;

*Formulation en u grad p
      Si (Exist rvp 'MATG');
      cpre = 'KMF' (rvp . 'MATG') PPI                    'TRAN' ;
*Formulation en p div u
      Sinon;
      cpre = 'KMF' (rvp . 'MATC') PPI                    'TRAN' ;
      Finsi ;

      Si TVNP ;
      cxre = 'KMF' (rvp . 'MBTR') PPI                    'TRAN' ;
      cpre = cpre et cxre ;
      Finsi ;

*******PROJ MDM1
      Si TMDM1 ;
* consistence selon (Gresho)
      gradpres = 'KOPS' MDM1 '*' cpre ;
      Sinon ;
* consistence selon (Guermond)
      gradpres =                 cpre     ;
      Finsi ;

      oublier cpre ;
      rv.'INCO'.'GRADPRES' = 'NOMC' lcu lc gradpres;

*Formulation en u grad p
      Si (Exist rvp 'MATG');
      s2 = sf - (rv.'INCO'.'GRADPRES') + st ;
      Sinon ;
*Formulation en p div u
      s2 = sf + (rv.'INCO'.'GRADPRES') + st     ;
      Finsi ;

  'FINSI' ;

*******PROJ                  t  n-1                             FIN*****
*******PROJ        MDM1 '*' C  P        -> rv.'INCO'.'GRADPRES' FIN*****
*******PROJ                                          t  n-1     FIN*****
*******PROJ         si (p div u)        -> S2 = F + C  P  + st  FIN*****
*******PROJ                                          t  n-1     FIN*****
*******PROJ         si (u grad p)       -> S2 = F - C  P  + st  FIN*****
*******PROJ/////////////////////////////////////////////////////////////
'FINSI' ;

**********  FIN Traitement PRESSION ************************
*   mess ' FIN Traitement PRESSION' ;

    'SI' ('EXISTE' rv 'CLIM') ;
      s1 = rv . 'CLIM' ;
    'SINON' ;
      s1=' ' ;
    'FINSI'  ;
    rv . 'S2' = s2 ;

*******************************************************************
* Résolution hors QDM méthode de projection
 'SI'((NON testprj) ou (EGA mdfdt 0));
*mess '* Résolution hors QDM méthode de projection ';

    rv . 'METHINV' . 'XINIT'  = resmn ;
    
   'SI' ('EXISTE' rv 'GPROJ') ;
      res = GRESP ma1 s1 s2 rv ;
   'SINON' ;      
       res = KRESP  ma1 'TYPI' (rv . 'METHINV')
                        'CLIM' s1
                        'SMBR' s2
                        'IMPR' IMPKRES ;
   'FINSI' ;                         

 'FINSI' ;

'SI' (testprj et (non(EGA mdfdt 0)));
*******PROJ*************************************************************
*******PROJ  Résolution QDM                                     DEBUT*** *******
*******PROJ  (On éclate les résolutions)                        DEBUT*** *******

* Résolution QDM méthode de projection (On éclate les résolutions)
*mess ' Résolution QDM méthode de projection';

     lpart  = KOPS 'EXTRCOUP' ma1  ;
     nbpart = dime lpart           ;

    $ma1  =table 'ESCLAVE';
    $s1   =table 'ESCLAVE';
    $s2   =table 'ESCLAVE';
    $resmn=table 'ESCLAVE';
    $tab1 =table 'ESCLAVE';

   repeter Bclcom nbpart ;
    nmc=lpart.&Bclcom;
*   mess ' Liste des composantes';
*   list nmc ;
    MA1i              = kops 'EXTRINCO' nmc nmc ma1 ;
    $ma1  .&Bclcom    = MA1i                        ;
    $s1   .&Bclcom    = exco s1 nmc 'NOID'          ;
    $s2   .&Bclcom    = exco s2 nmc 'NOID'          ;
    resmni            = exco resmn nmc 'NOID'       ;
    TAB1   = copie rv . 'METHINV'                   ;
*   TAB1   =       rv . 'METHINV'                   ;
    TAB1 . 'XINIT'  = resmni                        ;
    $tab1 .&Bclcom    = TAB1                        ;

*-------------------------------------------------------------------------
*mess ' Mise à jour du préconditionnement';
'SI' ('>EG' TAB1 . 'TYPINV' 2) ; 
  lword1 lword2 = 'KOPS' 'EXTRNINC' MA1i ;
  word1 = 'EXTR' lword1 1 ;

  'SI' ('NON' ('EXIS' rv 'TABRES')) ;
    rv . 'TABRES' = 'TABLE' ;
  'FINSI' ;
  TABRES=rv.'TABRES';

  'SI' (TAB1 . 'CALPREC') ;
    'SI' ('NON' ('EXIS' TABRES word1)) ;
      TABRES . word1 = 'TABLE' ;
      TABRES . word1 . 'MATASS' = MA1i ;
      TABRES . word1 . 'MAPREC' = MA1i ;
      'MESS' 'On recalcule le preconditionneur';
    'SINON' ;
      TABRES . word1 . 'MATASS' = MA1i ;
      TABRES . word1 . 'MAPREC' = MA1i ;
      'MESS' 'On recalcule le preconditionneur';
    'FINSI' ;
  'FINSI' ;
  TAB1 . 'MATASS' = TABRES.word1.'MATASS' ;
  TAB1 . 'MAPREC' = TABRES.word1.'MAPREC' ;
'FINSI' ;
*-------------------------------------------------------------------------
   Fin Bclcom ;
 
Si(EXIST RV 'KPR');
KPR=rv.'KPR';
Sinon;
KPR=FAUX;
Finsi;

Si KPR ;
mess ' On parallélise les résolutions ';
    res = ASSI TOUS KRES   $ma1 'TYPI' $tab1
                                'CLIM' $s1
                                'SMBR' $s2
                                'IMPR' IMPKRES ;
    res=et res ;
Sinon ;
mess ' Les résolutions sont traitées séquentiellement';
    res mm1= kops 'MATRIK'        ;
    repeter Bclcom nbpart ;
    resi= KRES   ($ma1.&bclcom) 'TYPI' ($tab1.&bclcom)
                                'CLIM' ($s1.&bclcom)
                                'SMBR' ($s2.&bclcom)
                                'IMPR' IMPKRES ;
    res=res et resi ;
    Fin Bclcom ;
Finsi;

*mess ' Fin Résolution QDM méthode de projection';
* Fin Résolution QDM méthode de projection
*-------------------------------------------------------------------------------
*******PROJ  Résolution QDM                              -> res FIN*****
*******PROJ  (On éclate les résolutions)                        FIN*****
*******PROJ/////////////////////////////////////////////////////////////
'FINSI' ;

'SI' testprj ;
*******PROJ*************************************************************
*******PROJ  ETAPE DE PROJECTION                                DEBUT*** *******
*******PROJ  on calcul cun (alias c U tilde)                    DEBUT*** *******
*     mess 'ETAPE DE PROJECTION ' ;
*   _n
* C U
      cun = ('KMF' (rvp . 'MATC') res) ;
      Si TVNP ;
      cxn = ('KMF' (rvp . 'MBTR') res) ;
      cun = cun et cxn ;
      Finsi ;

*    mess ' On calcule les seconds membres de l équation de pression ' ;
*    mess ' s ils existent (opérateurs FIMP)                         ' ;

    'REPETER' blocpj ('DIME' (rvp . 'LISTOPER')) ;
      nomper = 'EXTRAIRE' &blocpj (rvp . 'LISTOPER') ;
      notable= 'MOT' ('TEXTE' ('CHAINE' &blocpj nomper)) ;
      si (ega nomper 'FIMP    ');
*   mess 'Seconds membres de l équation de pression Opérateur ' nomper ;
      msi mai=('TEXTE' nomper) (rvp . notable) ;
      cun = cun et msi ;
      finsi ;
    'FIN' blocpj ;

      cun = cun * ROW * (-1./dt) ;

      rvp . 'METHINV' . 'XINIT'  = resmn2  ;

                    sl=  KRESP matpr 'TYPI' (rvp . 'METHINV')
                        'CLIM' (rvp . CLIM )
                        'SMBR' cun
                        'IMPR' IMPKRES  ;  oublier cun ;

      oublier resmn2 ; resmn2 = sl ;
      ctl = 'KMF' (rvp . 'MATC') sl 'TRAN' ;
      Si TVNP ;
      cxl = 'KMF' (rvp . 'MBTR') sl 'TRAN' ;
      ctl = ctl + cxl ;
      Finsi ;

* elimination of end of step velocity (Si exist PNM2) (on ne corrige pas)
     'SI'(NON TPNM2);
     'SI' (NON ('EXIST' (rv.'INCO') 'PNM2')) ;
      a=nomc lcu lc ('KOPS' IDigv '*' ctl) ;
      res  =  res + ((1./ROW)*a*dt) ;
     'FINSI' ;
     'FINSI' ;

      oublier ctl ;

*******PROJ                                -> (alias c U tilde) FIN*****
*******PROJ                                -> sl (alias lambda) FIN*****
*******PROJ                                -> on corrige u      FIN*****
*******PROJ/////////////////////////////////////////////////////////////
'FINSI' ;

*****  Avancement en temps
*?  Si (testprj et ('EGA' mdfdt 0)) ; IMPTCRR=RV.'IMPR' ;  finsi ;
 Si ('EGA' mdfdt 0) ; IMPTCRR=RV.'IMPR' ;  finsi ;
    IMPKRES=0 ;
*   mess 'On passe dans TCRR ' ;
    eps = 'TCRR' res omeg (rv . 'INCO') 'IMPR' IMPTCRR ;

    'OUBLIER' resmn ;
    rv.'METHINV'.'CALPREC' = FAUX ;
    'SI' (testprj) ;
      rvp.'METHINV'.'CALPREC'= FAUX ;
    'FINSI' ;
    resmn = res ;
    'MENAGE' ;
*
   'SI' (rv . 'STOPITER') ;
      rv . 'STOPITER' = FAUX ;
      'QUITTER' bloci ;
   'FINSI' ;
  'FIN' bloci ;
  irt=0 ;

* mess ' ITMA= ' (rv . 'ITMA') ' mdfdt= ' mdfdt ;
  'SI' ('EGA' (rv . 'ITMA') 0) ;
*   mess ' ON PASSE PLUS DS TCNM>>>>>>>>>>>>>>' ;
     irt = 'TCNM' rv 'NOUP';
  'SINON'                                         ;
*   mess ' ON PASSE DS TCNM>>>>' (rv . 'ITMA') mdfdt;
     irt = 'TCNM' rv ;
  'FINSI' ;

*****  Avancement en temps  algorithme de projection
*******PROJ*************************************************************
*******PROJ  P(n+1) = P(n) + sl                                 DEBUT*** *******
* Algorithme de projection  P(n+1) = P(n) + sl
   'SI' testprj ;

     'SI' (NON ('EXIST' (rv.'INCO') 'PRESSION'));
         rv.'INCO'.'PRESSION' = sl;
         Si TPNM2 ; rv.'INCO'.'PNM2' = sl ; Finsi ;
     'SINON' ;
         PNM1=rv.'INCO'.'PRESSION';
        si ( EGA ISCHT 1) ;
         PN=PNM1 + (1.5*sl) ;
        sinon ;
         PN=PNM1 + sl ;
        finsi ;

       Si TPNM2 ; rv.'INCO'.'PNM2'=PNM1 ;  Finsi ;
       rv.'INCO'.'PRESSION'=PN   ;
     'FINSI' ;

     Si ((EGA TYPROJ 'PSCT') ou (EGA TYPROJ 'PENA')) ;
     rv.'INCO'.'PRESSION' = sl;
     Finsi ;

       'OUBLIER' sl ;
   'FINSI' ;
*******PROJ  P(n+1) = P(n) + sl                                 FIN*****
*******PROJ/////////////////////////////////////////////////////////////
*****  Avancement en temps  FIN


  'SI' testran ;
    rv . 'CO' . 'VITESSE' ='KOPS' (rv . 'INCO' . nomvi)
                             '-'  (rv . 'SEDIM') ;
    k = 'ABS' (rv . 'INCO' . 'KN') ; 
    e = 'ABS' (rv . 'INCO' . 'EN') ;
    k = 'KOPS' ('KOPS' k '*' k) '/' e ;
    dif = 'KOPS' ('KOPS' k '*' 0.09) '+' (rv . 'COEF') ;
    rv . 'CO' . 'DIFFU'   = 'NOEL' rvpd dif ;
    rv . 'CO' . 'TEMPERA' = rv . 'INCO' . 'CN' ;
  'FINSI' ;
  rv.'METHINV'.'CALPREC' = FAUX ;
  'SI' (testprj) ;
    rvp.'METHINV'.'CALPREC'= FAUX ;
  'FINSI' ;
  'MENAGE' ;

  'SI' ('EGA' irt 1) ;
    'MESSAGE' ' Temps final atteint : '
              (rv . 'PASDETPS' . 'TPS') ;
    'QUITTER' bloc1 ;
  'FINSI' ;
  'SI' (rv . 'STOPPDT') ;
    rv . 'STOPPDT' = FAUX ;
    'MESSAGE' ' Arret demandé au temps : ' (rv . 'PASDETPS' . 'TPS') ;
    'QUITTER' bloc1 ;
  'FINSI' ;
'FIN' bloc1 ;

'SI' testpr ;
  rvp . 'PRESSION' = 'KCHT' rvpd 'SCAL' 'CENTRE'
                            (rvp . 'PRESSION') ;
  rvp . 'PN'       = 'ELNO' rvpd
                            (rvp . 'PRESSION') ;
'FINSI' ;
************************    E X E C ************************************
'FINPROC' ;