C PRDILU SOURCE PV 20/09/26 21:19:20 10724 SUBROUTINE PRDILU(KMORS,KISA,MATRIK,IMPR,IRET) IMPLICIT INTEGER(I-N) IMPLICIT REAL*8 (A-H,O-Z) C*********************************************************************** C NOM : PRDILU C DESCRIPTION : C Calcul du préconditionneur D-ILU d'une matrice Morse. C D-ILU : Diagonal Incomplete LU factorization C C On stocke l'inverse des pivots du préconditionneur 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 Normalement, aucun indice strictement nul n'est stocké C (appel de clmors dans melim) donc il suffirait de C tester l'existence de l'indice diagonal C On met toutefois un warning si l'indice est très petit. C C Ce sous-programme est en fait une interface à : C mdilus.eso : dans le cas symétrique C mdilun.eso : dans le cas non-symétrique C qui sont en Fortran presque pur (pour raison de rapidité) C et effectuent la construction proprement dite du C préconditionneur. C C ATTENTION : pour une matrice A quelconque, la factorisation C --------- D-ILU 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/TTMF) 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*********************************************************************** C APPELES : - C*********************************************************************** C ENTREES : MATRIK, IMPR C ENTREES/SORTIES : - C SORTIES : INVPIV (IDIAG dans KIDMAT(1) dans MATRIK) C 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 INVPIV : pointeur sur segment IZA de l'include SMMATRIK C vecteur contenant l'inverse de la factorisation C D-ILU de la matrice morse pointée par MATRIK C (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*********************************************************************** 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 SMMATRIK POINTEUR KMORS.PMORS POINTEUR KISA.IZA POINTEUR INVPIV.IZA * * .. Variables locales C*** IRET=0 IF (IMPR.GT.5) WRITE(IOIMP,*) 'Entrée dans prdilu' C On récupère les segments utiles SEGACT MATRIK NTTT =KNTTT IDMAT=KIDMAT(1) KSIM =KSYM 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 les lignes de la matrice C comme il sied au stockage Morse C NBVA=NTTT SEGINI INVPIV SEGACT KMORS SEGACT KISA NNZ=KISA.A(/1) IF (KSIM.EQ.0) THEN CALL MDILUS(NTTT,NNZ, $ KMORS.IA,KMORS.JA,KISA.A, $ INVPIV.A, $ IRET) IF (IRET.NE.0) GOTO 9999 ELSEIF (KSIM.EQ.2) THEN CALL MDILUN(NTTT,NNZ, $ KMORS.IA,KMORS.JA,KISA.A, $ INVPIV.A, $ IRET) IF (IRET.NE.0) GOTO 9999 ELSE WRITE(IOIMP,*) 'Valeur incorrecte KSYM=',KSIM WRITE(IOIMP,*) 'dans le MATRIK',MATRIK IRET=-2 GOTO 9999 ENDIF IF (IMPR.GT.6) THEN WRITE(IOIMP,*) 'création du pointeur INVPIV=' $ ,INVPIV IF (IMPR.GT.7) THEN WRITE(IOIMP,*) 'INVPIV(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 prdilu.eso' RETURN * * End of PRDILU * END