Télécharger kres11.eso

Retour à la liste

Numérotation des lignes :

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

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