Télécharger prdilu.eso

Retour à la liste

Numérotation des lignes :

prdilu
  1. C PRDILU SOURCE PV 20/09/26 21:19:20 10724
  2. SUBROUTINE PRDILU(KMORS,KISA,MATRIK,IMPR,IRET)
  3. IMPLICIT INTEGER(I-N)
  4. IMPLICIT REAL*8 (A-H,O-Z)
  5. C***********************************************************************
  6. C NOM : PRDILU
  7. C DESCRIPTION :
  8. C Calcul du préconditionneur D-ILU d'une matrice Morse.
  9. C D-ILU : Diagonal Incomplete LU factorization
  10. C
  11. C On stocke l'inverse des pivots du préconditionneur
  12. C dans un segment de type
  13. C IZA pointé par IDIAG du segment IDMAT
  14. C pointé par KIDMAT(1) du segment MATRIK.
  15. C (Toujours la réutilisation de l'existant...)
  16. C Normalement, aucun indice strictement nul n'est stocké
  17. C (appel de clmors dans melim) donc il suffirait de
  18. C tester l'existence de l'indice diagonal
  19. C On met toutefois un warning si l'indice est très petit.
  20. C
  21. C Ce sous-programme est en fait une interface à :
  22. C mdilus.eso : dans le cas symétrique
  23. C mdilun.eso : dans le cas non-symétrique
  24. C qui sont en Fortran presque pur (pour raison de rapidité)
  25. C et effectuent la construction proprement dite du
  26. C préconditionneur.
  27. C
  28. C ATTENTION : pour une matrice A quelconque, la factorisation
  29. C --------- D-ILU peut ne pas exister (pivot nul) ou avoir
  30. C des pivots négatifs MEME SI la factorisation
  31. C complète de A existe et n'a que des pivots
  32. C positifs.
  33. C
  34. C
  35. C LANGAGE : ESOPE
  36. C AUTEUR : Stéphane GOUNAND (CEA/DRN/DMT/SEMT/TTMF)
  37. C mél : gounand@semt2.smts.cea.fr
  38. C REFERENCE (bibtex-like) :
  39. C @BOOK{templates,
  40. C AUTHOR={R.Barrett, M.Berry, T.F.Chan, J.Demmel, J.Donato,
  41. C J.Dongarra, V.Eijkhout, R.Pozo, C.Romine,
  42. C H. Van der Vorst},
  43. C TITLE={Templates for the Solution of Linear Systems :
  44. C Building Blocks for Iterative Methods},
  45. C PUBLISHER={SIAM}, YEAR={1994}, ADDRESS={Philadelphia,PA} }
  46. C -> URL : http://www.netlib.org/templates/Templates.html
  47. C***********************************************************************
  48. C APPELES : -
  49. C***********************************************************************
  50. C ENTREES : MATRIK, IMPR
  51. C ENTREES/SORTIES : -
  52. C SORTIES : INVPIV (IDIAG dans KIDMAT(1) dans MATRIK)
  53. C IRET
  54. C CODE RETOUR (IRET) : 0 si ok
  55. C <0 si problème
  56. C MATRIK : pointeur sur segment MATRIK de l'include SMMATRIK
  57. C on pioche dedans les informations nécessaires
  58. C (différents pointeurs, nb. de ddl...)
  59. C IMPR : niveau d'impression
  60. C INVPIV : pointeur sur segment IZA de l'include SMMATRIK
  61. C vecteur contenant l'inverse de la factorisation
  62. C D-ILU de la matrice morse pointée par MATRIK
  63. C (KIDMAT(4-5))
  64. C***********************************************************************
  65. C VERSION : v1, 01/04/98, version initiale
  66. C HISTORIQUE : v1, 01/04/98, création
  67. C HISTORIQUE : 09/02/98, on ne construit pas le préconditionneur s'il
  68. C existe déjà.
  69. C HISTORIQUE :
  70. C***********************************************************************
  71. C Prière de PRENDRE LE TEMPS de compléter les commentaires
  72. C en cas de modification de ce sous-programme afin de faciliter
  73. C la maintenance !
  74. C***********************************************************************
  75.  
  76. -INC PPARAM
  77. -INC CCOPTIO
  78. POINTEUR KMORS.PMORS
  79. POINTEUR KISA.IZA
  80. POINTEUR INVPIV.IZA
  81. *
  82. * .. Variables locales
  83. C***
  84. IRET=0
  85. IF (IMPR.GT.5) WRITE(IOIMP,*) 'Entrée dans prdilu'
  86. C On récupère les segments utiles
  87. SEGACT MATRIK
  88. NTTT =KNTTT
  89. IDMAT=KIDMAT(1)
  90. KSIM =KSYM
  91. SEGACT IDMAT
  92. INVPIV=IDIAG
  93. SEGDES IDMAT
  94. SEGDES MATRIK
  95. C Le préconditionneur est-il déjà construit ?
  96. C!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
  97. C En fait, on surcharge tout le temps INVPIV car IDIAG peut ne pas
  98. C etre nul mais contenir autre chose que le préconditionneur
  99. C!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
  100. INVPIV=0
  101. C!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
  102. IF (INVPIV.EQ.0) THEN
  103. C
  104. C On parcourt les lignes de la matrice
  105. C comme il sied au stockage Morse
  106. C
  107. NBVA=NTTT
  108. SEGINI INVPIV
  109. SEGACT KMORS
  110. SEGACT KISA
  111. NNZ=KISA.A(/1)
  112. IF (KSIM.EQ.0) THEN
  113. CALL MDILUS(NTTT,NNZ,
  114. $ KMORS.IA,KMORS.JA,KISA.A,
  115. $ INVPIV.A,
  116. $ IRET)
  117. IF (IRET.NE.0) GOTO 9999
  118.  
  119. ELSEIF (KSIM.EQ.2) THEN
  120. CALL MDILUN(NTTT,NNZ,
  121. $ KMORS.IA,KMORS.JA,KISA.A,
  122. $ INVPIV.A,
  123. $ IRET)
  124. IF (IRET.NE.0) GOTO 9999
  125. ELSE
  126. WRITE(IOIMP,*) 'Valeur incorrecte KSYM=',KSIM
  127. WRITE(IOIMP,*) 'dans le MATRIK',MATRIK
  128. IRET=-2
  129. GOTO 9999
  130. ENDIF
  131.  
  132. IF (IMPR.GT.6) THEN
  133. WRITE(IOIMP,*) 'création du pointeur INVPIV='
  134. $ ,INVPIV
  135. IF (IMPR.GT.7) THEN
  136. WRITE(IOIMP,*) 'INVPIV(1..',NBVA,')= '
  137. WRITE(IOIMP,1002)(INVPIV.A(II),II=1,NBVA)
  138. ENDIF
  139. ENDIF
  140.  
  141. SEGDES KISA
  142. SEGDES KMORS
  143. SEGDES INVPIV
  144. C
  145. C On stocke l'inverse de la diagonale obtenue
  146. C
  147. SEGACT IDMAT*MOD
  148. IDIAG=INVPIV
  149. SEGDES IDMAT
  150. ELSE
  151. IF (IMPR.GT.6) THEN
  152. WRITE(IOIMP,*) 'Le préconditionneur est déjà construit :'
  153. WRITE(IOIMP,*) 'INVPIV=',INVPIV
  154. IF (IMPR.GT.7) THEN
  155. SEGACT INVPIV
  156. WRITE(IOIMP,*) 'INVPIV.A(1..',NBVA,')= '
  157. WRITE(IOIMP,1002)(INVPIV.A(II),II=1,NBVA)
  158. SEGDES INVPIV
  159. ENDIF
  160. ENDIF
  161. ENDIF
  162. *
  163. * Normal termination
  164. *
  165. RETURN
  166. *
  167. * Format handling
  168. *
  169. 1002 FORMAT(10(1X,1PE11.4))
  170. *
  171. * Error handling
  172. *
  173. 9999 CONTINUE
  174. WRITE(IOIMP,*) 'An error was detected in prdilu.eso'
  175. RETURN
  176. *
  177. * End of PRDILU
  178. *
  179. END
  180.  
  181.  
  182.  
  183.  
  184.  
  185.  
  186.  
  187.  
  188.  
  189.  
  190.  
  191.  
  192.  
  193.  
  194.  
  195.  

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