C PRILTP    SOURCE    GOUNAND   25/04/30    21:15:32     12258          
      SUBROUTINE PRILTP(AMORS,AISA,MATRIK,XLFIL,XDTOL,XSPIV,
     $     IVARI,
     $     INCX,NORMP,LRACOU,
     $     IMPR,IRET)
      IMPLICIT REAL*8 (A-H,O-Z)
      IMPLICIT INTEGER (I-N)
C***********************************************************************
C NOM         : PRILTP
C DESCRIPTION :
C     Calcul du préconditionneur ILUTP d'une matrice Morse
C     Le préconditionneur est une matrice stockée
C     au format MSR (Modified Sparse Row, stockage de l'inverse de la
C     diagonale) de meme profil que la matrice Morse (format CSR) qu'il
C     préconditionne.
C     Le profil et les valeurs du préconditionneur sont
C     stockés dans KIDMAT(6 et 7) (réutilisation de l'existant).
C
C     Ce sous-programme est en fait une interface à :
C       meiltp
C     qui est en Fortran presque pur (pour raison de rapidité)
C     et effectue la construction proprement dite du
C     préconditionneur.
C
C  ATTENTION : pour une matrice A quelconque, la factorisation
C  ---------   ILUTP peut ne pas exister (pivot nul) ou avoir
C              des pivots négatifs MEME SI la factorisation
C              complète de A existe et n'a que des pivots
C              positifs.
C
C
C LANGAGE     : ESOPE
C AUTEUR      : Stéphane GOUNAND (CEA/DRN/DMT/SEMT/LTMF)
C               mél : gounand@semt2.smts.cea.fr
C REFERENCE (bibtex-like) :
C  @BOOK{templates,
C        AUTHOR={R.Barrett, M.Berry, T.F.Chan, J.Demmel, J.Donato,
C                J.Dongarra, V.Eijkhout, R.Pozo, C.Romine,
C                H. Van der Vorst},
C        TITLE={Templates for the Solution of Linear Systems :
C               Building Blocks for Iterative Methods},
C        PUBLISHER={SIAM}, YEAR={1994}, ADDRESS={Philadelphia,PA} }
C  -> URL : http://www.netlib.org/templates/Templates.html
C  Sparskit : a basic tool kit for sparse matrix computations
C  Version 2 (Youcef Saad)
C  -> URL : http://www.cs.umn.edu/Research/arpa/SPARSKIT/sparskit.html
C REMARQUES :
C On pourrait effectuer l'ordonnancement des colonnes du préconditionneur
C (peut-être que la montée-descente va plus vite après, il faudrait
C  tester)
C***********************************************************************
C APPELES       : MEILUT
C APPELE PAR    : KRES2
C***********************************************************************
C ENTREES            : MATRIK, XLFIL, XDTOL, IMPR
C SORTIES            : ILUM, ILUI (KIDMAT(6-7) dans MATRIK), IRET
C CODE RETOUR (IRET) : 0 si ok
C                     <0 si problème
C     MATRIK : pointeur sur segment MATRIK de l'include SMMATRIK
C              on pioche dedans les informations nécessaires
C              (différents pointeurs, nb. de ddl...)
C     IMPR   : niveau d'impression
C***********************************************************************
C VERSION    : v1, 08/04/2004, version initiale
C HISTORIQUE : v1, 08/04/2004, création
C HISTORIQUE :
C***********************************************************************
C Prière de PRENDRE LE TEMPS de compléter les commentaires
C en cas de modification de ce sous-programme afin de faciliter
C la maintenance !
C***********************************************************************

-INC PPARAM
-INC CCOPTIO
-INC CCREEL
-INC SMMATRIK
      INTEGER NTT,NJA,NBVA
      POINTEUR AMORS.PMORS
      POINTEUR AISA.IZA
      POINTEUR KMORS.PMORS
      POINTEUR KISA.IZA
      POINTEUR ILUM.PMORS
      POINTEUR ILUI.IZA
      POINTEUR IDMATP.IDMAT
      POINTEUR IDMATD.IDMAT
      POINTEUR INCX.IZA
      POINTEUR INCX2.IZA
-INC SMLENTI
      INTEGER JG
      POINTEUR IWORK.MLENTI
      POINTEUR JWORK.MLENTI
      POINTEUR KWORK.MLENTI
      POINTEUR IPERM.MLENTI
      POINTEUR JPERM.MLENTI
-INC SMLREEL
      POINTEUR RWORK.MLREEL
      POINTEUR QWORK.MLREEL
      POINTEUR NORMP.MLREEL
      POINTEUR NORMP2.MLREEL
* -INC SMLLOGI
      SEGMENT MLLOGI
      LOGICAL LOGI(JG)
      ENDSEGMENT
      POINTEUR LWORK.MLLOGI
*
      REAL*8 XLFIL,XDTOL
      INTEGER IMPR,IRET
*
      INTEGER ILFIL,NNZA,NNZLU,NNZMLU,NTTDDL
      REAL*8 XNTLIG
      LOGICAL LRACOU
*
      IF (IMPR.GT.5) WRITE(IOIMP,*) 'Entrée dans priltp.eso'
C     On récupère les segments utiles
      SEGACT MATRIK
      ILUM =KIDMAT(6)
      ILUI =KIDMAT(7)
C     Le préconditionneur est-il déjà construit ?
      IF ((ILUM.EQ.0).OR.(ILUI.EQ.0)) THEN
C     La matrice Morse et son préconditionneur ont le meme
C     profil
*! Attention, on modifie la matrice dans meiltp
         SEGACT AMORS*MOD
         NTTDDL=AMORS.IA(/1)-1
         NNZA=AMORS.JA(/1)
         SEGACT AISA*MOD
*
* Ceci bugge avec ILUTPGO2 !!!!!
*
*
*         XNTLIG=DBLE(NNZA)/DBLE(NTTDDL)
*         ILFIL=INT(XNTLIG*XLFIL/2.D0)
*         NNZMLU=NTTDDL*(2*ILFIL+1)
*
         NNZMLU = INT(DBLE(NNZA)*XLFIL) + 1
*
         NTT=NTTDDL-1
         NJA=NNZMLU+1
         SEGINI ILUM
         NBVA=NNZMLU+1
         SEGINI ILUI
         JG=NTTDDL
         SEGINI IWORK
         JG=NTTDDL
         SEGINI JWORK
         JG=NTTDDL
         SEGINI RWORK
         JG=NTTDDL
         SEGINI IPERM
         JG=NTTDDL
         SEGINI JPERM
*         WRITE(IOIMP,*) 'NTTDDL=',NTTDDL,'NNZA=',NNZA,
*     $        'NNZMLU=',NNZMLU
*         WRITE(IOIMP,*) 'ILFIL=',ILFIL,'XLFIL=',XLFIL,'XDTOL=',XDTOL
*         IF (XSPIV.LT.XZERO.OR.XSPIV.GT.1.D0) THEN
*            WRITE(IOIMP,*) 'Error : XSPIV=',XSPIV
*            GOTO 9999
*         ENDIF
         IF (IVARI.EQ.0) THEN
C     Les boucles sont en Fortran pur
            CALL MEILTP(NTTDDL,NNZA,AISA.A,AMORS.JA,AMORS.IA,
     $           XLFIL,XDTOL,XSPIV,NNZMLU,
     $           ILUI.A,ILUM.JA,ILUM.IA,NNZLU,
     $           IWORK.LECT,JWORK.LECT,RWORK.PROG,
     $           IPERM.LECT,JPERM.LECT,
     $           IMPR,IRET)
            IF (IRET.NE.0) GOTO 9999
         ELSEIF (IVARI.EQ.1) THEN
            JG=NTTDDL
            SEGINI LWORK
            CALL MILTP2(NTTDDL,NNZA,AISA.A,AMORS.JA,AMORS.IA,
     $           XLFIL,XDTOL,XSPIV,NNZMLU,
     $           ILUI.A,ILUM.JA,ILUM.IA,NNZLU,
     $           IWORK.LECT,JWORK.LECT,RWORK.PROG,LWORK.LOGI,
     $           IPERM.LECT,JPERM.LECT,
     $           IMPR,IRET)
            IF (IRET.NE.0) GOTO 9999
            SEGSUP LWORK
         ELSEIF (IVARI.EQ.2) THEN
            JG=NTTDDL
            SEGINI LWORK
            JG=NTTDDL
            SEGINI KWORK
            JG=NTTDDL
            SEGINI QWORK
            CALL MILTP3(NTTDDL,NNZA,AISA.A,AMORS.JA,AMORS.IA,
     $           XLFIL,XDTOL,XSPIV,NNZMLU,
     $           ILUI.A,ILUM.JA,ILUM.IA,NNZLU,
     $           IWORK.LECT,JWORK.LECT,RWORK.PROG,LWORK.LOGI,
     $           KWORK.LECT,QWORK.PROG,
     $           IPERM.LECT,JPERM.LECT,
     $           IMPR,IRET)
            IF (IRET.NE.0) GOTO 9999
            SEGSUP QWORK
            SEGSUP KWORK
            SEGSUP LWORK
         ELSE
            WRITE(IOIMP,*) 'IVARI=',IVARI,' non prevu'
         ENDIF
         IF (IMPR.GT.1) THEN
            WRITE(IOIMP,*) 'Préconditionneur : nb.termesstockés=',
     $           NNZLU-1
         ENDIF
*
*         SEGPRT,JPERM
         SEGSUP IPERM
         SEGSUP RWORK
         SEGSUP JWORK
         SEGSUP IWORK
         NBVA=NNZLU
         SEGADJ,ILUI
         NTT=NTTDDL-1
         NJA=NNZLU
         SEGADJ,ILUM
*         SEGPRT,ILUM
*         SEGPRT,ILUI
*
* On s'occupe des permutations
*
         SEGACT MATRIK*MOD
         IF (.NOT.LRACOU) THEN
            KMORS=MATRIK.KIDMAT(4)
            SEGACT KMORS*MOD
            DO IJA=KMORS.IA(1),KMORS.IA(NTTDDL+1)-1
               KMORS.JA(IJA)=JPERM.LECT(KMORS.JA(IJA))
            ENDDO
         ENDIF
         IDMATP=MATRIK.KIDMAT(1)
         IDMATD=MATRIK.KIDMAT(2)
         IF (IDMATD.EQ.0) THEN
            WRITE(IOIMP,*) 'Erreur hyper grave'
         ENDIF
         IF (IDMATP.EQ.IDMATD) THEN
            SEGINI,IDMATP=IDMATD
         ELSE
            SEGACT IDMATP*MOD
         ENDIF
* Copie de NUAN dans NUNA
         DO ITTDDL=1,NTTDDL
            IDMATP.NUNA(ITTDDL)=IDMATP.NUAN(ITTDDL)
         ENDDO
* Mise à jour de NUAN
         DO ITTDDL=1,NTTDDL
            IDMATP.NUAN(ITTDDL)=JPERM.LECT(IDMATP.NUNA(ITTDDL))
         ENDDO
* Mise à jour de NUNA
         DO ITTDDL=1,NTTDDL
            IDMATP.NUNA(IDMATP.NUAN(ITTDDL))=ITTDDL
         ENDDO
         INCX2=INCX
         SEGACT INCX2
         NBVA=NTTDDL
         SEGINI,INCX
         DO ITTDDL=1,NTTDDL
            INCX.A(JPERM.LECT(ITTDDL))=INCX2.A(ITTDDL)
         ENDDO
         SEGSUP INCX2
         IF (NORMP.NE.0) THEN
            NORMP2=NORMP
            SEGACT NORMP2
            JG=NTTDDL
            SEGINI,NORMP
            DO ITTDDL=1,NTTDDL
               NORMP.PROG(JPERM.LECT(ITTDDL))=NORMP2.PROG(ITTDDL)
            ENDDO
            SEGSUP NORMP2
         ENDIF
         SEGSUP JPERM
         MATRIK.KIDMAT(1)=IDMATP
C
C     On stocke la factorisation obtenue du préconditionneur
C
         SEGACT MATRIK*MOD
         KIDMAT(6)=ILUM
         KIDMAT(7)=ILUI
C
         IF (IMPR.GT.6) THEN
            WRITE(IOIMP,*) 'création du préconditionneur Morse',
     $           ' de pointeurs',ILUM,'et',ILUI
            IF (IMPR.GT.8) THEN
               CALL ECMORS(ILUM,ILUI,(IMPR-1))
            ENDIF
         ENDIF
      ELSE
         IF (IMPR.GT.1) THEN
            WRITE(IOIMP,*)
     $           'Préconditionneur deja construit.'
         ENDIF
         IF (IMPR.GT.6) THEN
            WRITE(IOIMP,*) 'Le préconditionneur est déjà construit :',
     $           'ILUM=',ILUM,' et ILUI=',ILUI
            IF (IMPR.GT.8) THEN
               CALL ECMORS(ILUM,ILUI,(IMPR-1))
            ENDIF
         ENDIF
      ENDIF
*
*     Normal termination
*
      RETURN
*
* Format handling
*
*
* Error handling
*
 9999 CONTINUE
      WRITE(IOIMP,*) 'An error was detected in priltp.eso'
      IRET=1
      RETURN
*
*     End of PRILTP
*
      END