Télécharger kres11.eso

Retour à la liste

Numérotation des lignes :

  1. C KRES11 SOURCE PV 16/11/17 22:00:21 9180
  2. SUBROUTINE KRES11(KMORS,KIZA,KTYP,ISMBR,
  3. $ KTYPI,ITER,RESID,ICALRS,IRSTRT,LBCG,BRTOL,IDDOT,IMVEC,
  4. $ KPREC,RXMILU,RXILUP,XLFIL,XDTOL,XSPIV,
  5. $ KTIME,LTIME,
  6. $ INCX,LRES,LNMV,ICVG,IMPR,INODET)
  7. IMPLICIT REAL*8 (A-H,O-Z)
  8. IMPLICIT INTEGER (I-N)
  9. C***********************************************************************
  10. C NOM : KRES11
  11. C DESCRIPTION : - Construction du préconditionneur
  12. C - Appel des solveurs itératifs
  13. C
  14. C
  15. C LANGAGE : ESOPE
  16. C AUTEUR : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  17. C mél : gounand@semt2.smts.cea.fr
  18. C***********************************************************************
  19. C VERSION : v1, 29/08/2011, version initiale
  20. C HISTORIQUE : v1, 29/08/2011, création
  21. C HISTORIQUE :
  22. C HISTORIQUE :
  23. C***********************************************************************
  24. -INC CCOPTIO
  25. -INC SMLENTI
  26. POINTEUR LNMV.MLENTI
  27. POINTEUR ATYP.MLENTI
  28. -INC SMLREEL
  29. POINTEUR LRES.MLREEL
  30. POINTEUR MAPREC.MATRIK
  31. POINTEUR AMORS.PMORS
  32. POINTEUR AISA.IZA
  33. -INC SMVECTD
  34. POINTEUR ISMBR.MVECTD
  35. POINTEUR INCX.MVECTD
  36. C Tableau de correspondance : entier <-> mot
  37. C pour le type d'inversion
  38. INTEGER NBTYPI,LNTYPI
  39. PARAMETER (NBTYPI=16)
  40. PARAMETER (LNTYPI=18)
  41. CHARACTER*(LNTYPI) LTYPI(NBTYPI)
  42. C Tableau de correspondance : entier <-> mot
  43. C pour le type de préconditionnement (cas d'une méthode itérative)
  44. INTEGER NBPREC,LNPREC
  45. PARAMETER (NBPREC=11)
  46. PARAMETER (LNPREC=8)
  47. CHARACTER*(LNPREC) LPREC(NBPREC)
  48. -INC SMTABLE
  49. POINTEUR KTIME.MTABLE
  50. DIMENSION ITTIME(4)
  51. CHARACTER*8 CHARI
  52. CHARACTER*1 CCOMP
  53. LOGICAL LTIME,LOGII
  54. C
  55. C Initialisation des tableaux
  56. C
  57. DATA LTYPI/'Choleski ',
  58. $ 'Conjugate Gradient',
  59. $ 'BiCGSTAB G ',
  60. $ 'BiCGSTAB(l) G ',
  61. $ 'GMRES(m) ',
  62. $ 'CGS-Neumaier ',
  63. $ 'Multigrid FCG ',
  64. $ 'Multigrid GCR(m) ',
  65. $ 'Bi-CG ',
  66. $ 'CG old ',
  67. $ 'BiCGSTAB old ',
  68. $ 'GMRES(m) old ',
  69. $ 'CGS ',
  70. $ 'BiCGSTAB ',
  71. $ 'BiCGSTAB(l) ',
  72. $ 'BiCGSTAB(l)F '/
  73. DATA LPREC/'None ',
  74. $ 'Jacobi ',
  75. $ 'D-ILU ',
  76. $ 'ILU(0) ',
  77. $ 'MILU(0) ',
  78. $ 'ILUT ',
  79. $ 'ILUT2 ',
  80. $ 'ILUTP ',
  81. $ 'ILUTP+0 ',
  82. $ 'ILU0-PV ',
  83. $ 'ILU0-PVf'/
  84.  
  85. IVALI=0
  86. XVALI=0.D0
  87. LOGII=.FALSE.
  88. IRETI=0
  89. XVALR=0.D0
  90. IOBRE=0
  91. IRETR=0
  92. *
  93. * Executable statements
  94. *
  95. * WRITE(IOIMP,*) 'Entrée dans kres10.eso'
  96. ICVG=0
  97. *Debug IF (KTYPI.EQ.1) KTYPI=3
  98. IF (KTYPI.LT.2.OR.KTYPI.GT.NBTYPI) THEN
  99. WRITE(IOIMP,*) 'KTYPI incorrect (2..',NBTYPI,') :',KTYPI
  100. GOTO 9999
  101. ENDIF
  102. IF (KPREC.LT.0.OR.KPREC.GT.NBPREC-1) THEN
  103. WRITE(IOIMP,*) 'PRECOND ',
  104. $ 'incorrect (0..',NBPREC-1,') :',KPREC
  105. GOTO 9999
  106. ENDIF
  107. C
  108. C Impressions
  109. C
  110. IF (IMPR.GT.1) THEN
  111. IF (IDDOT.EQ.0) CCOMP=' '
  112. IF (IDDOT.EQ.1) CCOMP='c'
  113. IF (KTYPI.EQ.5.OR.KTYPI.EQ.8) THEN
  114. WRITE(IOIMP,*) CCOMP,LTYPI(KTYPI),' (m=',IRSTRT,')'
  115. ELSEIF (KTYPI.EQ.11.OR.KTYPI.EQ.12) THEN
  116. WRITE(IOIMP,*) CCOMP,LTYPI(KTYPI),' (l=',LBCG,')'
  117. ELSE
  118. WRITE(IOIMP,*) CCOMP,LTYPI(KTYPI)
  119. ENDIF
  120. IF (KPREC.EQ.4) THEN
  121. WRITE(IOIMP,110) ' ',LPREC(KPREC+1),' (r=',RXMILU,')'
  122. 110 FORMAT (3A,D9.2,A)
  123. ELSEIF (KPREC.GE.5.AND.KPREC.LE.6) THEN
  124. WRITE(IOIMP,111) ' ',LPREC(KPREC+1),' (lfil=',XLFIL,
  125. $ ' droptol=',XDTOL,')'
  126. 111 FORMAT (3A,D9.2,A,D9.2,A)
  127. ELSEIF (KPREC.GE.7.AND.KPREC.LE.8) THEN
  128. WRITE(IOIMP,112) ' ',LPREC(KPREC+1),' (lfil=',XLFIL,
  129. $ ' droptol=',XDTOL,' pivtol=',XSPIV,
  130. $ ')'
  131. 112 FORMAT (3A,D9.2,A,D9.2,A,D9.2,A,I4,A)
  132. ELSEIF (KPREC.EQ.10) THEN
  133. WRITE(IOIMP,110) ' ',LPREC(KPREC+1),' (r=',RXILUP,')'
  134. ELSE
  135. WRITE(IOIMP,*) LPREC(KPREC+1)
  136. ENDIF
  137. ENDIF
  138. IF (LTIME) THEN
  139. CALL TIMESPV(ITTIME)
  140. ITI1=(ITTIME(1)+ITTIME(2))/10
  141. ENDIF
  142. C
  143. C - Construction du préconditionneur (repris sur kres5)
  144. C
  145. SEGACT ISMBR
  146. INC=ISMBR.VECTBB(/1)
  147. NRIGE=0
  148. NMATRI=0
  149. NKID=9
  150. NKMT=7
  151. SEGINI MAPREC
  152. MAPREC.KNTTT=INC
  153. MAPREC.KSYM=2
  154. NTT=0
  155. NPT=0
  156. NBLK=0
  157. SEGINI IDMAT
  158. MAPREC.KIDMAT(1)=IDMAT
  159. *
  160. IF (KPREC.EQ.1) THEN
  161. CALL MEIDIA(KMORS,KIZA,MAPREC,IMPR,IRET)
  162. IF (IRET.NE.0) GOTO 9999
  163. ELSEIF (KPREC.EQ.2) THEN
  164. CALL PRDILU(KMORS,KIZA,MAPREC,IMPR,IRET)
  165. IF (IRET.NE.0) GOTO 9999
  166. ELSEIF (KPREC.EQ.3) THEN
  167. CALL PRILU0(KMORS,KIZA,MAPREC,IMPR,IRET)
  168. IF (IRET.NE.0) GOTO 9999
  169. ELSEIF (KPREC.EQ.4) THEN
  170. CALL PRMILU(KMORS,KIZA,MAPREC,RXMILU,IMPR,IRET)
  171. IF (IRET.NE.0) GOTO 9999
  172. ELSEIF (KPREC.GE.5.AND.KPREC.LE.6) THEN
  173. IVARI=KPREC-5
  174. CALL PRILUT(KMORS,KIZA,MAPREC,XLFIL,XDTOL,IVARI,
  175. $ IMPR,IRET)
  176. IF (IRET.NE.0) GOTO 9999
  177. * WRITE(IOIMP,*) 'PRILUT !'
  178. * ELSEIF (KPREC.GE.7.AND.KPREC.LE.9) THEN
  179. ELSEIF (KPREC.GE.7.AND.KPREC.LE.8) THEN
  180. WRITE(IOIMP,*) 'KPREC=',KPREC, 'non dispo'
  181. GOTO 9999
  182. ELSEIF (KPREC.EQ.9.OR.KPREC.EQ.10) THEN
  183. CALL PRILUP(KMORS,KIZA,MAPREC,RXILUP,KPREC-9,IMPR,IRET)
  184. IF (IRET.NE.0) GOTO 9999
  185. ELSEIF (KPREC.NE.0) THEN
  186. WRITE(IOIMP,*) 'Erreur de programmation'
  187. GOTO 9999
  188. ENDIF
  189. IF (LTIME) THEN
  190. CALL TIMESPV(ITTIME)
  191. ITI2=(ITTIME(1)+ITTIME(2))/10
  192. ENDIF
  193. C MATRIK=MAPREC
  194. C SEGACT MATRIK
  195. C nkid=kidmat(/1)
  196. C WRITE(IOIMP,*) 'Apres constr precon'
  197. C WRITE(IOIMP,*) ' MAPREC KIDMAT NKID=',nkid
  198. C WRITE(IOIMP,2020) (KIDMAT(II),II=1,KIDMAT(/1))
  199. C 2020 FORMAT (20(2X,I12) )
  200. C
  201. C - Appel des solveurs itératifs
  202. C Apparemment, le segment MATRIK ne sert qu'à vérifier la symétrie
  203. C et à construire INCTYP pour le multigrille
  204. C
  205. NRIGE=0
  206. NMATRI=0
  207. NKID=0
  208. NKMT=0
  209. SEGINI MATRIK
  210. IF (KTYPI.NE.2) KSYM=2
  211. SEGACT ISMBR
  212. INC=ISMBR.VECTBB(/1)
  213. SEGINI INCX
  214. KS2B=ISMBR
  215. AMORS=KMORS
  216. AISA=KIZA
  217. ATYP=KTYP
  218. *
  219. IF (KTYPI.EQ.2) THEN
  220. CALL PRCG2(AMORS,AISA,KS2B,MATRIK,MAPREC,LRES,LNMV,INCX,
  221. $ ITER,INMV,RESID,KPREC,BRTOL,ICALRS,IDDOT,IMVEC,
  222. $ IMPR,IRET)
  223. ELSEIF (KTYPI.EQ.3) THEN
  224. LBCG=1
  225. CALL PRBCGG(AMORS,AISA,KS2B,MATRIK,MAPREC,LRES,LNMV,INCX,
  226. $ ITER,INMV,RESID,KPREC,BRTOL,LBCG,ICALRS,IDDOT,IMVEC,
  227. $ IMPR,IRET)
  228. ELSEIF (KTYPI.EQ.4) THEN
  229. CALL PRBCGG(AMORS,AISA,KS2B,MATRIK,MAPREC,LRES,LNMV,INCX,
  230. $ ITER,INMV,RESID,KPREC,BRTOL,LBCG,ICALRS,IDDOT,IMVEC,
  231. $ IMPR,IRET)
  232. ELSEIF (KTYPI.EQ.5) THEN
  233. C WRITE(IOIMP,*) 'ITER=',ITER
  234. C WRITE(IOIMP,*) 'INMV=',INMV
  235. C WRITE(IOIMP,*) 'RESID=',RESID
  236. C WRITE(IOIMP,*) 'KPREC=',KPREC
  237. C WRITE(IOIMP,*) 'IRSTRT=',IRSTRT
  238. C WRITE(IOIMP,*) 'ICALRS=',ICALRS
  239. C WRITE(IOIMP,*) 'IDDOT=',IDDOT
  240. CALL PRGMR2(AMORS,AISA,KS2B,MATRIK,MAPREC,LRES,LNMV,INCX,
  241. $ ITER,INMV,RESID,KPREC,IRSTRT,ICALRS,IDDOT,IMVEC,
  242. $ IMPR,IRET)
  243. ELSEIF (KTYPI.EQ.6) THEN
  244. CALL PRCGSN(AMORS,AISA,KS2B,MATRIK,MAPREC,LRES,LNMV,INCX,
  245. $ ITER,INMV,RESID,KPREC,BRTOL,ICALRS,IDDOT,IMVEC,
  246. $ IMPR,IRET)
  247. ELSEIF (KTYPI.EQ.7) THEN
  248. * WRITE(IOIMP,*) 'Appel de pragmg'
  249. CALL PRAGMG(AMORS,ATYP,AISA,KS2B,MATRIK,MAPREC,LRES,LNMV,
  250. $ INCX,
  251. $ ITER,INMV,RESID,KPREC,1,ICALRS,IDDOT,IMVEC,KTIME,LTIME,
  252. $ IMPR,IRET)
  253. SEGSUP ATYP
  254. ELSEIF (KTYPI.EQ.8) THEN
  255. CALL PRAGMG(AMORS,ATYP,AISA,KS2B,MATRIK,MAPREC,LRES,LNMV,
  256. $ INCX,
  257. $ ITER,INMV,RESID,KPREC,IRSTRT,ICALRS,IDDOT,IMVEC,KTIME,
  258. $ LTIME,
  259. $ IMPR,IRET)
  260. SEGSUP ATYP
  261. ELSEIF (KTYPI.EQ.9) THEN
  262. CALL PRBCG(AMORS,AISA,KS2B,MATRIK,MAPREC,LRES,LNMV,INCX,
  263. $ ITER,INMV,RESID,KPREC,BRTOL,ICALRS,IDDOT,IMVEC,
  264. $ IMPR,IRET)
  265. ELSE
  266. WRITE(IOIMP,*) 'KTYPI=',KTYPI,' invalide (1..',NBTYPI,')'
  267. GOTO 9999
  268. ENDIF
  269. C IRET<0 : 'vrai erreur' ou breakdown (dans le cas de BiCGSTAB)
  270. C IRET>0 : l'itération n'a pas convergé mais on veut
  271. C quand meme la solution calculée
  272. ICVG=IRET
  273. IF (IRET.GT.0) THEN
  274. WRITE(IOIMP,*) 'Convergence to tolerance not achieved : ',
  275. $ 'ITER=',ITER,' INMV=',INMV,' RESID=',RESID
  276. ELSEIF (IRET.EQ.0) THEN
  277. IF (IMPR.GT.0) THEN
  278. WRITE(IOIMP,*) 'ITER=',ITER,' INMV=',INMV,' RESID=',RESID
  279. ENDIF
  280. ELSEIF (IRET.LT.0) THEN
  281. WRITE(IOIMP,*) 'Error or breakdown in iterative method'
  282. GOTO 9999
  283. ELSE
  284. WRITE(IOIMP,*) 'KTYPI=',KTYPI,' invalide (1..',NBTYPI,')'
  285. GOTO 9999
  286. ENDIF
  287. SEGSUP MATRIK
  288. IF (LTIME) THEN
  289. CALL TIMESPV(ITTIME)
  290. ITI3=(ITTIME(1)+ITTIME(2))/10
  291. ITP=ITI2-ITI1
  292. ITR=ITI3-ITI2
  293. CHARI='PRECONDI'
  294. CALL ECCTAB(KTIME,'MOT ',IVALI,XVALI,CHARI,LOGII,IRETI,
  295. $ 'ENTIER ',ITP,XVALR,CHARR,LOGIR,IRETR)
  296. CHARI='RESOLUTI'
  297. CALL ECCTAB(KTIME,'MOT ',IVALI,XVALI,CHARI,LOGII,IRETI,
  298. $ 'ENTIER ',ITR,XVALR,CHARR,LOGIR,IRETR)
  299. ENDIF
  300. C
  301. C Destruction du préconditionneur
  302. C
  303. SEGACT MAPREC
  304. IDMAT=MAPREC.KIDMAT(1)
  305. IF (IDMAT.NE.0) THEN
  306. SEGACT IDMAT
  307. IZA=IDIAG
  308. SEGSUP IZA
  309. SEGSUP IDMAT
  310. ENDIF
  311. PMORS=MAPREC.KIDMAT(6)
  312. IF (PMORS.NE.0) SEGSUP PMORS
  313. IZA=MAPREC.KIDMAT(7)
  314. IF (IZA.NE.0) SEGSUP IZA
  315. SEGSUP MAPREC
  316. C
  317. C Destruction de la matrice Morse
  318. C
  319. IF (INODET.EQ.0) THEN
  320. SEGSUP AMORS
  321. SEGSUP AISA
  322. ENDIF
  323. *
  324. * Normal termination
  325. *
  326. RETURN
  327. *
  328. * Format handling
  329. *
  330. *
  331. * Error handling
  332. *
  333. 9999 CONTINUE
  334. MOTERR(1:8)='KRES11 '
  335. CALL ERREUR(349)
  336. RETURN
  337. *
  338. * End of subroutine KRES11
  339. *
  340. END
  341.  
  342.  
  343.  
  344.  
  345.  
  346.  
  347.  

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