Télécharger cgresi.eso

Retour à la liste

Numérotation des lignes :

  1. C CGRESI SOURCE CHAT 05/01/12 21:53:39 5004
  2. SUBROUTINE CGRESI(KSTO,
  3. 1 NL,ILG,
  4. 2 IMAT,IPRC,IA,JA,KA,
  5. 3 IZB,IZP,ICOLD,
  6. 4 NPT,NPITE,NEFF,ICONV,EPI,IPOU,VPOU,
  7. 5 NIMPR,IPAT)
  8. C
  9. CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
  10. C C
  11. C RESOLUTION ITERATIVE D'UN SYSTEME C
  12. C SYMETRIQUE CREUX PAR LA METHODE DU C
  13. C GRADIENT CONJUGUE PRECONDITIONNE PAR C
  14. C UN CHOLESKI INCOMPLET. C
  15. C C
  16. C***** INFOS SUR LA METHODE (LIGNE 0) C
  17. C C
  18. C DEUX MODES DE STOCKAGE SONT PROPOSES C
  19. C MORSE OU COMPRESSE. C
  20. C C
  21. C MORSE (KSTO=0) C
  22. C COMPRESSE (KSTO=1) C
  23. C C
  24. C POUR LE CAS ICG, ON RESOUD L'EQUATION SUIVANTE C
  25. C C
  26. C (K-1)*A*x = (K-1)*b C
  27. C C
  28. C AVEC K = LE CHOLESKI INCOMPLET DE LA MATRICE A C
  29. C C
  30. C***** INFOS SUR LES DIMENSIONS (LIGNE 1) C
  31. C C
  32. C NL (L'ORDRE DU SYSTEME) C
  33. C C
  34. C ILG (LONGUEUR LA POUR LE MORSE) C
  35. C (NOMBRE MAXI DE VALEURS PAR LIGNE DE LA C
  36. C MATRICE POUR LE COMPRESSE) C
  37. C C
  38. C***** INFOS SUR LA MATRICE (LIGNE 2) C
  39. C C
  40. C IMAT (LA MATRICE DANS UN LISTREEL) C
  41. C C
  42. C IPRC (LE PRECONDITIONNEMENT DANS UN LISTREEL) C
  43. C C
  44. C DCG : IPRC EST LE POINTEUR SUR UN LISTREEL C
  45. C CONTENANT D**(-1/2) CE LISTREEL A POUR C
  46. C DIMENSION NL AU MOINS C
  47. C C
  48. C ICG : IPRC EST LE POINTEUR SUR UN LISTREEL C
  49. C CONTENANT LA MATRICE DE C
  50. C PRECONDITIONNEMENT FACTORISEE C
  51. C LE TABLEAU KA D'ADRESSAGE EST LE MEME C
  52. C QUE CELUI DE LA MATRICE DE BASE. C
  53. C C
  54. C IA (TABLEAU DE POINTEURS DU STOCKAGE MORSE) C
  55. C C
  56. C IL DOIT ETRE DIMENSIONNE A NL+1 AU MINIMUM C
  57. C C
  58. C JA (TABLEAU DE CONNECTIVITES DE LA MATRICE) C
  59. C C
  60. C IL EST STOCKE EN MORSE ET DE LONGUEUR LA C
  61. C C
  62. C KA (TABLEAU DE CONNECTIVITES EN MODE COMPRESSE) C
  63. C C
  64. C IL EST DE LA FORME KA(NL,NNZ) AVEC NNZ C
  65. C LE NOMBRE MAXI DE TERMES NON NULS POUR UNE C
  66. C LIGNE DE LA MATRICE C
  67. C C
  68. C NB : EN MORSE, KA N'EST PAS UTILISE C
  69. C EN COMPRESSE, IA ET JA NE SONT PAS UTILISES C
  70. C C
  71. C IA,JA ET KA SONT DES LISTENTI. C
  72. C C
  73. C***** INFOS SUR LE SECOND MEMBRE, LA SOLUTION INITIALE C
  74. C***** ET LA SOLUTION FINALE. (LIGNE 3) C
  75. C C
  76. C IZB SECOND MEMBRE EN ENTREE C
  77. C INCHANGE EN SORTIE. C
  78. C C
  79. C IZP ESTIMATION DE LA SOLUTION EN ENTREE C
  80. C LA SOLUTION EN SORTIE, IL FAUT LA C
  81. C CONSERVER POUR LA RESOLUTION SUIVANTE C
  82. C C
  83. C IZB ET IZP SONT DES LISTREEL DIMENSIONNES A NL C
  84. C C
  85. C ICOLD =1 (DEMARRAGE FROID) C
  86. C ICOLD!=1 (DEMARRAGE CHAUD A PARTIR DE LA C
  87. C SOLUTION INITIALE FOURNIE DANS IZP) C
  88. C C
  89. C***** INFOS SUR LA CONVERGENCE (LIGNE 4) C
  90. C C
  91. C NPT NOMBRE MAXI D'ITERATIONS EN ENTREE C
  92. C SI IL EST < 10, IL EST MIS A 10 C
  93. C NPITE FREQUENCE DES TESTS DE CONVERGENCE C
  94. C NEFF NOMBRE REEL EN SORTIE C
  95. C ICONV 0 NON CONVERGE, 1 CONVERGE C
  96. C EPI PRECISION DEMANDEE C
  97. C IPOU=1 LA PRECISION VA ETRE POUSSE D'UN FACTEUR C
  98. C VPOU. CECI PEUT S'AVERER UTILE LORS D'UN C
  99. C CALCUL ITERATIF OU ON VEUT DES LE DEBUT C
  100. C UNE SOLUTION PRECISE POUR EVITER DE PARTIR C
  101. C SUR UNE MAUVAISE PISTE. VPOU A UNE VALEUR C
  102. C PAR DEFAUT DE 100. C
  103. C C
  104. C***** INFOS SUR LE CONTROLE DES IMPRESSIONS (LIGNE 5) C
  105. C C
  106. C NIMPR =1 (INFOS SUR LA CONVERGENCE) C
  107. C NIMPR!=1 (SILENCIEUX SAUF EN CAS DE DIVERGENCE) C
  108. C C
  109. C IPAT = (FACULTATIF, PAR EXEMPLE UN PAS DE TEMPS)C
  110. C C
  111. C RSETD(B,A) B <- A C
  112. C ADIVEC(A,B,C,alpha) C <- A + alpha*B C
  113. C PIMDV(A,B,C) C <- A*B C
  114. C C
  115. CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
  116. C
  117. IMPLICIT INTEGER(I-N)
  118. IMPLICIT REAL*8 (A-H,O-Z)
  119. -INC SMLENTI
  120. POINTEUR IA.MLENTI
  121. POINTEUR JA.MLENTI
  122. POINTEUR KA.MLENTI
  123. -INC SMLREEL
  124. POINTEUR IMAT.MLREEL
  125. POINTEUR IPRC.MLREEL
  126. POINTEUR IZB.MLREEL
  127. POINTEUR IZBB.MLREEL
  128. POINTEUR IZP.MLREEL
  129. POINTEUR IZPS.MLREEL
  130. POINTEUR IZZPR.MLREEL
  131. POINTEUR IZPMV.MLREEL
  132. POINTEUR IZTMP.MLREEL
  133. POINTEUR IZDIR.MLREEL
  134.  
  135. DATA VPOU0/100./
  136.  
  137. IF(KSTO.EQ.0) LA=ILG
  138. IF(KSTO.EQ.1) NNZ=ILG
  139. IF(IPOU.EQ.1) THEN
  140. VPOUE=VPOU0
  141. IF(VPOU.NE.VPOU0) VPOUE=VPOU
  142. ENDIF
  143.  
  144. CALL OOOMRU(1)
  145.  
  146. C On sauve le second membre en le dupliquant dans IZBB
  147. C qui lui sera ecrase.
  148. SEGACT IZB
  149. JG=IZB.PROG(/1)
  150. SEGINI IZBB
  151. CALL RSETD(IZBB.PROG,IZB.PROG,JG)
  152. SEGDES IZB
  153.  
  154. SEGACT IZP
  155.  
  156. JG=NL
  157. SEGINI IZPS
  158. SEGINI IZZPR
  159. SEGINI IZPMV
  160. SEGINI IZTMP
  161. SEGINI IZDIR
  162.  
  163. SEGACT IMAT
  164. SEGACT IPRC
  165. IF(KSTO.EQ.0) SEGACT IA
  166. IF(KSTO.EQ.0) SEGACT JA
  167. IF(KSTO.EQ.1) SEGACT KA
  168.  
  169. IF(NPT.LT.10) NPT=10
  170. IF(NPT.GT.NL) NPT=NL
  171.  
  172. IF(ICOLD.EQ.1) CALL INITD (IZP.PROG,NL,0.D0)
  173. CALL RSETD(IZPS.PROG,IZP.PROG,NL)
  174. IF(KSTO.EQ.0) THEN
  175. C CALL PMVM(A,IA,JA,X,Y,NL,LA)
  176. CALL PMVM(IMAT.PROG,IA.LECT,JA.LECT,
  177. 1 IZP.PROG,IZPMV.PROG,NL,LA)
  178. ELSE
  179. C CALL PMVC(A,KA,X,Y,NL,NNZ)
  180. CALL PMVC(IMAT.PROG,KA.LECT,
  181. 1 IZP.PROG,IZPMV.PROG,NL,NNZ)
  182. ENDIF
  183. CALL ADIVEC (IZBB.PROG,IZPMV.PROG,IZBB.PROG,-1.D0,NL)
  184. ANORB=SQRT(PISCAL(IZBB.PROG,IZBB.PROG,NL))
  185. IF(KSTO.EQ.0) THEN
  186. CALL ICGPRM(IZZPR.PROG,IZBB.PROG,IPRC.PROG,IA.LECT,JA.LECT,NL)
  187. ELSE
  188. CALL ICGPRC(IZZPR.PROG,IZBB.PROG,IPRC.PROG,KA.LECT,NL,NNZ)
  189. ENDIF
  190. CALL RSETD(IZDIR.PROG,IZZPR.PROG,NL)
  191. SCAR1=PISCAL(IZBB.PROG,IZDIR.PROG,NL)
  192. C
  193. C===========DEBUT DE LA BOUCLE ITERATIVE=====600========
  194. C
  195. NEFF=0
  196. NCOMP=0
  197. ICONV=0
  198. DO 600 I=1,NPT
  199. NEFF=NEFF+1
  200. NCOMP=NCOMP+1
  201. IF(KSTO.EQ.0) THEN
  202. C CALL PMVM(A,IA,JA,X,Y,NL,LA)
  203. CALL PMVM(IMAT.PROG,IA.LECT,JA.LECT,
  204. 1 IZDIR.PROG,IZPMV.PROG,NL,LA)
  205. ELSE
  206. C CALL PMVC(A,KA,X,Y,NL,NNZ)
  207. CALL PMVC(IMAT.PROG,KA.LECT,
  208. 1 IZDIR.PROG,IZPMV.PROG,NL,NNZ)
  209. ENDIF
  210. SCADIR=PISCAL(IZDIR.PROG,IZPMV.PROG,NL)
  211. ALPHA=SCAR1/SCADIR
  212. CALL ADIVEC(IZP.PROG,IZDIR.PROG,IZP.PROG,ALPHA,NL)
  213. ALP=-1.D0*ALPHA
  214. CALL ADIVEC (IZBB.PROG,IZPMV.PROG,IZBB.PROG,ALP,NL)
  215. IF(KSTO.EQ.0) THEN
  216. CALL ICGPRM(IZZPR.PROG,IZBB.PROG,IPRC.PROG,IA.LECT,JA.LECT,NL)
  217. ELSE
  218. CALL ICGPRC(IZZPR.PROG,IZBB.PROG,IPRC.PROG,KA.LECT,NL,NNZ)
  219. ENDIF
  220. SCAR2=PISCAL(IZBB.PROG,IZZPR.PROG,NL)
  221. BETA=SCAR2/SCAR1
  222. CALL ADIVEC (IZZPR.PROG,IZDIR.PROG,IZDIR.PROG,BETA,NL)
  223. IF (NCOMP.EQ.NPITE.OR.I.EQ.NPT) THEN
  224. IF (I.EQ.NPT) INDIC=1
  225. C On pousse la precision si necessaire
  226. ZEPS=EPI
  227. IF(IPOU.EQ.1) ZEPS=EPI/VPOUE
  228. CALL COTEST (IZP.PROG,IZPS.PROG,
  229. 1 ZEPS,DELTA,DELTAP,
  230. 2 INDIC,ICONV,NL)
  231. SEPA=DELTA
  232. RSEPA=DELTA/DELTAP
  233. CALL RSETD(IZPS.PROG,IZP.PROG,NL)
  234. NCOMP=0
  235. ENDIF
  236. IF (ICONV.EQ.1) GOTO 601
  237. SCAR1=SCAR2
  238. 600 CONTINUE
  239. 601 CONTINUE
  240. CALL MTABL(IZBB.PROG,RES,IMA,NL)
  241. RRES=RES/ANORB
  242. IF (ICONV.NE.1) THEN
  243. WRITE(6,1111)IPAT,NPT,RSEPA,DELTAP,RES,IZBB.PROG(IMA)
  244. ELSE
  245. IF(NIMPR.EQ.1) WRITE(6,2222)IPAT,EPI,
  246. 1 RES,SEPA,I,
  247. 2 RRES,RSEPA,IMA
  248. ENDIF
  249.  
  250. SEGDES IZP
  251.  
  252. SEGSUP IZBB
  253. SEGSUP IZPS
  254. SEGSUP IZZPR
  255. SEGSUP IZDIR
  256. SEGSUP IZTMP
  257. SEGSUP IZPMV
  258.  
  259. SEGDES IMAT*(NOMOD,MRU)
  260. SEGDES IPRC*(NOMOD,MRU)
  261. IF(KSTO.EQ.0) SEGDES IA*(NOMOD,MRU)
  262. IF(KSTO.EQ.0) SEGDES JA*(NOMOD,MRU)
  263. IF(KSTO.EQ.1) SEGDES KA*(NOMOD,MRU)
  264.  
  265. CALL OOOMRU(0)
  266.  
  267. C
  268. 1111 FORMAT(1X,'*CGRESI* ',I6,' NON CONV EN ',I5,
  269. 1' ITERA, RSEPA = ',D8.2,1X,'DELTAP = ',D8.2,1X,
  270. 2' RES = ',D8.2,' P = ',D8.2)
  271. 2222 FORMAT(1X,'*CGRESI* ',I6,' *PREC ',D8.2,
  272. 1 ' * RESIDU-SEPARA @ NITE ',D8.2,1X,D8.2,1X,I5,
  273. 2 ' * RRES-RSEP ',D8.2,1X,D8.2,
  274. 3 ' * ELEMENT ',I7)
  275. C
  276. C FIN RESOLUTION
  277. C
  278. RETURN
  279. END
  280.  
  281.  
  282.  

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