Télécharger kres2.eso

Retour à la liste

Numérotation des lignes :

  1. C KRES2 SOURCE CB215821 19/07/30 21:17:08 10273
  2. SUBROUTINE KRES2()
  3. IMPLICIT REAL*8 (A-H,O-Z)
  4. IMPLICIT INTEGER (I-N)
  5. C*************************************************************************
  6. C Operateur KRES
  7. C
  8. C OBJET : Resoud l'equation AX=B par différentes méthodes
  9. C * directe (factorisation Choleski)
  10. C * itérative (Gradient conjugué, BiCGSTAB,
  11. C GMRES(m)) avec différents préconditionneurs
  12. C diagonal (Jacobi), D-ILU, ILU(0) (Choleski
  13. C incomplet), MILU(0)
  14. C SYNTAXE : CHPO3 = KRES MA1 'TYPI' TAB1
  15. C 'CLIM' CHPO1
  16. C 'SMBR' CHPO2
  17. C 'IMPR' VAL1 ;
  18. C Voir la notice pour plus de précisions.
  19. C
  20. C***********************************************************************
  21. C APPELES : KRES3, KRES4, KRES5
  22. C APPELES (E/S) : LIROBJ, ECROBJ, ERREUR, LIRMOT, LIRENT, LIRTAB,
  23. C ACME, ACMO, ACMM, ACMF, ECMO
  24. C APPELES (UTIL.) : QUETYP
  25. C APPELE PAR : KOPS
  26. C***********************************************************************
  27. C***********************************************************************
  28. C HISTORIQUE : 26/10/98 : prise en compte du cas particulier
  29. C où A est diagonale. C'est en fait assez pénible
  30. C car on utilise le CHPOINT comme structure de
  31. C données pour la matrice A et les vecteurs X,B,CLIM
  32. C HISTORIQUE : 09/02/99 : on peut utiliser le préconditionnement d'une
  33. C autre matrice que celle que l'on inverse
  34. C HISTORIQUE : 01/06/99 : on modifie la partie résolution itérative
  35. C en effet, lors de l'imposition des CL. de
  36. C Dirichlet, on changeait les valeurs de la
  37. C matrice Morse.
  38. C Ceci n'est pas bon lorsqu'on veut utiliser la
  39. C matrice assemblée pour plusieurs pas de temps.
  40. C On travaille maintenant sur une copie.
  41. C HISTORIQUE : 20/12/99 : reprogrammation de l'assemblage
  42. C Le nouvel assemblage n'est, pour l'instant effectif que
  43. C lorsqu'une méthode itérative est sélectionnée (-> fabrication
  44. C d'une matrice Morse). Le nouvel assemblage est plus performant
  45. C (temps de calcul, utilisation de la mémoire) et plus général (cas
  46. C particuliers à peu près supprimés) que le précédent.
  47. C HISTORIQUE : 05/01/00 : On ne supprime plus les 0.D0 de la matrice
  48. C assemblée (cf. clmors appelé par melim). Ceci pour ne plus avoir
  49. C perte éventuelle de symétrie de la matrice assemblée. Cela permet
  50. C aussi de ne plus dupliquer la matrice assemblée.
  51. C HISTORIQUE : 13/01/00 : Mise en conformité du solveur direct avec le
  52. C nouvel assemblage.
  53. C HISTORIQUE : 22/03/00 : Rajout des préconditionneurs ILUT
  54. C HISTORIQUE : 13/04/00 : Séparation Lecture des données
  55. C Ecriture des résultats (kres2)
  56. C Assemblage kres3
  57. C Méthode directe kres4
  58. C Méthodes itératives kres5
  59. C HISTORIQUE : 06/04/04 : Renumérotation (double mult.)
  60. C HISTORIQUE : 06/04/04 : Scaling
  61. C HISTORIQUE : 08/04/04 : ajout ILUTP
  62. C HISTORIQUE : 09/02/06 : ajout nb prod matrice vecteur (NMATVEC)
  63. C simplification du code
  64. C HISTORIQUE : 22/02/06 : Gros nettoyage au niveau de l'entrée des
  65. C arguments (Nouvelle syntaxe)
  66. C HISTORIQUE : 08/2011 : En vue de la suppression de l'objet MATRIK
  67. C on utilise l'assemblage de RESOU.
  68. C HISTORIQUE : 04/2019 : remplacement de NOEL par NELIM
  69. C Idéalement, il faudrait reprendre ce que Pierre a fait dans
  70. C RESOU dans les fiches 10[0,1]?? et avec MREM.En vue de la
  71. C suppression de l'objet MATRIK
  72. C
  73. C***********************************************************************
  74. C Prière de PRENDRE LE TEMPS de compléter les commentaires
  75. C en cas de modification de ce sous-programme afin de faciliter
  76. C la maintenance !
  77. C***********************************************************************
  78.  
  79. -INC PPARAM
  80. -INC CCOPTIO
  81. -INC SMLREEL
  82. POINTEUR LRES.MLREEL
  83. -INC SMLENTI
  84. POINTEUR LNMV.MLENTI
  85. -INC SMCHPOI
  86. POINTEUR KCLIM.MCHPOI
  87. POINTEUR KSMBR.MCHPOI
  88. POINTEUR MCHINI.MCHPOI
  89. POINTEUR MCHSOL.MCHPOI
  90. -INC SMTABLE
  91. POINTEUR MTINV.MTABLE
  92. POINTEUR KTIME.MTABLE
  93. DIMENSION ITTIME(4)
  94. CHARACTER*8 CHARI
  95. * MATRIK est la matrice à inverser
  96. * MAPREC est la matrice dont on utilise le préconditionnement
  97. * MATASS est la matrice dont on utilise l'assemblage
  98. * pour préconditionner celui de MATRIK
  99. POINTEUR MAPREC.MATRIK
  100. POINTEUR MATASS.MATRIK
  101. *STAT -INC SMSTAT
  102. C
  103. CHARACTER*8 TYPE
  104. * Paramètre m du GMRES(m)
  105. INTEGER RESTRT
  106. INTEGER ITER
  107. REAL*8 BRTOL,RESID,RXMILU,RXILUP
  108. *
  109. REAL*8 XLFIL,XDTOL
  110. INTEGER KPREC
  111. INTEGER NMATRI
  112. INTEGER IP,KTYPI
  113. INTEGER IMPINV,IIMPR
  114. CHARACTER*4 MRENU,MMULAG
  115. LOGICAL LRIG,LTIME,LDETR,LDEPE,LASRIG,LDMULT,LOGII
  116. INTEGER IMPR,IRET
  117. C
  118. C Lecture des arguments et mise à défaut des optionnels ()
  119. C
  120. C MATRIK : La matrice lue en entrée au format MATRIK
  121. C MTINV : L'éventuelle table d'inversion (obsolète)
  122. C IMPR : Niveau d'impression solveur direct
  123. C KCLIM : Chpoint éventuel de conditions aux limites de Dirichlet
  124. C KSMBR : Chpoint second membre
  125. C KTYPI : Type de méthode de résolution
  126. C MATASS : Matrice utilisée pour préconditionner l'assemblage
  127. C MAPREC : Matrice utilisée pour préconditionner la construction du
  128. C préconditionneur
  129. C MRENU : Type de renumérotation
  130. C MMULAG : Placement des multiplicateurs de Lagrange
  131. C ISCAL : Scaling de la matrice
  132. C IOUBL : Oubli des matrices élémentaires ?
  133. C IMPINV : Niveau d'impression solveur itératif
  134. C MCHINI : Chpoint estimation de l'inconnue
  135. C ITER : Nombre maxi d'itérations à effectuer
  136. C RESID : Norme L2 maxi du résidu
  137. C BRTOL : Breakdown tolerance pour les méthodes de type Bi-CG
  138. C IRSTRT : Paramètre m de redémarrage pour GMRES
  139. C KPREC : Type du préconditionneur
  140. C RXMILU : Paramètre de relaxation pour MILU(0)
  141. C RXILUP : Paramètre de filtre pour ILU(0)-PV
  142. C XLFIL : Paramètre de remplissage pour les préconditionneurs ILUT
  143. C XDTOL : Drop tolerance pour les préconditionneurs ILUT
  144. C XSPIV : Sensibilité du pivoting pour les préconditionneurs ILUTP
  145. C LBCG : le l dans BicgStab(l)
  146. C ICALRS : façon de calculer le résidu
  147. C METASS : Méthode d'assemblage
  148. C LTIME : construit une table avec des statistiques temporelles
  149. C si vrai
  150. C LDEPE : élimine les dépendances si VRAI
  151. C et matrice d'entrée RIGIDITE
  152. C IDDOT : 0 => utilisation du produit scalaire normal dans les
  153. C méthodes itératives
  154. C 1 => utilisation du produit scalaire compensé
  155. * IMPR=4
  156.  
  157. IVALI=0
  158. XVALI=0.D0
  159. LOGII=.FALSE.
  160. IRETI=0
  161. XVALR=0.D0
  162. *inutile IOBRE=0
  163. IRETR=0
  164.  
  165. IMPR=0
  166. * WRITE(IOIMP,*) 'coucou kres2'
  167. *
  168. *STAT CALL INMSTA(MSTAT,IMPR)
  169. *
  170. CALL PRKRES(MATRIK,MTINV,IMPR,KCLIM,KSMBR,KTYPI,MATASS,MAPREC,
  171. $ MRENU,MMULAG,ISCAL,INORMU,IOUBL,IMPINV,MCHINI,ITER,RESID
  172. $ ,BRTOL,IRSTRT,KPREC,RXMILU,RXILUP,XLFIL,XDTOL,XSPIV,LBCG
  173. $ ,ICALRS,METASS,LTIME,LDEPE,MRIGID,IDDOT,IMVEC,IRET)
  174. IF (IRET.NE.0) GOTO 9999
  175. IMPR=MAX(IMPR,IMPINV)
  176. *
  177. * Préparation de la matrice et du second membre
  178. * suivant les cas
  179. *
  180. * LASRIG=.TRUE. on utilise l'assemblage de RESOU
  181. * LDMULT=.TRUE. on dédouble les multiplicateurs de Lagrange
  182. * LDEPE=.TRUE. On fait l'élimination des relations
  183.  
  184. LASRIG=(METASS.EQ.6)
  185. * Pour l'instant, il faut toujours dédoubler les multiplicateurs
  186. * quand on assemble avec l'assemblage de RESO car le traitement des
  187. * multiplicateurs dans ldmt1 l'impose (simple multiplicateur non
  188. * prévu)
  189. LDMULT=LASRIG
  190. MRIGI0=MRIGID
  191. * Nouveau separm gère la non élimination avec NOEL
  192. IF (LDEPE) then
  193. * noel=0
  194. nelim=1
  195. ELSE
  196. * noel=1
  197. nelim=0
  198. ENDIF
  199. *dbg write(ioimp,*) 'LASRIG=',LASRIG
  200. *dbg write(ioimp,*) 'LDMULT=',LDMULT
  201. *dbg write(ioimp,*) 'LDEPE=',LDEPE,' NELIM=',NELIM
  202.  
  203. IF (MRIGID.NE.0) THEN
  204. *old IF (LDEPE) THEN
  205. KSMBR0=KSMBR
  206. *old CALL KRES6(MRIGID,KSMBR,IDEPE,
  207. *old $ MRIGIC,KSMBRC,KSMBR0,KSMBR1)
  208. CALL KRES6(MRIGID,KSMBR,LDMULT,NELIM,
  209. $ MRIGIC,KSMBRC,KSMBR1)
  210. IF (IERR.NE.0) RETURN
  211. MRIGID=MRIGIC
  212. KSMBR=KSMBRC
  213. *old ENDIF
  214. IF (LASRIG) THEN
  215. * Gestion de la normalisation
  216. NORICO=NORINC
  217. NORIDO=NORIND
  218. IF (ISCAL.EQ.0) THEN
  219. NORINC=0
  220. NORIND=0
  221. ELSE
  222. NORINC=-1
  223. NORIND=0
  224. ENDIF
  225. * Gestion de la renumérotation
  226. NUCROO=NUCROU
  227. IF (MRENU.EQ.'RIEN') THEN
  228. NUCROU=-1
  229. * La renumérotation sera en fait Reverse Cuthill-McKee dans NUMOPT
  230. ELSEIF (MRENU.EQ.'SLOA') THEN
  231. NUCROU=1
  232. * La renumérotation sera en fait Nested Dissection dans NUMOPT
  233. ELSEIF (MRENU.EQ.'GIPR'.OR.MRENU.EQ.'GIBA') THEN
  234. NUCROU=2
  235. ELSE
  236. WRITE(IOIMP,*) 'MRENU=',MRENU
  237. CALL ERREUR(5)
  238. RETURN
  239. ENDIF
  240. CALL KRES8(MRIGID,KSMBR,INORMU,
  241. $ KTYPI,ITER,RESID,ICALRS,IRSTRT,LBCG,BRTOL,IDDOT,IMVEC,
  242. $ KPREC,RXMILU,RXILUP,XLFIL,XDTOL,XSPIV,
  243. $ KTIME,LTIME,
  244. $ MCHSOL,LRES,LNMV,ICVG,IMPR)
  245. IF (IERR.NE.0) RETURN
  246. * Gestion de la normalisation
  247. NORINC=NORICO
  248. NORIND=NORIDO
  249. NUCROU=NUCROO
  250. IF (LTIME) CALL ECROBJ('TABLE ',KTIME)
  251. IF (MTINV.NE.0) THEN
  252. CALL ECMO(MTINV,'CVGOK','ENTIER',ICVG)
  253. IF (LRES.NE.0) CALL ECMO(MTINV,'CONVINV','LISTREEL',LRES)
  254. IF (LNMV.NE.0) CALL ECMO(MTINV,'NMATVEC','LISTENTI',LNMV)
  255. ENDIF
  256. * On décondense si nécessaire
  257. *
  258. * write (6,*) ' resou - mchsol '
  259. * call ecchpo(mchsol,0)
  260. * call mucpri(mchsol,mrigid,iresi)
  261. * write (6,*) ' kres - iresi '
  262. * call ecchpo(iresi,0)
  263. * WRITE(IOIMP,*) 'Avant KRES7'
  264. *old IF (MRIGI0.NE.0.AND.LDEPE) THEN
  265. IF (MRIGI0.NE.0) THEN
  266. MSOLC=MCHSOL
  267. IDTARG=1
  268. * On détruit MSOLC dans KRES7
  269. CALL KRES7(MSOLC,MRIGI0,KSMBR0,KSMBR1,IDTARG,
  270. $ MCHSOL)
  271. IF (IERR.NE.0) RETURN
  272. ENDIF
  273. CALL ECROBJ('CHPOINT ',MCHSOL)
  274. RETURN
  275. ELSE
  276. CALL ECROBJ('RIGIDITE',MRIGID)
  277. CALL RIMA
  278. IF (IERR.NE.0) GOTO 9999
  279. CALL MACHI2(1)
  280. IF (IERR.NE.0) GOTO 9999
  281. CALL LIROBJ('MATRIK',MATRIK,1,IRET)
  282. IF(IRET.EQ.0) GOTO 9999
  283.  
  284. * Changement des noms d'inconnues du second membre
  285. IF (KSMBR.NE.0) THEN
  286. CALL ECROBJ('CHPOINT ',KSMBR)
  287. CALL MACHI2(1)
  288. CALL LIROBJ('CHPOINT ',KSMBR,1,IRET)
  289. IF (IRET.EQ.0) GOTO 9999
  290. ENDIF
  291. ENDIF
  292. * write (6,*) ' le vecteur 2'
  293. * call ecchpo(ksmbr,0)
  294. * write (6,*) ' la matrice 2'
  295. * call ecrent(5)
  296. * call ecmatk(matrik)
  297. ENDIF
  298. *
  299. SEGACT MATRIK
  300. NMATRI=IRIGEL(/2)
  301. IF(NMATRI.EQ.0)THEN
  302. C% Résolution impossible : la matrice de RIGIDITE est vide
  303. CALL ERREUR(727)
  304. RETURN
  305. ENDIF
  306. SEGDES MATRIK
  307. IF (MATASS.EQ.0) MATASS=MATRIK
  308. IF (MAPREC.EQ.0) MAPREC=MATRIK
  309. * WRITE(IOIMP,*) 'Sortie de prkres'
  310. * WRITE(IOIMP,*) 'IOUBL=',IOUBL
  311. C
  312. IF (LTIME) THEN
  313. CALL CRTABL(KTIME)
  314. call timespv(ittime,oothrd)
  315. ITI1=(ITTIME(1)+ITTIME(2))/10
  316. ELSE
  317. KTIME=0
  318. ENDIF
  319. *STAT CALL PRMSTA('Lectures',MSTAT,IMPR)
  320. *
  321. C
  322. C Assemblage proprement dit
  323. C
  324. IIMPR=0
  325. CALL KRES3(MATRIK,MATASS,MRENU,MMULAG,METASS,
  326. $ KTIME,LTIME,
  327. $ IIMPR,IRET)
  328. IF (IRET.NE.0) GOTO 9999
  329. *! WRITE(IOIMP,*) 'Aprés assemblage'
  330. *! CALL ECRENT(5)
  331. *! CALL ECROBJ('MATRIK ',MATRIK)
  332. *! CALL PRLIST
  333. IF (LTIME) THEN
  334. call timespv(ittime,oothrd)
  335. ITI2=(ITTIME(1)+ITTIME(2))/10
  336. ENDIF
  337. *STAT CALL PRMSTA('Assemblage',MSTAT,IMPR)
  338. *
  339. * "Oubli" des valeurs des matrice élémentaires
  340. * On met les tableaux de LIZAFM à 0 => à MENAGE de les supprimmer
  341. * si besoin est.
  342. *
  343. IOUBD=MOD(IOUBL,10)
  344. *! WRITE(IOIMP,*) 'IOUBD=',IOUBD
  345. IF (IOUBD.EQ.1) THEN
  346. CALL OUBIMA(MATRIK,IMPR,IRET)
  347. IF (IRET.NE.0) GOTO 9999
  348. IF (IMPR.GT.2) THEN
  349. WRITE(IOIMP,*) 'Oubli des mat. elem.'
  350. ENDIF
  351. ELSEIF (IOUBD.EQ.2) THEN
  352. call ooohor(0)
  353. SEGACT MATRIK*MOD
  354. LDETR=.FALSE.
  355. NMATE=IRIGEL(/2)
  356. DO IMATE=1,NMATE
  357. IMATRI=IRIGEL(4,IMATE)
  358. SEGACT IMATRI*MOD
  359. NSOUM =LIZAFM(/1)
  360. NTOTIN=LIZAFM(/2)
  361. DO ITOTIN=1,NTOTIN
  362. DO ISOUM=1,NSOUM
  363. IZAFM=LIZAFM(ISOUM,ITOTIN)
  364. IF (IZAFM.NE.0) THEN
  365. LDETR=.TRUE.
  366. SEGSUP IZAFM
  367. LIZAFM(ISOUM,ITOTIN)=0
  368. ENDIF
  369. ENDDO
  370. ENDDO
  371. SEGDES IMATRI
  372. ENDDO
  373. IF (IMPR.GT.2.AND.LDETR) THEN
  374. WRITE(IOIMP,*) 'Destruction des mat. elem.'
  375. ENDIF
  376. ELSEIF (IOUBD.NE.0) THEN
  377. WRITE(IOIMP,*) 'IOUBL=',IOUBL, ' non prevu'
  378. GOTO 9999
  379. ENDIF
  380. *STAT CALL PRMSTA('Oubli',MSTAT,IMPR)
  381. *! WRITE(IOIMP,*) 'Aprés oubli'
  382. C
  383. C Méthode directe
  384. C
  385. IF (KTYPI.EQ.1) THEN
  386. CALL KRES4(MATRIK,KCLIM,KSMBR,
  387. $ ISCAL,
  388. $ MCHSOL,
  389. $ IMPR,IRET)
  390. IF (IRET.NE.0) GOTO 9999
  391. *STAT CALL PRMSTA('Methode directe',MSTAT,IMPR)
  392. C
  393. C Méthodes itératives
  394. C
  395. ELSE
  396. C
  397. CALL KRES5(MATRIK,KCLIM,KSMBR,KTYPI,
  398. $ MCHINI,ITER,RESID,
  399. $ BRTOL,IRSTRT,LBCG,ICALRS,
  400. $ MAPREC,KPREC,
  401. $ RXMILU,RXILUP,XLFIL,XDTOL,XSPIV,
  402. $ ISCAL,
  403. $ KTIME,LTIME,IDDOT,IMVEC,
  404. $ MCHSOL,LRES,LNMV,ICVG,
  405. $ IMPR,IRET)
  406. IF (IRET.NE.0) GOTO 9999
  407. *STAT CALL PRMSTA('Methode iterative',MSTAT,IMPR)
  408. IF (MTINV.NE.0) THEN
  409. CALL ECMO(MTINV,'CVGOK','ENTIER',ICVG)
  410. CALL ECMO(MTINV,'CONVINV','LISTREEL',LRES)
  411. CALL ECMO(MTINV,'NMATVEC','LISTENTI',LNMV)
  412. ENDIF
  413. ENDIF
  414. IF (LTIME) THEN
  415. call timespv(ittime,oothrd)
  416. ITI3=(ITTIME(1)+ITTIME(2))/10
  417. CHARI='ASS+RENU'
  418. CALL ECCTAB(KTIME,'MOT ',IVALI,XVALI,CHARI,LOGII,IRETI,
  419. $ 'ENTIER ',ITI2-ITI1,XVALR,CHARR,LOGIR,IRETR)
  420. CHARI='PRE+RESO'
  421. CALL ECCTAB(KTIME,'MOT ',IVALI,XVALI,CHARI,LOGII,IRETI,
  422. $ 'ENTIER ',ITI3-ITI2,XVALR,CHARR,LOGIR,IRETR)
  423. CHARI='TOTAL '
  424. CALL ECCTAB(KTIME,'MOT ',IVALI,XVALI,CHARI,LOGII,IRETI,
  425. $ 'ENTIER ',ITI3-ITI1,XVALR,CHARR,LOGIR,IRETR)
  426. SEGDES KTIME
  427. CALL ECROBJ('TABLE ',KTIME)
  428. ENDIF
  429. IOUBE=IOUBL/10
  430. *! WRITE(IOIMP,*) 'IOUBE=',IOUBE
  431. IF (IOUBE.GE.1) THEN
  432. call ooohor(0)
  433. SEGACT MATRIK*MOD
  434. IF (IOUBE.EQ.2) THEN
  435. PMORS=KIDMAT(4)
  436. IF (PMORS.NE.0) THEN
  437. IF (IMPR.GT.2) THEN
  438. WRITE(IOIMP,*) 'Destruction du profil morse'
  439. ENDIF
  440. SEGSUP PMORS
  441. KIDMAT(4)=0
  442. ENDIF
  443. ENDIF
  444. IZA=KIDMAT(5)
  445. IF (IZA.NE.0) THEN
  446. IF (IMPR.GT.2) THEN
  447. WRITE(IOIMP,*) 'Destruction des valeurs morses'
  448. ENDIF
  449. SEGSUP IZA
  450. KIDMAT(5)=0
  451. ENDIF
  452. PMORS=KIDMAT(6)
  453. IF (PMORS.NE.0) THEN
  454. IF (IMPR.GT.2) THEN
  455. WRITE(IOIMP,*) 'Destruction du profil du precon'
  456. ENDIF
  457. SEGSUP PMORS
  458. KIDMAT(6)=0
  459. ENDIF
  460. IZA=KIDMAT(7)
  461. IF (IZA.NE.0) THEN
  462. IF (IMPR.GT.2) THEN
  463. WRITE(IOIMP,*) 'Destruction des valeurs du precon'
  464. ENDIF
  465. SEGSUP IZA
  466. KIDMAT(7)=0
  467. ENDIF
  468. SEGDES MATRIK
  469. ELSEIF (IOUBE.NE.0) THEN
  470. WRITE(IOIMP,*) 'IOUBL=',IOUBL, ' non prevu'
  471. GOTO 9999
  472. ENDIF
  473. *
  474. * On décondense si nécessaire
  475. *
  476. * write (6,*) ' resou - mchsol '
  477. * call ecchpo(mchsol,0)
  478. * call mucpri(mchsol,mrigid,iresi)
  479. * write (6,*) ' kres - iresi '
  480. * call ecchpo(iresi,0)
  481. * WRITE(IOIMP,*) 'Avant KRES7'
  482. IF (MRIGI0.NE.0.AND.LDEPE) THEN
  483. MSOLC=MCHSOL
  484. * CALL KRES7(MRIGID,IDEPE,KSMBR0,KSMBR1,
  485. * $ MSOLC,
  486. * $ MCHSOL)
  487. CALL KRES7(MSOLC,MRIGI0,KSMBR0,KSMBR1,1,
  488. $ MCHSOL)
  489. IF (IERR.NE.0) RETURN
  490. ENDIF
  491. CALL ECROBJ('CHPOINT ',MCHSOL)
  492. *STAT CALL SUMSTA(MSTAT,IMPR)
  493. *
  494. * Normal termination
  495. *
  496. RETURN
  497. *
  498. * Format handling
  499. *
  500. *
  501. * Error handling
  502. *
  503. 9999 CONTINUE
  504. WRITE(IOIMP,*) 'An error was detected in kres2.eso'
  505. * 153 2
  506. * Opération illicite dans ce contexte
  507. CALL ERREUR(153)
  508. RETURN
  509. *
  510. * End of KRES2
  511. *
  512. END
  513.  
  514.  
  515.  

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