Sources Esope | Cast3M

meidia
C MEIDIA    SOURCE    PV        20/09/26    21:18:46     10724          
      SUBROUTINE MEIDIA(KMORS,KISA,MATRIK,IMPR,IRET)
      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8 (A-H,O-Z)
C***********************************************************************
C NOM         : MEIDIA
C DESCRIPTION :
C     Calcul du préconditionneur "Jacobi" d'une matrice Morse.
C     Son calcul est simple : c'est la diagonale de la matrice !
C
C     On stocke l'inverse de la diagonale de la matrice
C     dans un segment de type
C     IZA pointé par IDIAG du segment IDMAT
C     pointé par KIDMAT(1) du segment MATRIK.
C     (Toujours la réutilisation de l'existant...)
C     Si on rencontre une diagonale nulle dans la matrice Morse,
C     on affiche un Warning et on stocke 1.D0 dans l'inverse
C     du préconditionneur.
C     De meme, on affiche un Warning si on rencontre des indices
C     très "petits"...
C
C
C LANGAGE     : ESOPE
C AUTEUR      : Stéphane GOUNAND (CEA/DRN/DMT/SEMT/TTMF)
C               mél : gounand@semt2.smts.cea.fr
C***********************************************************************
C APPELES : -
C***********************************************************************
C ENTREES            : MATRIK, IMPR, IRET
C ENTREES/SORTIES    : -
C SORTIES            : INVPIV (IDIAG dans KIDMAT(1) dans MATRIK)
C CODE RETOUR (IRET) : 0 si ok
C                     &lt;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     INVPIV : pointeur sur segment IZA de l'include SMMATRIK
C              vecteur contenant l'inverse de la diagonale
C              de la matrice morse pointée par MATRIK (KIDMAT(4-5))
C***********************************************************************
C VERSION    : v1, 01/04/98, version initiale
C HISTORIQUE : v1, 01/04/98, création
C HISTORIQUE : 09/02/98, on ne construit pas le préconditionneur s'il
C              existe déjà.
C HISTORIQUE :
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
      POINTEUR KMORS.PMORS
      POINTEUR KISA.IZA
      POINTEUR INVPIV.IZA
*
*     .. Variables locales
*     .. Parameters
      REAL*8     ONE      ,ZERO
      PARAMETER (ONE=1.0D0,ZERO=0.0D0)
*     ..
C     Nombre de pivots nul
      INTEGER NBPIVN
C     Nombre de pivots petits
      INTEGER NBPIVP
C     Nombre de pivots inférieurs à 0
      INTEGER NBPIVI
      REAL*8 VALPIV
C***
      IRET=0
      IF (IMPR.GT.5) WRITE(IOIMP,*) 'Entrée dans meidia'
      NBPIVP=0
C     On récupère les segments utiles
      SEGACT MATRIK
      NTTT =KNTTT
      IDMAT=KIDMAT(1)
      SEGACT IDMAT
      INVPIV=IDIAG
      SEGDES IDMAT
      SEGDES MATRIK
C     Le préconditionneur est-il déjà construit ?
C!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
C En fait, on surcharge tout le temps INVPIV car IDIAG peut ne pas
C etre nul mais contenir autre chose que le préconditionneur
C!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      INVPIV=0
C!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      IF (INVPIV.EQ.0) THEN
C
C     On parcourt la matrice a la recherche des indices diagonaux
C
         NBVA=NTTT
         SEGINI INVPIV
         SEGACT KMORS
         SEGACT KISA
         DO 1 INTTT=1,NTTT
            IDEB=KMORS.IA(INTTT)
            IFIN=KMORS.IA(INTTT+1)-1
            ICOL=IDEB
            IF (IDEB.LE.IFIN) THEN
C     On cherche le terme Aii
 11            CONTINUE
               IF (KMORS.JA(ICOL).LT.INTTT.AND.ICOL.LT.IFIN) THEN
                  ICOL=ICOL+1
                  GOTO 11
               ENDIF
C     On ne l'a pas trouvé
               IF (KMORS.JA(ICOL).NE.INTTT) THEN
                  WRITE(IOIMP,*) 'diag.',INTTT,'inexistante'
                  WRITE(IOIMP,*) 'le préconditionnement par la diag.'
                  WRITE(IOIMP,*) 'est impossible.'
                  IRET=-1
                  GOTO 9999
               ELSE
C     On l'a trouvé
                  VALPIV=KISA.A(ICOL)
                  IF (VALPIV.LT.XPETIT) THEN
                     NBPIVP=NBPIVP+1
                     VALPIV=ONE
                  ENDIF
                  INVPIV.A(INTTT)=ONE/VALPIV
               ENDIF
            ELSE
               WRITE(IOIMP,*) 'Ligne',INTTT,'vide'
               WRITE(IOIMP,*) 'le préconditionnement par la diag.'
               WRITE(IOIMP,*) 'est impossible.'
               IRET=-3
               GOTO 9999
            ENDIF
 1       CONTINUE
*
*     Warning(s)
*
         IF (NBPIVP.GT.0) THEN
            IF (IMPR.GT.1) THEN
               WRITE(IOIMP,*) 'WARNING !'
               WRITE(IOIMP,*) NBPIVP,' |diag.|&lt;',XPETIT
            ENDIF
         ENDIF
 
         IF (IMPR.GT.6) THEN
            WRITE(IOIMP,*) 'création du pointeur INVPIV=',INVPIV
            IF (IMPR.GT.7) THEN
               WRITE(IOIMP,*) 'INVPIV.A(1..',NBVA,')= '
               WRITE(IOIMP,1002)(INVPIV.A(II),II=1,NBVA)
            ENDIF
         ENDIF
 
         SEGDES KISA
         SEGDES KMORS
         SEGDES INVPIV
C
C     On stocke l'inverse de la diagonale obtenue
C
         SEGACT IDMAT*MOD
         IDIAG=INVPIV
         SEGDES IDMAT
      ELSE
         IF (IMPR.GT.6) THEN
            WRITE(IOIMP,*) 'Le préconditionneur est déjà construit :'
            WRITE(IOIMP,*) 'INVPIV=',INVPIV
            IF (IMPR.GT.7) THEN
               SEGACT INVPIV
               WRITE(IOIMP,*) 'INVPIV.A(1..',NBVA,')= '
               WRITE(IOIMP,1002)(INVPIV.A(II),II=1,NBVA)
               SEGDES INVPIV
            ENDIF
         ENDIF
      ENDIF
*
*     Normal termination
*
      RETURN
*
* Format handling
*
 1002 FORMAT(10(1X,1PE11.4))
*
* Error handling
*
 9999 CONTINUE
      WRITE(IOIMP,*) 'An error was detected in meidia.eso'
      RETURN
*
*     End of MEIDIA
*
      END