Télécharger kres4.eso

Retour à la liste

Numérotation des lignes :

  1. C KRES4 SOURCE PV 16/11/17 22:00:24 9180
  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 (IRET.NE.0) GOTO 9999
  126. C
  127. C Obtention de la solution (montée-descente)
  128. C
  129. CALL REZOLU(IDMAT,
  130. $ KS2B,
  131. $ IMPR,IRET)
  132. IF (IRET.NE.0) GOTO 9999
  133. INCX=KS2B
  134. IF (ISCAL.GT.0) THEN
  135. C
  136. C On dénorme le vecteur solution : attention modification...
  137. C
  138. CALL NORMV2(INCX,NORMP,IMPR,IRET)
  139. IF (IRET.NE.0) GOTO 9999
  140. ENDIF
  141. C
  142. C Transformation du vecteur-solution en chpoint
  143. C
  144. CALL XDISP(MATRIK,INCX,MCHSOL,IMPR,IRET)
  145. IF (IRET.NE.0) GOTO 9999
  146. C
  147. C Suppression des objets temporaires
  148. C
  149. SEGSUP,AMORS
  150. SEGSUP,AISA
  151. SEGSUP INCX
  152. *
  153. * Normal termination
  154. *
  155. IRET=0
  156. RETURN
  157. *
  158. * Error handling
  159. *
  160. 9999 CONTINUE
  161. IRET=1
  162. WRITE(IOIMP,*) 'An error was detected in kres4.eso'
  163. RETURN
  164. *
  165. * End of KRES4
  166. *
  167. END
  168.  
  169.  
  170.  
  171.  
  172.  
  173.  
  174.  
  175.  
  176.  
  177.  
  178.  
  179.  
  180.  
  181.  
  182.  
  183.  
  184.  
  185.  
  186.  

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