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. -INC CCOPTIO
  39. -INC SMCHPOI
  40. POINTEUR KCLIM.MCHPOI
  41. POINTEUR KSMBR.MCHPOI
  42. POINTEUR MCHSOL.MCHPOI
  43. POINTEUR INCX.IZA
  44. POINTEUR KS2B.IZA
  45. POINTEUR KMORS.PMORS
  46. POINTEUR KISA.IZA
  47. POINTEUR AMORS.PMORS
  48. POINTEUR AISA.IZA
  49. POINTEUR NORMP.IZA
  50. POINTEUR NORMD.IZA
  51. *
  52. INTEGER IMPR,IRET
  53. *
  54. * Executable statements
  55. *
  56. IF (IMPR.GT.5) WRITE(IOIMP,*) 'Entrée dans kres4.eso'
  57. C
  58. C On vérifie que la matrice est correctement assemblée
  59. C
  60. CALL VERMAT(MATRIK,IMPR,IRET)
  61. IF (IRET.NE.0) GOTO 9999
  62. C Pas besoin d'estimation de départ pour obtenir la solution
  63. INCX=0
  64. C
  65. C On transforme le chpoint second membre en vecteur second membre
  66. C In MESMBR : SEGINI KS2B
  67. CALL MESMBR(MATRIK,KSMBR,
  68. $ KS2B,
  69. $ IMPR,IRET)
  70. IF (IRET.NE.0) GOTO 9999
  71. C
  72. C On applique les conditions aux limites
  73. C
  74. C In MELIM : SEGINI AMORS
  75. C SEGINI AISA
  76. CALL MELIM(MATRIK,KCLIM,
  77. $ INCX,KS2B,
  78. $ AMORS,AISA,
  79. $ IMPR,IRET)
  80. IF (IRET.NE.0) GOTO 9999
  81. C Eventuellement, on élimine les 0.d0 de la matrice
  82. C On ne le fait que pour les méthodes directes (où le gain en temps est
  83. C important) .
  84. C Il faudrait régler les problèmes de précision pour les termes
  85. C proches de zéro pouvant engendrer des pertes de symétrie
  86. C de la matrice...
  87. C Il y a aussi des problèmes lorsqu'un
  88. C terme vaut zéro à une itération et autre chose aux suivantes.
  89. C L'idéal serait de ne virer que les 0.D0 issus des conditions
  90. C aux limites...
  91. CALL CLMORS(AMORS,AISA,IMPR,IRET)
  92. IF (IRET.NE.0) GOTO 9999
  93. IF (ISCAL.GT.0) THEN
  94. C
  95. C Calcul des normes primales (colonnes) et duales (lignes)
  96. C de la matrice. Norme = norme L2, soit :
  97. C {\sum_{i ou j} a_{ij}^2}^{1/2}
  98. C
  99. CALL NORMAT(AMORS,AISA,ISCAL,NORMP,NORMD,IMPR,IRET)
  100. IF (IRET.NE.0) GOTO 9999
  101. C
  102. C On norme la matrice : attention modification...
  103. C
  104. CALL NORMAM(AMORS,AISA,NORMP,NORMD,IMPR,IRET)
  105. IF (IRET.NE.0) GOTO 9999
  106. C
  107. C On norme le second membre : attention modification...
  108. C
  109. CALL NORMV2(KS2B,NORMD,IMPR,IRET)
  110. IF (IRET.NE.0) GOTO 9999
  111. ENDIF
  112. C
  113. C On donne des infos sur la matrice
  114. C
  115. CALL INFMAT(AMORS,AISA,IMPR,IRET)
  116. IF (IRET.NE.0) GOTO 9999
  117. C
  118. C Factorisation LDU de la matrice
  119. C
  120. CALL TRIALU(MATRIK,AMORS,AISA,
  121. $ IDMAT,
  122. $ IMPR,IRET)
  123. IF (IRET.NE.0) GOTO 9999
  124. C
  125. C Obtention de la solution (montée-descente)
  126. C
  127. CALL REZOLU(IDMAT,
  128. $ KS2B,
  129. $ IMPR,IRET)
  130. IF (IRET.NE.0) GOTO 9999
  131. INCX=KS2B
  132. IF (ISCAL.GT.0) THEN
  133. C
  134. C On dénorme le vecteur solution : attention modification...
  135. C
  136. CALL NORMV2(INCX,NORMP,IMPR,IRET)
  137. IF (IRET.NE.0) GOTO 9999
  138. ENDIF
  139. C
  140. C Transformation du vecteur-solution en chpoint
  141. C
  142. CALL XDISP(MATRIK,INCX,MCHSOL,IMPR,IRET)
  143. IF (IRET.NE.0) GOTO 9999
  144. C
  145. C Suppression des objets temporaires
  146. C
  147. SEGSUP,AMORS
  148. SEGSUP,AISA
  149. SEGSUP INCX
  150. *
  151. * Normal termination
  152. *
  153. IRET=0
  154. RETURN
  155. *
  156. * Error handling
  157. *
  158. 9999 CONTINUE
  159. IRET=1
  160. WRITE(IOIMP,*) 'An error was detected in kres4.eso'
  161. RETURN
  162. *
  163. * End of KRES4
  164. *
  165. END
  166.  
  167.  
  168.  
  169.  
  170.  
  171.  
  172.  
  173.  
  174.  
  175.  
  176.  
  177.  
  178.  
  179.  
  180.  
  181.  
  182.  
  183.  
  184.  

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