Télécharger kres4.eso

Retour à la liste

Numérotation des lignes :

kres4
  1. C KRES4 SOURCE GOUNAND 22/08/25 21:15:06 11434
  2. SUBROUTINE KRES4(MATRIK,KCLIM,KSMBR,
  3. $ ISCAL,
  4. $ MCHSOL,
  5. $ IMPR,IRET)
  6. IMPLICIT REAL*8 (A-H,O-Z)
  7. IMPLICIT INTEGER (I-N)
  8. C***********************************************************************
  9. C NOM : KRES4
  10. C DESCRIPTION : Résolution d'un système par une méthode directe
  11. C (Factorisation LDU sans pivoting).
  12. C
  13. C
  14. C LANGAGE : ESOPE
  15. C AUTEUR : Stéphane GOUNAND (CEA/DRN/DMT/SEMT/LTMF)
  16. C mél : gounand@semt2.smts.cea.fr
  17. C***********************************************************************
  18. C APPELES : VERMAT, MESMBR, MELIM, CLMORS, TRIALU, REZOLU,
  19. C XDISP
  20. C APPELES (E/S) : INFMAT
  21. C APPELE PAR : KRES2
  22. C***********************************************************************
  23. C ENTREES : MATRIK, KCLIM, KSMBR
  24. C ENTREES/SORTIES : -
  25. C SORTIES : MCHSOL
  26. C CODE RETOUR (IRET) : = 0 si tout s'est bien passé
  27. C***********************************************************************
  28. C VERSION : v1, 14/04/2000, version initiale
  29. C HISTORIQUE : v1, 14/04/2000, création
  30. C HISTORIQUE : 06/04/04 : Scaling
  31. C HISTORIQUE :
  32. C HISTORIQUE :
  33. C***********************************************************************
  34. C Prière de PRENDRE LE TEMPS de compléter les commentaires
  35. C en cas de modification de ce sous-programme afin de faciliter
  36. C la maintenance !
  37. C***********************************************************************
  38.  
  39. -INC PPARAM
  40. -INC CCOPTIO
  41. -INC SMCHPOI
  42. POINTEUR KCLIM.MCHPOI
  43. POINTEUR KSMBR.MCHPOI
  44. POINTEUR MCHSOL.MCHPOI
  45. POINTEUR INCX.IZA
  46. POINTEUR KS2B.IZA
  47. POINTEUR KMORS.PMORS
  48. POINTEUR KISA.IZA
  49. POINTEUR AMORS.PMORS
  50. POINTEUR AISA.IZA
  51. POINTEUR NORMP.IZA
  52. POINTEUR NORMD.IZA
  53. *
  54. INTEGER IMPR,IRET
  55. *
  56. * Executable statements
  57. *
  58. IF (IMPR.GT.5) WRITE(IOIMP,*) 'Entrée dans kres4.eso'
  59. C
  60. C On vérifie que la matrice est correctement assemblée
  61. C
  62. CALL VERMAT(MATRIK,IMPR,IRET)
  63. IF (IRET.NE.0) GOTO 9999
  64. C Pas besoin d'estimation de départ pour obtenir la solution
  65. INCX=0
  66. C
  67. C On transforme le chpoint second membre en vecteur second membre
  68. C In MESMBR : SEGINI KS2B
  69. CALL MESMBR(MATRIK,KSMBR,
  70. $ KS2B,
  71. $ IMPR,IRET)
  72. IF (IRET.NE.0) GOTO 9999
  73. C
  74. C On applique les conditions aux limites
  75. C
  76. C In MELIM : SEGINI AMORS
  77. C SEGINI AISA
  78. CALL MELIM(MATRIK,KCLIM,
  79. $ INCX,KS2B,
  80. $ AMORS,AISA,
  81. $ IMPR,IRET)
  82. IF (IRET.NE.0) GOTO 9999
  83. C Eventuellement, on élimine les 0.d0 de la matrice
  84. C On ne le fait que pour les méthodes directes (où le gain en temps est
  85. C important) .
  86. C Il faudrait régler les problèmes de précision pour les termes
  87. C proches de zéro pouvant engendrer des pertes de symétrie
  88. C de la matrice...
  89. C Il y a aussi des problèmes lorsqu'un
  90. C terme vaut zéro à une itération et autre chose aux suivantes.
  91. C L'idéal serait de ne virer que les 0.D0 issus des conditions
  92. C aux limites...
  93. CALL CLMORS(AMORS,AISA,IMPR,IRET)
  94. IF (IRET.NE.0) GOTO 9999
  95. IF (ISCAL.GT.0) THEN
  96. C
  97. C Calcul des normes primales (colonnes) et duales (lignes)
  98. C de la matrice. Norme = norme L2, soit :
  99. C {\sum_{i ou j} a_{ij}^2}^{1/2}
  100. C
  101. CALL NORMAT(AMORS,AISA,ISCAL,NORMP,NORMD,IMPR,IRET)
  102. IF (IRET.NE.0) GOTO 9999
  103. C
  104. C On norme la matrice : attention modification...
  105. C
  106. CALL NORMAM(AMORS,AISA,NORMP,NORMD,IMPR,IRET)
  107. IF (IRET.NE.0) GOTO 9999
  108. C
  109. C On norme le second membre : attention modification...
  110. C
  111. CALL NORMV2(KS2B,NORMD,IMPR,IRET)
  112. IF (IRET.NE.0) GOTO 9999
  113. ENDIF
  114. C
  115. C On donne des infos sur la matrice
  116. C
  117. CALL INFMAT(AMORS,AISA,IMPR,IRET)
  118. IF (IRET.NE.0) GOTO 9999
  119. C
  120. C Factorisation LDU de la matrice
  121. C
  122. CALL TRIALU(MATRIK,AMORS,AISA,
  123. $ IDMAT,
  124. $ IMPR,IRET)
  125. if (ierr.ne.0) return
  126. IF (IRET.NE.0) GOTO 9999
  127. C
  128. C Obtention de la solution (montée-descente)
  129. C
  130. CALL REZOLU(IDMAT,
  131. $ KS2B,
  132. $ IMPR,IRET)
  133. if (ierr.ne.0) return
  134. IF (IRET.NE.0) GOTO 9999
  135. INCX=KS2B
  136. IF (ISCAL.GT.0) THEN
  137. C
  138. C On dénorme le vecteur solution : attention modification...
  139. C
  140. CALL NORMV2(INCX,NORMP,IMPR,IRET)
  141. IF (IRET.NE.0) GOTO 9999
  142. ENDIF
  143. C
  144. C Transformation du vecteur-solution en chpoint
  145. C
  146. CALL XDISP(MATRIK,INCX,MCHSOL,IMPR,IRET)
  147. IF (IRET.NE.0) GOTO 9999
  148. C
  149. C Suppression des objets temporaires
  150. C
  151. SEGSUP,AMORS
  152. SEGSUP,AISA
  153. SEGSUP INCX
  154. *
  155. * Normal termination
  156. *
  157. IRET=0
  158. RETURN
  159. *
  160. * Error handling
  161. *
  162. 9999 CONTINUE
  163. IRET=1
  164. WRITE(IOIMP,*) 'An error was detected in kres4.eso'
  165. RETURN
  166. *
  167. * End of KRES4
  168. *
  169. END
  170.  
  171.  

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