Télécharger kres5.eso

Retour à la liste

Numérotation des lignes :

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

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