Télécharger ellp58.eso

Retour à la liste

Numérotation des lignes :

  1. C ELLP58 SOURCE CHAT 05/01/12 23:35:55 5004
  2. SUBROUTINE ELLP58(ZA,ZB,ZU,N,N2,NNPOI,ICHAR,NIT,
  3. * A,B,U,R,R1,P,P1,CH,EPS)
  4. C
  5. IMPLICIT INTEGER(I-N)
  6. IMPLICIT REAL*8(A-H,O-Y)
  7. IMPLICIT COMPLEX*16(Z)
  8. C
  9. C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
  10. C
  11. C OPERATEUR ...
  12. C
  13. C RESOLUTION D'UN SYSTEME LINEAIRE COMPLEXE ZA1 * ZX = ZV
  14. C PAR LA METHODE ITERATIVE DU GRADIANT CONJUGUE
  15. C
  16. C PARAMETRES :
  17. C
  18. C ZA(N*N) : MATRICE (N*N) REPRESENTANT LE SYSTEME LINEAIRE
  19. C ZB(N) : SECOND MEMBRE
  20. C ZU(N) : TABLEAU SOLUTION
  21. C
  22. C INTERMEDIAIRES :
  23. C
  24. C A(2N*2N) :
  25. C B(2N) :---> MATRICES REELLES DU SYSTEME REEL CORRESPONDANT
  26. C U(2N) :
  27. C
  28. C AUTRES MATRICES : MATRICES DE TRAVAIL
  29. C
  30. C SORTIE :
  31. C
  32. C ZX : SOLUTION
  33. C
  34. C AUTEUR : SAINT - DIZIER
  35. C DATE : 16 FEV 1990
  36. C
  37. C ====================================================================
  38. C
  39. COMPLEX*16 ZA(N,*),ZB(*),ZU(*)
  40. C COMPLEX*16 ZR(*),ZDELTA(*),ZP(*),ZH(*),ZH1(*)
  41. C
  42. REAL*8 A(N2,*)
  43. REAL*8 B(*)
  44. REAL*8 U(*)
  45. REAL*8 R(*)
  46. REAL*8 R1(*)
  47. REAL*8 P(*)
  48. REAL*8 P1(*)
  49. REAL*8 CH(*)
  50. C
  51. ZI = (0.D0 , -1.D0)
  52. N2 = N*2
  53. C
  54. DO 1 I = 1 , N
  55. B(I ) = ZB(I)
  56. B(I+N) = ZB(I)*ZI
  57. U(I ) = ZU(I)
  58. U(I+N) = ZU(I)*ZI
  59. DO 2 J = 1 , N
  60. A(I ,J ) = ZA(I,J)
  61. A(I+N,J+N) = A(I,J)
  62. A(I+N,J ) = ZA(I,J)*ZI
  63. A(I ,J+N) = -A(I+N,J)
  64. 2 CONTINUE
  65. 1 CONTINUE
  66. C
  67. C
  68. C EPS = 1.D-8
  69. C
  70. C
  71. C ............................... R = B - A*U
  72. C
  73. DO 10 I = 1 , N2
  74. R(I) = B(I)
  75. DO 11 J = 1 , N2
  76. R(I) = R(I) - A(I,J)*U(J)
  77. 11 CONTINUE
  78. 10 CONTINUE
  79. C
  80. C ............................... P1 = P = R1 = R
  81. C
  82. DO 12 I = 1 , N2
  83. R1(I) = R(I)
  84. P (I) = R(I)
  85. P1(I) = R(I)
  86. 12 CONTINUE
  87. C
  88. NIT = 0
  89. 100 NIT = NIT + 1
  90. C ............................... CH = (A) * P
  91. C
  92. DO 20 I = 1 , N2
  93. CH (I) = 0.D0
  94. DO 21 J = 1 , N2
  95. CH (I) = CH(I) + A(I,J)*P(J)
  96. 21 CONTINUE
  97. 20 CONTINUE
  98. C
  99. C ............................... DEL1 = (R)T * R1
  100. C
  101. DEL1 = 0.D0
  102. DO 30 I = 1 , N2
  103. DEL1 = DEL1 + R(I)*R1(I)
  104. 30 CONTINUE
  105. C
  106. C ............................... ALFA = DEL1/(CH)T*P1
  107. C
  108. D = 0.D0
  109. DO 40 I = 1 , N2
  110. D = D + CH(I)*P1(I)
  111. 40 CONTINUE
  112. C
  113. ALFA = DEL1 / D
  114. C
  115. C ............................... ZU = ZU + ZLANDA*ZP
  116. C ............................... ZDELTA = ZDELTA + ZLANDA*ZH1
  117. C ............................... DELTA1 = (ZDELTA)T * ZDELTA
  118. C
  119. DEL2 = 0.D0
  120. DO 80 I = 1 , N2
  121. U(I) = U(I) + ALFA * P(I)
  122. R(I) = R(I) - ALFA * CH(I)
  123. 80 CONTINUE
  124. C
  125. DO 81 I = 1 , N2
  126. CH (I) = 0.D0
  127. DO 82 J = 1 , N2
  128. CH (I) = CH(I) + A(J,I)*P1(J)
  129. 82 CONTINUE
  130. C
  131. R1(I) = R1(I) - ALFA*CH(I)
  132. DEL2 = DEL2 + R(I)*R1(I)
  133. C
  134. 81 CONTINUE
  135. C
  136. BETA = DEL2 / DEL1
  137. C
  138. XNU = 0.D0
  139. XDU = 0.D0
  140. XNR = 0.D0
  141. XDR = 0.D0
  142. XNF = 0.D0
  143. XDF = 0.D0
  144. XNM = 0.D0
  145. XDM = 0.D0
  146. ND12 = N / 12
  147. DO 85 I = 1 , ND12
  148. II = (I-1)*12
  149. DO 86 J = 1 , 3
  150. XNU = XNU + ALFA*ALFA*(P(II+J )**2 + P(II+J +N)**2)
  151. XNR = XNR + ALFA*ALFA*(P(II+J+3)**2 + P(II+J+3+N)**2)
  152. XNF = XNF + ALFA*ALFA*(P(II+J+6)**2 + P(II+J+6+N)**2)
  153. XNM = XNM + ALFA*ALFA*(P(II+J+9)**2 + P(II+J+9+N)**2)
  154. C
  155. XDU = XDU + U(II+J )**2 + U(II+J +N)**2
  156. XDR = XDR + U(II+J+3)**2 + U(II+J+3+N)**2
  157. XDF = XDF + U(II+J+6)**2 + U(II+J+6+N)**2
  158. XDM = XDM + U(II+J+9)**2 + U(II+J+9+N)**2
  159. 86 CONTINUE
  160. 85 CONTINUE
  161. C
  162. C NN = (NNPOI-1)*12 + ICHAR
  163. C XNUM = P(NN)**2 +P(NN+N)**2
  164. C XDEN = U(NN)**2 +U(NN+N)**2
  165. C
  166. XST1 = SQRT(XNU/XDU)
  167. XST2 = SQRT(XNR/XDR)
  168. XST3 = SQRT(XNF/XDF)
  169. XST4 = SQRT(XNM/XDM)
  170. XSTOP = SQRT(XST1*XST1+XST2*XST2+XST3*XST3+XST4*XST4)
  171. C WRITE(6,*)'XST1 :',XST1
  172. C WRITE(6,*)'XST2 :',XST2
  173. C WRITE(6,*)'XST3 :',XST3
  174. C WRITE(6,*)'XST4 :',XST4
  175. C
  176. IF (XSTOP.LE.EPS) THEN
  177. DO 88 I = 1 , N
  178. ZU(I) = U(I) - ZI*U(I+N)
  179. 88 CONTINUE
  180. WRITE(6,*) 'VECTEUR SOLUTION :',NNPOI,ICHAR
  181. C DO 89 I = 1 , N
  182. C WRITE(6,*) I , ZU(I)
  183. C9 CONTINUE
  184. C
  185. RETURN
  186. END IF
  187. C
  188. C
  189. C ************************************************************
  190. C
  191. DO 90 I = 1 , N2
  192. P (I) = BETA*P (I) + R (I)
  193. P1(I) = BETA*P1(I) + R1(I)
  194. 90 CONTINUE
  195. C
  196. GO TO 100
  197. C
  198. END
  199.  
  200.  
  201.  

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