Télécharger kres5.eso

Retour à la liste

Numérotation des lignes :

kres5
  1. C KRES5 SOURCE GOUNAND 22/08/25 21:15:07 11434
  2. SUBROUTINE KRES5(MATRIK,KCLIM,KSMBR,KTYPI,
  3. $ MCHINI,ITER,RESID,
  4. $ BRTOL,RESTRT,LBCG,ICALRS,
  5. $ MAPREC,KPREC,
  6. $ RXMILU,RXILUP,XLFIL,XDTOL,XSPIV,
  7. $ ISCAL,
  8. $ KTIME,LTIME,IDDOT,IMVEC,
  9. $ MCHSOL,LRES,LNMV,ICVG,
  10. $ IMPR,IRET)
  11. IMPLICIT REAL*8 (A-H,O-Z)
  12. IMPLICIT INTEGER (I-N)
  13. C***********************************************************************
  14. C NOM : KRES5
  15. C DESCRIPTION : Résolution d'un système par une méthode itérative
  16. C (type gradient conjugué).
  17. C
  18. C
  19. C LANGAGE : ESOPE
  20. C AUTEUR : Stéphane GOUNAND (CEA/DRN/DMT/SEMT/LTMF)
  21. C mél : gounand@semt2.smts.cea.fr
  22. C***********************************************************************
  23. C APPELES : VERMAT, MEXINI, MESMBR, MELIM, INFMAT,
  24. C MEIDIA, PRDILU, PRILU, PRMILU, PRILUT,
  25. C PRCG, PRBCGS, PRBCS2, PRGMR,
  26. C XDISP
  27. C APPELE PAR : KRES2
  28. C***********************************************************************
  29. C ENTREES : MATRIK, KCLIM, KSMBR, KTYPI, MCHINI, ITER, RESID,
  30. C BRTOL, RESTRT, MAPREC, KPREC, RXMILU, XLFIL,
  31. C XDTOL
  32. C ENTREES/SORTIES : -
  33. C SORTIES : MCHSOL, LRES
  34. C CODE RETOUR (IRET) : = 0 si tout s'est bien passé
  35. C***********************************************************************
  36. C VERSION : v1, 14/04/2000, version initiale
  37. C HISTORIQUE : v1, 14/04/2000, création
  38. C HISTORIQUE : 06/04/04 : Scaling
  39. C HISTORIQUE : 08/04/04 : ajout ILUTP
  40. C HISTORIQUE :
  41. C***********************************************************************
  42. C Prière de PRENDRE LE TEMPS de compléter les commentaires
  43. C en cas de modification de ce sous-programme afin de faciliter
  44. C la maintenance !
  45. C***********************************************************************
  46.  
  47. -INC PPARAM
  48. -INC CCOPTIO
  49. -INC SMLENTI
  50. POINTEUR LNMV.MLENTI
  51. POINTEUR ATYP.MLENTI
  52. -INC SMLREEL
  53. POINTEUR LRES.MLREEL
  54. -INC SMCHPOI
  55. POINTEUR KCLIM.MCHPOI
  56. POINTEUR KSMBR.MCHPOI
  57. POINTEUR MCHINI.MCHPOI
  58. POINTEUR MCHSOL.MCHPOI
  59. POINTEUR MAPREC.MATRIK
  60. POINTEUR INCX.IZA
  61. POINTEUR KS2B.IZA
  62. POINTEUR KMORS.PMORS
  63. POINTEUR KISA.IZA
  64. POINTEUR AMORS.PMORS
  65. POINTEUR AISA.IZA
  66. POINTEUR NORMP.IZA
  67. POINTEUR NORMD.IZA
  68. *
  69. * .. Parameters
  70. REAL*8 ZERO ,ONE
  71. PARAMETER (ZERO=0.0D0,ONE=1.0D0)
  72. * Paramètre m du GMRES(m)
  73. INTEGER RESTRT
  74. INTEGER ITER,IVARI
  75. REAL*8 BRTOL,RESID,RXMILU,RXILUP
  76. *
  77. REAL*8 XLFIL,XDTOL
  78. INTEGER KPREC
  79. INTEGER KTYPI
  80. C LOGICAL LCLZER
  81. LOGICAL LRACOU
  82. LOGICAL LFIRST
  83. C
  84. INTEGER IMPR,IRET
  85. C Tableau de correspondance : entier <-> mot
  86. C pour le type d'inversion
  87. INTEGER NBTYPI,LNTYPI
  88. PARAMETER (NBTYPI=9)
  89. PARAMETER (LNTYPI=18)
  90. CHARACTER*(LNTYPI) LTYPI(NBTYPI)
  91. C Tableau de correspondance : entier <-> mot
  92. C pour le type de préconditionnement (cas d'une méthode itérative)
  93. INTEGER NBPREC,LNPREC
  94. PARAMETER (NBPREC=11)
  95. PARAMETER (LNPREC=8)
  96. CHARACTER*(LNPREC) LPREC(NBPREC)
  97. -INC SMTABLE
  98. POINTEUR KTIME.MTABLE
  99. DIMENSION ITTIME(4)
  100. CHARACTER*8 CHARI
  101. CHARACTER*1 CCOMP
  102. LOGICAL LTIME,LOGII
  103. C
  104. C Initialisation des tableaux
  105. C
  106. DATA LTYPI/'Choleski ',
  107. $ 'Conjugate Gradient',
  108. $ 'BiCGSTAB G ',
  109. $ 'BiCGSTAB(l) G ',
  110. $ 'GMRES(m) ',
  111. $ 'CGS-Neumaier ',
  112. $ 'Multigrid FCG ',
  113. $ 'Multigrid GCR(m) ',
  114. $ 'Bi-CG '/
  115. DATA LPREC/'None ',
  116. $ 'Jacobi ',
  117. $ 'D-ILU ',
  118. $ 'ILU(0) ',
  119. $ 'MILU(0) ',
  120. $ 'ILUT ',
  121. $ 'ILUT2 ',
  122. $ 'ILUTP ',
  123. $ 'ILUTP+0 ',
  124. $ 'ILU0-PV ',
  125. $ 'ILU0-PVf'/
  126.  
  127. IVALI=0
  128. XVALI=0.D0
  129. LOGII=.FALSE.
  130. IRETI=0
  131. XVALR=0.D0
  132. IOBRE=0
  133. IRETR=0
  134. *
  135. * Executable statements
  136. *
  137. ICVG=0
  138. IF (IMPR.GT.5) WRITE(IOIMP,*) 'Entrée dans kres5.eso'
  139. C Batterie de tests
  140. C
  141. C KTYPI : Méthode d'inversion du système (cf. LTYPI)
  142. C 1 : résolution directe (Choleski)
  143. C 2 : Gradient Conjugué
  144. C 3 : Bi-Gradient Conjugué Stabilisé (BiCGSTAB)
  145. C 4 : BiCGSTAB(2)
  146. C 5 : GMRES(m)
  147. C 6 : CGS
  148. C 7 : Multigrille symétrique
  149. C 8 : Multigrille non symétrique
  150. C 9 : BiCG
  151. IF (KTYPI.LT.2.OR.KTYPI.GT.NBTYPI) THEN
  152. WRITE(IOIMP,*) 'KTYPI incorrect (2..',NBTYPI,') :',KTYPI
  153. GOTO 9999
  154. ENDIF
  155. * Pas de préconditionnement en cas de multigrille
  156. C - Type de préconditionnement : (cf. LPREC)
  157. C 0 : pas de préconditionnement
  158. C 1 : préconditionnement par la diagonale
  159. C 2 : préconditionnement D-ILU
  160. C 3 : préconditionnement ILU(0) (Choleski)
  161. C 4 : préconditionnement MILU(0) (Choleski modifié)
  162. C 5 : préconditionnement ILUT (dual truncation)
  163. C 6 : préconditionnement ILUT2 (une variante du
  164. C précédent qui remplit mieux la mémoire et
  165. C fonctionne mieux quelquefois)
  166. C 7 : préconditionnement ILUTP
  167. C 8 : préconditionnement ILUTP version Goon
  168. C 9 : préconditionnement ILU-PV
  169. C 10 : préconditionnement ILU-PVf
  170. IF (KPREC.LT.0.OR.KPREC.GT.NBPREC-1) THEN
  171. WRITE(IOIMP,*) 'PRECOND ',
  172. $ 'incorrect (0..',NBPREC-1,') :',KPREC
  173. GOTO 9999
  174. ENDIF
  175. C
  176. C On vérifie que la matrice est correctement assemblée
  177. C
  178. CALL VERMAT(MATRIK,IMPR,IRET)
  179. IF (IRET.NE.0) GOTO 9999
  180. C
  181. C On transforme le chpoint estimation initiale de l'inconnue en vecteur
  182. C estimation initiale de l'inconnue
  183. C In MEXINI : SEGINI INCX
  184. CALL MEXINI(MATRIK,MCHINI,
  185. $ INCX,
  186. $ IMPR,IRET)
  187. IF (IRET.NE.0) GOTO 9999
  188. C
  189. C On transforme le chpoint second membre en vecteur second membre
  190. C In MESMBR : SEGINI KS2B
  191. CALL MESMBR(MATRIK,KSMBR,
  192. $ KS2B,
  193. $ IMPR,IRET)
  194. IF (IRET.NE.0) GOTO 9999
  195. C
  196. C On applique les conditions aux limites
  197. C
  198. LRACOU=(KCLIM.EQ.0)
  199. * LRACOU=.FALSE.
  200. * WRITE(IOIMP,*) 'LRACOU=',LRACOU
  201. IF (LRACOU) THEN
  202. SEGACT MATRIK
  203. AMORS=MATRIK.KIDMAT(4)
  204. AISA=MATRIK.KIDMAT(5)
  205. SEGDES MATRIK
  206. ELSE
  207. C
  208. C In MELIM : SEGINI AMORS
  209. C SEGINI AISA
  210. CALL MELIM(MATRIK,KCLIM,
  211. $ INCX,KS2B,
  212. $ AMORS,AISA,
  213. $ IMPR,IRET)
  214. IF (IRET.NE.0) GOTO 9999
  215. ENDIF
  216. C
  217. LFIRST=.FALSE.
  218. IF (ISCAL.GT.0) THEN
  219. IF (LRACOU) THEN
  220. SEGACT MATRIK
  221. NORMP=MATRIK.KKMMT(4)
  222. NORMD=MATRIK.KKMMT(5)
  223. SEGDES MATRIK
  224. IF (NORMP.EQ.0.AND.NORMD.EQ.0) THEN
  225. LFIRST=.TRUE.
  226. ENDIF
  227. ELSE
  228. NORMP=0
  229. NORMD=0
  230. ENDIF
  231. IF (NORMP.EQ.0.OR.NORMD.EQ.0) THEN
  232. C Calcul des normes primales (colonnes) et duales (lignes)
  233. C de la matrice. Norme = norme L2, soit :
  234. C {\sum_{i ou j} a_{ij}^2}^{1/2}
  235. C
  236. CALL NORMAT(AMORS,AISA,ISCAL,NORMP,NORMD,IMPR,IRET)
  237. IF (IRET.NE.0) GOTO 9999
  238. C
  239. C On norme la matrice : attention modification...
  240. C
  241. CALL NORMAM(AMORS,AISA,NORMP,NORMD,IMPR,IRET)
  242. IF (IRET.NE.0) GOTO 9999
  243. ENDIF
  244. C
  245. C On norme le second membre : attention modification...
  246. C
  247. CALL NORMV2(KS2B,NORMD,IMPR,IRET)
  248. IF (IRET.NE.0) GOTO 9999
  249. C
  250. C On dénorme l'inconnue X : attention modification...
  251. C
  252. CALL NORMV1(INCX,NORMP,IMPR,IRET)
  253. IF (IRET.NE.0) GOTO 9999
  254. ELSE
  255. NORMP=0
  256. NORMD=0
  257. ENDIF
  258. C
  259. C On donne des infos sur la matrice
  260. C
  261. CALL INFMAT(AMORS,AISA,IMPR,IRET)
  262. IF (IRET.NE.0) GOTO 9999
  263. C
  264. C ON CHOISIT LE TYPE D'INVERSION
  265. C
  266. C Méthodes de gradient conjugué
  267. C
  268. C KTYPI = 2 : Gradient conjugué (CG)
  269. C KTYPI = 3 : Bi-Gradient conjugué stabilisé (BiCGSTAB)
  270. C KTYPI = 4 : BiCGSTAB(2)
  271. C KTYPI = 5 : GMRES(m)
  272. C KTYPI = 6 : CGS
  273. C KTYPI = 7 : Algebraic Multigrid FCG
  274. C KTYPI = 8 : Algebraic Multigrid GCR(m)
  275. C KTYPI = 9 : BiCG
  276. C
  277. C KPREC = 0 : pas de préconditionnement
  278. C KPREC = 1 : préconditionnement par la diagonale
  279. C KPREC = 2 : préconditionnement par factorisation incomplete
  280. C D-ILU
  281. C KPREC = 3 : préconditionnement par factorisation incomplete
  282. C ILU(0) (i.e. Crout)
  283. C KPREC = 4 : préconditionnement par factorisation incomplete
  284. C modifiée MILU(0) (i.e. Crout)
  285. C KPREC = 5 : préconditionnement par factorisation incomplete
  286. C ILUT (dual truncation strategy)
  287. C KPREC = 6 : préconditionnement ILUT2 (une variante du
  288. C précédent qui remplit mieux la mémoire et
  289. C fonctionne mieux quelquefois)
  290. C KPREC = 7 : préconditionnement ILUTP
  291. C KPREC = 8 : préconditionnement ILUTPGOO
  292. C KPREC = 9 : préconditionnement ILU-PV
  293. C KPREC = 10 : préconditionnement ILU-PV filtré
  294. IF(KTYPI.GE.2.AND.KTYPI.LE.NBTYPI)THEN
  295. IF (IMPR.GT.1) THEN
  296. IF (IDDOT.EQ.0) CCOMP=' '
  297. IF (IDDOT.EQ.1) CCOMP='c'
  298. IF (KTYPI.EQ.5.OR.KTYPI.EQ.8) THEN
  299. WRITE(IOIMP,*) CCOMP,LTYPI(KTYPI),' (m=',RESTRT,')'
  300. ELSEIF (KTYPI.EQ.11.OR.KTYPI.EQ.12) THEN
  301. WRITE(IOIMP,*) CCOMP,LTYPI(KTYPI),' (l=',LBCG,')'
  302. ELSE
  303. WRITE(IOIMP,*) CCOMP,LTYPI(KTYPI)
  304. ENDIF
  305. IF (KPREC.EQ.4) THEN
  306. WRITE(IOIMP,110) ' ',LPREC(KPREC+1),' (r=',RXMILU,')'
  307. 110 FORMAT (3A,D9.2,A)
  308. ELSEIF (KPREC.GE.5.AND.KPREC.LE.6) THEN
  309. WRITE(IOIMP,111) ' ',LPREC(KPREC+1),' (lfil=',XLFIL,
  310. $ ' droptol=',XDTOL,')'
  311. 111 FORMAT (3A,D9.2,A,D9.2,A)
  312. ELSEIF (KPREC.GE.7.AND.KPREC.LE.8) THEN
  313. WRITE(IOIMP,112) ' ',LPREC(KPREC+1),' (lfil=',XLFIL,
  314. $ ' droptol=',XDTOL,' pivtol=',XSPIV,
  315. $ ')'
  316. 112 FORMAT (3A,D9.2,A,D9.2,A,D9.2,A,I4,A)
  317. ELSEIF (KPREC.EQ.10) THEN
  318. WRITE(IOIMP,110) ' ',LPREC(KPREC+1),' (r=',RXILUP,')'
  319. ELSE
  320. WRITE(IOIMP,*) LPREC(KPREC+1)
  321. ENDIF
  322. ENDIF
  323. IF (LTIME) THEN
  324. call timespv(ittime,oothrd)
  325. ITI1=(ITTIME(1)+ITTIME(2))/10
  326. ENDIF
  327. C
  328. C Construction (éventuelle) du préconditionneur
  329. C
  330. * Gestion du CTRL-C
  331. if (ierr.NE.0) return
  332. IF (KPREC.EQ.1) THEN
  333. CALL MEIDIA(AMORS,AISA,MAPREC,IMPR,IRET)
  334. IF (IRET.NE.0) GOTO 9999
  335. ELSEIF (KPREC.EQ.2) THEN
  336. CALL PRDILU(AMORS,AISA,MAPREC,IMPR,IRET)
  337. IF (IRET.NE.0) GOTO 9999
  338. ELSEIF (KPREC.EQ.3) THEN
  339. CALL PRILU0(AMORS,AISA,MAPREC,IMPR,IRET)
  340. IF (IRET.NE.0) GOTO 9999
  341. ELSEIF (KPREC.EQ.4) THEN
  342. CALL PRMILU(AMORS,AISA,MAPREC,RXMILU,IMPR,IRET)
  343. IF (IRET.NE.0) GOTO 9999
  344. ELSEIF (KPREC.GE.5.AND.KPREC.LE.6) THEN
  345. IVARI=KPREC-5
  346. CALL PRILUT(AMORS,AISA,MAPREC,XLFIL,XDTOL,IVARI,
  347. $ IMPR,IRET)
  348. IF (IRET.NE.0) GOTO 9999
  349. * WRITE(IOIMP,*) 'PRILUT !'
  350. * ELSEIF (KPREC.GE.7.AND.KPREC.LE.9) THEN
  351. ELSEIF (KPREC.GE.7.AND.KPREC.LE.8) THEN
  352. * Ici, on pivote les colonnes lors de la construction du
  353. * préconditionneur...
  354. * On modifiera donc IDMATP et INCX et NORMP
  355. IF (KPREC.EQ.7) THEN
  356. IVARI=0
  357. ELSEIF (KPREC.EQ.8) THEN
  358. IVARI=2
  359. ENDIF
  360. CALL PRILTP(AMORS,AISA,MAPREC,XLFIL,XDTOL,XSPIV,
  361. $ IVARI,
  362. $ INCX,NORMP,LRACOU,
  363. $ IMPR,IRET)
  364. IF (IRET.NE.0) GOTO 9999
  365. ELSEIF (KPREC.EQ.9.OR.KPREC.EQ.10) THEN
  366. CALL PRILUP(AMORS,AISA,MAPREC,RXILUP,KPREC-9,IMPR,IRET)
  367. IF (IRET.NE.0) GOTO 9999
  368. ELSEIF (KPREC.NE.0) THEN
  369. WRITE(IOIMP,*) 'Erreur de programmation'
  370. GOTO 9999
  371. ENDIF
  372. IF (LTIME) THEN
  373. call timespv(ittime,oothrd)
  374. ITI2=(ITTIME(1)+ITTIME(2))/10
  375. ENDIF
  376. C
  377. C Nettoyage des zéros stockés (utilité ?)
  378. C
  379. IF (LFIRST) THEN
  380. * IIMPR=IMPR
  381. * IMPR=6
  382. CALL CLMORS(AMORS,AISA,IMPR,IRET)
  383. IF (IRET.NE.0) GOTO 9999
  384. * IMPR=IIMPR
  385. ENDIF
  386. * SEGPRT, AMORS
  387. * SEGPRT, AISA
  388. C
  389. C Résolution itérative proprement dite
  390. C
  391. * Gestion du CTRL-C
  392. if (ierr.NE.0) return
  393. IF (LTIME) THEN
  394. call timespv(ittime,oothrd)
  395. ITI3=(ITTIME(1)+ITTIME(2))/10
  396. ENDIF
  397. INMV=0
  398. IF (KTYPI.EQ.2) THEN
  399. CALL PRCG2(AMORS,AISA,KS2B,MATRIK,MAPREC,LRES,LNMV,INCX,
  400. $ ITER,INMV,RESID,KPREC,BRTOL,ICALRS,IDDOT,IMVEC,
  401. $ IMPR,IRET)
  402. ELSEIF (KTYPI.EQ.3) THEN
  403. LBCG=1
  404. CALL PRBCGG(AMORS,AISA,KS2B,MATRIK,MAPREC,LRES,LNMV,INCX,
  405. $ ITER,INMV,RESID,KPREC,BRTOL,LBCG,ICALRS,IDDOT,IMVEC,
  406. $ IMPR,IRET)
  407. ELSEIF (KTYPI.EQ.4) THEN
  408. CALL PRBCGG(AMORS,AISA,KS2B,MATRIK,MAPREC,LRES,LNMV,INCX,
  409. $ ITER,INMV,RESID,KPREC,BRTOL,LBCG,ICALRS,IDDOT,IMVEC,
  410. $ IMPR,IRET)
  411. ELSEIF (KTYPI.EQ.5) THEN
  412. CALL PRGMR2(AMORS,AISA,KS2B,MATRIK,MAPREC,LRES,LNMV,INCX,
  413. $ ITER,INMV,RESID,KPREC,RESTRT,ICALRS,IDDOT,IMVEC,
  414. $ IMPR,IRET)
  415. ELSEIF (KTYPI.EQ.6) THEN
  416. CALL PRCGSN(AMORS,AISA,KS2B,MATRIK,MAPREC,LRES,LNMV,INCX,
  417. $ ITER,INMV,RESID,KPREC,BRTOL,ICALRS,IDDOT,IMVEC,
  418. $ IMPR,IRET)
  419. ELSEIF (KTYPI.EQ.7) THEN
  420. CALL INCTYP(MATRIK,
  421. $ ATYP,
  422. $ IMPR,IRET)
  423. IF (IRET.NE.0) GOTO 9999
  424. * WRITE(IOIMP,*) 'Appel de pragmg'
  425. CALL PRAGMG(AMORS,ATYP,AISA,KS2B,MATRIK,MAPREC,LRES,LNMV,
  426. $ INCX,
  427. $ ITER,INMV,RESID,KPREC,1,ICALRS,IDDOT,IMVEC,KTIME,LTIME,
  428. $ IMPR,IRET)
  429. SEGSUP ATYP
  430. ELSEIF (KTYPI.EQ.8) THEN
  431. CALL INCTYP(MATRIK,
  432. $ ATYP,
  433. $ IMPR,IRET)
  434. IF (IRET.NE.0) GOTO 9999
  435. CALL PRAGMG(AMORS,ATYP,AISA,KS2B,MATRIK,MAPREC,LRES,LNMV,
  436. $ INCX,
  437. $ ITER,INMV,RESID,KPREC,RESTRT,ICALRS,IDDOT,IMVEC,
  438. $ KTIME,LTIME,
  439. $ IMPR,IRET)
  440. SEGSUP ATYP
  441. ELSEIF (KTYPI.EQ.9) THEN
  442. CALL PRBCG(AMORS,AISA,KS2B,MATRIK,MAPREC,LRES,LNMV,INCX,
  443. $ ITER,INMV,RESID,KPREC,BRTOL,ICALRS,IDDOT,IMVEC,
  444. $ IMPR,IRET)
  445. ELSE
  446. WRITE(IOIMP,*) 'KTYPI=',KTYPI,' incorrect'
  447. GOTO 9999
  448. ENDIF
  449. * Gestion du CTRL-C
  450. if (ierr.NE.0) return
  451. C IRET<0 : 'vrai erreur' ou breakdown (dans le cas de BiCGSTAB)
  452. C IRET>0 : l'itération n'a pas convergé mais on veut
  453. C quand meme la solution calculée
  454. ICVG=IRET
  455. IF (IRET.GT.0) THEN
  456. WRITE(IOIMP,*) 'Convergence to tolerance not achieved : ',
  457. $ 'ITER=',ITER,' INMV=',INMV,' RESID=',RESID
  458. ELSEIF (IRET.EQ.0) THEN
  459. IF (IMPR.GT.0) THEN
  460. WRITE(IOIMP,*) 'ITER=',ITER,' INMV=',INMV,' RESID=',RESID
  461. ENDIF
  462. ELSEIF (IRET.LT.0) THEN
  463. WRITE(IOIMP,*) 'Error or breakdown in iterative method'
  464. GOTO 9999
  465. ENDIF
  466. IF (LTIME) THEN
  467. call timespv(ittime,oothrd)
  468. ITI4=(ITTIME(1)+ITTIME(2))/10
  469. ITP=ITI2-ITI1
  470. ITR=ITI4-ITI3
  471. CHARI='PRECONDI'
  472. CALL ECCTAB(KTIME,'MOT ',IVALI,XVALI,CHARI,LOGII,IRETI,
  473. $ 'ENTIER ',ITP,XVALR,CHARR,LOGIR,IRETR)
  474. CHARI='RESOLUTI'
  475. CALL ECCTAB(KTIME,'MOT ',IVALI,XVALI,CHARI,LOGII,IRETI,
  476. $ 'ENTIER ',ITR,XVALR,CHARR,LOGIR,IRETR)
  477. ENDIF
  478. ELSE
  479. WRITE(IOIMP,*) 'KTYPI=',KTYPI,' invalide (1..',NBTYPI,')'
  480. GOTO 9999
  481. ENDIF
  482. IF (ISCAL.GT.0) THEN
  483. C
  484. C On dénorme le vecteur solution : attention modification...
  485. C
  486. CALL NORMV2(INCX,NORMP,IMPR,IRET)
  487. IF (IRET.NE.0) GOTO 9999
  488. IF (LRACOU) THEN
  489. SEGACT MATRIK*MOD
  490. MATRIK.KKMMT(4)=NORMP
  491. MATRIK.KKMMT(5)=NORMD
  492. SEGDES MATRIK
  493. ENDIF
  494. ENDIF
  495. C
  496. C Transformation du vecteur-solution en chpoint
  497. C
  498. CALL XDISP(MATRIK,INCX,MCHSOL,IMPR,IRET)
  499. IF (IRET.NE.0) GOTO 9999
  500. C
  501. C Suppression des objets temporaires
  502. C
  503. IF (LRACOU) THEN
  504. SEGACT MATRIK*MOD
  505. MATRIK.KIDMAT(4)=AMORS
  506. MATRIK.KIDMAT(5)=AISA
  507. SEGDES MATRIK
  508. ELSE
  509. SEGSUP,AMORS
  510. SEGSUP,AISA
  511. ENDIF
  512. SEGSUP INCX
  513. SEGSUP KS2B
  514. *
  515. * Normal termination
  516. *
  517. IRET=0
  518. RETURN
  519. *
  520. * Format handling
  521. *
  522. *
  523. * Error handling
  524. *
  525. 9999 CONTINUE
  526. IRET=1
  527. WRITE(IOIMP,*) 'An error was detected in kres5.eso'
  528. RETURN
  529. *
  530. * End of KRES5
  531. *
  532. END
  533.  
  534.  

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