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. -INC CCOPTIO
  79. -INC SMLREEL
  80. POINTEUR LRES.MLREEL
  81. -INC SMLENTI
  82. POINTEUR LNMV.MLENTI
  83. -INC SMCHPOI
  84. POINTEUR KCLIM.MCHPOI
  85. POINTEUR KSMBR.MCHPOI
  86. POINTEUR MCHINI.MCHPOI
  87. POINTEUR MCHSOL.MCHPOI
  88. -INC SMTABLE
  89. POINTEUR MTINV.MTABLE
  90. POINTEUR KTIME.MTABLE
  91. DIMENSION ITTIME(4)
  92. CHARACTER*8 CHARI
  93. * MATRIK est la matrice à inverser
  94. * MAPREC est la matrice dont on utilise le préconditionnement
  95. * MATASS est la matrice dont on utilise l'assemblage
  96. * pour préconditionner celui de MATRIK
  97. POINTEUR MAPREC.MATRIK
  98. POINTEUR MATASS.MATRIK
  99. *STAT -INC SMSTAT
  100. C
  101. CHARACTER*8 TYPE
  102. * Paramètre m du GMRES(m)
  103. INTEGER RESTRT
  104. INTEGER ITER
  105. REAL*8 BRTOL,RESID,RXMILU,RXILUP
  106. *
  107. REAL*8 XLFIL,XDTOL
  108. INTEGER KPREC
  109. INTEGER NMATRI
  110. INTEGER IP,KTYPI
  111. INTEGER IMPINV,IIMPR
  112. CHARACTER*4 MRENU,MMULAG
  113. LOGICAL LRIG,LTIME,LDETR,LDEPE,LASRIG,LDMULT,LOGII
  114. INTEGER IMPR,IRET
  115. C
  116. C Lecture des arguments et mise à défaut des optionnels ()
  117. C
  118. C MATRIK : La matrice lue en entrée au format MATRIK
  119. C MTINV : L'éventuelle table d'inversion (obsolète)
  120. C IMPR : Niveau d'impression solveur direct
  121. C KCLIM : Chpoint éventuel de conditions aux limites de Dirichlet
  122. C KSMBR : Chpoint second membre
  123. C KTYPI : Type de méthode de résolution
  124. C MATASS : Matrice utilisée pour préconditionner l'assemblage
  125. C MAPREC : Matrice utilisée pour préconditionner la construction du
  126. C préconditionneur
  127. C MRENU : Type de renumérotation
  128. C MMULAG : Placement des multiplicateurs de Lagrange
  129. C ISCAL : Scaling de la matrice
  130. C IOUBL : Oubli des matrices élémentaires ?
  131. C IMPINV : Niveau d'impression solveur itératif
  132. C MCHINI : Chpoint estimation de l'inconnue
  133. C ITER : Nombre maxi d'itérations à effectuer
  134. C RESID : Norme L2 maxi du résidu
  135. C BRTOL : Breakdown tolerance pour les méthodes de type Bi-CG
  136. C IRSTRT : Paramètre m de redémarrage pour GMRES
  137. C KPREC : Type du préconditionneur
  138. C RXMILU : Paramètre de relaxation pour MILU(0)
  139. C RXILUP : Paramètre de filtre pour ILU(0)-PV
  140. C XLFIL : Paramètre de remplissage pour les préconditionneurs ILUT
  141. C XDTOL : Drop tolerance pour les préconditionneurs ILUT
  142. C XSPIV : Sensibilité du pivoting pour les préconditionneurs ILUTP
  143. C LBCG : le l dans BicgStab(l)
  144. C ICALRS : façon de calculer le résidu
  145. C METASS : Méthode d'assemblage
  146. C LTIME : construit une table avec des statistiques temporelles
  147. C si vrai
  148. C LDEPE : élimine les dépendances si VRAI
  149. C et matrice d'entrée RIGIDITE
  150. C IDDOT : 0 => utilisation du produit scalaire normal dans les
  151. C méthodes itératives
  152. C 1 => utilisation du produit scalaire compensé
  153. * IMPR=4
  154.  
  155. IVALI=0
  156. XVALI=0.D0
  157. LOGII=.FALSE.
  158. IRETI=0
  159. XVALR=0.D0
  160. *inutile IOBRE=0
  161. IRETR=0
  162.  
  163. IMPR=0
  164. * WRITE(IOIMP,*) 'coucou kres2'
  165. *
  166. *STAT CALL INMSTA(MSTAT,IMPR)
  167. *
  168. CALL PRKRES(MATRIK,MTINV,IMPR,KCLIM,KSMBR,KTYPI,MATASS,MAPREC,
  169. $ MRENU,MMULAG,ISCAL,INORMU,IOUBL,IMPINV,MCHINI,ITER,RESID
  170. $ ,BRTOL,IRSTRT,KPREC,RXMILU,RXILUP,XLFIL,XDTOL,XSPIV,LBCG
  171. $ ,ICALRS,METASS,LTIME,LDEPE,MRIGID,IDDOT,IMVEC,IRET)
  172. IF (IRET.NE.0) GOTO 9999
  173. IMPR=MAX(IMPR,IMPINV)
  174. *
  175. * Préparation de la matrice et du second membre
  176. * suivant les cas
  177. *
  178. * LASRIG=.TRUE. on utilise l'assemblage de RESOU
  179. * LDMULT=.TRUE. on dédouble les multiplicateurs de Lagrange
  180. * LDEPE=.TRUE. On fait l'élimination des relations
  181.  
  182. LASRIG=(METASS.EQ.6)
  183. * Pour l'instant, il faut toujours dédoubler les multiplicateurs
  184. * quand on assemble avec l'assemblage de RESO car le traitement des
  185. * multiplicateurs dans ldmt1 l'impose (simple multiplicateur non
  186. * prévu)
  187. LDMULT=LASRIG
  188. MRIGI0=MRIGID
  189. * Nouveau separm gère la non élimination avec NOEL
  190. IF (LDEPE) then
  191. * noel=0
  192. nelim=1
  193. ELSE
  194. * noel=1
  195. nelim=0
  196. ENDIF
  197. *dbg write(ioimp,*) 'LASRIG=',LASRIG
  198. *dbg write(ioimp,*) 'LDMULT=',LDMULT
  199. *dbg write(ioimp,*) 'LDEPE=',LDEPE,' NELIM=',NELIM
  200.  
  201. IF (MRIGID.NE.0) THEN
  202. *old IF (LDEPE) THEN
  203. KSMBR0=KSMBR
  204. *old CALL KRES6(MRIGID,KSMBR,IDEPE,
  205. *old $ MRIGIC,KSMBRC,KSMBR0,KSMBR1)
  206. CALL KRES6(MRIGID,KSMBR,LDMULT,NELIM,
  207. $ MRIGIC,KSMBRC,KSMBR1)
  208. IF (IERR.NE.0) RETURN
  209. MRIGID=MRIGIC
  210. KSMBR=KSMBRC
  211. *old ENDIF
  212. IF (LASRIG) THEN
  213. * Gestion de la normalisation
  214. NORICO=NORINC
  215. NORIDO=NORIND
  216. IF (ISCAL.EQ.0) THEN
  217. NORINC=0
  218. NORIND=0
  219. ELSE
  220. NORINC=-1
  221. NORIND=0
  222. ENDIF
  223. * Gestion de la renumérotation
  224. NUCROO=NUCROU
  225. IF (MRENU.EQ.'RIEN') THEN
  226. NUCROU=-1
  227. * La renumérotation sera en fait Reverse Cuthill-McKee dans NUMOPT
  228. ELSEIF (MRENU.EQ.'SLOA') THEN
  229. NUCROU=1
  230. * La renumérotation sera en fait Nested Dissection dans NUMOPT
  231. ELSEIF (MRENU.EQ.'GIPR'.OR.MRENU.EQ.'GIBA') THEN
  232. NUCROU=2
  233. ELSE
  234. WRITE(IOIMP,*) 'MRENU=',MRENU
  235. CALL ERREUR(5)
  236. RETURN
  237. ENDIF
  238. CALL KRES8(MRIGID,KSMBR,INORMU,
  239. $ KTYPI,ITER,RESID,ICALRS,IRSTRT,LBCG,BRTOL,IDDOT,IMVEC,
  240. $ KPREC,RXMILU,RXILUP,XLFIL,XDTOL,XSPIV,
  241. $ KTIME,LTIME,
  242. $ MCHSOL,LRES,LNMV,ICVG,IMPR)
  243. IF (IERR.NE.0) RETURN
  244. * Gestion de la normalisation
  245. NORINC=NORICO
  246. NORIND=NORIDO
  247. NUCROU=NUCROO
  248. IF (LTIME) CALL ECROBJ('TABLE ',KTIME)
  249. IF (MTINV.NE.0) THEN
  250. CALL ECMO(MTINV,'CVGOK','ENTIER',ICVG)
  251. IF (LRES.NE.0) CALL ECMO(MTINV,'CONVINV','LISTREEL',LRES)
  252. IF (LNMV.NE.0) CALL ECMO(MTINV,'NMATVEC','LISTENTI',LNMV)
  253. ENDIF
  254. * On décondense si nécessaire
  255. *
  256. * write (6,*) ' resou - mchsol '
  257. * call ecchpo(mchsol,0)
  258. * call mucpri(mchsol,mrigid,iresi)
  259. * write (6,*) ' kres - iresi '
  260. * call ecchpo(iresi,0)
  261. * WRITE(IOIMP,*) 'Avant KRES7'
  262. *old IF (MRIGI0.NE.0.AND.LDEPE) THEN
  263. IF (MRIGI0.NE.0) THEN
  264. MSOLC=MCHSOL
  265. IDTARG=1
  266. * On détruit MSOLC dans KRES7
  267. CALL KRES7(MSOLC,MRIGI0,KSMBR0,KSMBR1,IDTARG,
  268. $ MCHSOL)
  269. IF (IERR.NE.0) RETURN
  270. ENDIF
  271. CALL ECROBJ('CHPOINT ',MCHSOL)
  272. RETURN
  273. ELSE
  274. CALL ECROBJ('RIGIDITE',MRIGID)
  275. CALL RIMA
  276. IF (IERR.NE.0) GOTO 9999
  277. CALL MACHI2(1)
  278. IF (IERR.NE.0) GOTO 9999
  279. CALL LIROBJ('MATRIK',MATRIK,1,IRET)
  280. IF(IRET.EQ.0) GOTO 9999
  281.  
  282. * Changement des noms d'inconnues du second membre
  283. IF (KSMBR.NE.0) THEN
  284. CALL ECROBJ('CHPOINT ',KSMBR)
  285. CALL MACHI2(1)
  286. CALL LIROBJ('CHPOINT ',KSMBR,1,IRET)
  287. IF (IRET.EQ.0) GOTO 9999
  288. ENDIF
  289. ENDIF
  290. * write (6,*) ' le vecteur 2'
  291. * call ecchpo(ksmbr,0)
  292. * write (6,*) ' la matrice 2'
  293. * call ecrent(5)
  294. * call ecmatk(matrik)
  295. ENDIF
  296. *
  297. SEGACT MATRIK
  298. NMATRI=IRIGEL(/2)
  299. IF(NMATRI.EQ.0)THEN
  300. C% Résolution impossible : la matrice de RIGIDITE est vide
  301. CALL ERREUR(727)
  302. RETURN
  303. ENDIF
  304. SEGDES MATRIK
  305. IF (MATASS.EQ.0) MATASS=MATRIK
  306. IF (MAPREC.EQ.0) MAPREC=MATRIK
  307. * WRITE(IOIMP,*) 'Sortie de prkres'
  308. * WRITE(IOIMP,*) 'IOUBL=',IOUBL
  309. C
  310. IF (LTIME) THEN
  311. CALL CRTABL(KTIME)
  312. call timespv(ittime,oothrd)
  313. ITI1=(ITTIME(1)+ITTIME(2))/10
  314. ELSE
  315. KTIME=0
  316. ENDIF
  317. *STAT CALL PRMSTA('Lectures',MSTAT,IMPR)
  318. *
  319. C
  320. C Assemblage proprement dit
  321. C
  322. IIMPR=0
  323. CALL KRES3(MATRIK,MATASS,MRENU,MMULAG,METASS,
  324. $ KTIME,LTIME,
  325. $ IIMPR,IRET)
  326. IF (IRET.NE.0) GOTO 9999
  327. *! WRITE(IOIMP,*) 'Aprés assemblage'
  328. *! CALL ECRENT(5)
  329. *! CALL ECROBJ('MATRIK ',MATRIK)
  330. *! CALL PRLIST
  331. IF (LTIME) THEN
  332. call timespv(ittime,oothrd)
  333. ITI2=(ITTIME(1)+ITTIME(2))/10
  334. ENDIF
  335. *STAT CALL PRMSTA('Assemblage',MSTAT,IMPR)
  336. *
  337. * "Oubli" des valeurs des matrice élémentaires
  338. * On met les tableaux de LIZAFM à 0 => à MENAGE de les supprimmer
  339. * si besoin est.
  340. *
  341. IOUBD=MOD(IOUBL,10)
  342. *! WRITE(IOIMP,*) 'IOUBD=',IOUBD
  343. IF (IOUBD.EQ.1) THEN
  344. CALL OUBIMA(MATRIK,IMPR,IRET)
  345. IF (IRET.NE.0) GOTO 9999
  346. IF (IMPR.GT.2) THEN
  347. WRITE(IOIMP,*) 'Oubli des mat. elem.'
  348. ENDIF
  349. ELSEIF (IOUBD.EQ.2) THEN
  350. call ooohor(0)
  351. SEGACT MATRIK*MOD
  352. LDETR=.FALSE.
  353. NMATE=IRIGEL(/2)
  354. DO IMATE=1,NMATE
  355. IMATRI=IRIGEL(4,IMATE)
  356. SEGACT IMATRI*MOD
  357. NSOUM =LIZAFM(/1)
  358. NTOTIN=LIZAFM(/2)
  359. DO ITOTIN=1,NTOTIN
  360. DO ISOUM=1,NSOUM
  361. IZAFM=LIZAFM(ISOUM,ITOTIN)
  362. IF (IZAFM.NE.0) THEN
  363. LDETR=.TRUE.
  364. SEGSUP IZAFM
  365. LIZAFM(ISOUM,ITOTIN)=0
  366. ENDIF
  367. ENDDO
  368. ENDDO
  369. SEGDES IMATRI
  370. ENDDO
  371. IF (IMPR.GT.2.AND.LDETR) THEN
  372. WRITE(IOIMP,*) 'Destruction des mat. elem.'
  373. ENDIF
  374. ELSEIF (IOUBD.NE.0) THEN
  375. WRITE(IOIMP,*) 'IOUBL=',IOUBL, ' non prevu'
  376. GOTO 9999
  377. ENDIF
  378. *STAT CALL PRMSTA('Oubli',MSTAT,IMPR)
  379. *! WRITE(IOIMP,*) 'Aprés oubli'
  380. C
  381. C Méthode directe
  382. C
  383. IF (KTYPI.EQ.1) THEN
  384. CALL KRES4(MATRIK,KCLIM,KSMBR,
  385. $ ISCAL,
  386. $ MCHSOL,
  387. $ IMPR,IRET)
  388. IF (IRET.NE.0) GOTO 9999
  389. *STAT CALL PRMSTA('Methode directe',MSTAT,IMPR)
  390. C
  391. C Méthodes itératives
  392. C
  393. ELSE
  394. C
  395. CALL KRES5(MATRIK,KCLIM,KSMBR,KTYPI,
  396. $ MCHINI,ITER,RESID,
  397. $ BRTOL,IRSTRT,LBCG,ICALRS,
  398. $ MAPREC,KPREC,
  399. $ RXMILU,RXILUP,XLFIL,XDTOL,XSPIV,
  400. $ ISCAL,
  401. $ KTIME,LTIME,IDDOT,IMVEC,
  402. $ MCHSOL,LRES,LNMV,ICVG,
  403. $ IMPR,IRET)
  404. IF (IRET.NE.0) GOTO 9999
  405. *STAT CALL PRMSTA('Methode iterative',MSTAT,IMPR)
  406. IF (MTINV.NE.0) THEN
  407. CALL ECMO(MTINV,'CVGOK','ENTIER',ICVG)
  408. CALL ECMO(MTINV,'CONVINV','LISTREEL',LRES)
  409. CALL ECMO(MTINV,'NMATVEC','LISTENTI',LNMV)
  410. ENDIF
  411. ENDIF
  412. IF (LTIME) THEN
  413. call timespv(ittime,oothrd)
  414. ITI3=(ITTIME(1)+ITTIME(2))/10
  415. CHARI='ASS+RENU'
  416. CALL ECCTAB(KTIME,'MOT ',IVALI,XVALI,CHARI,LOGII,IRETI,
  417. $ 'ENTIER ',ITI2-ITI1,XVALR,CHARR,LOGIR,IRETR)
  418. CHARI='PRE+RESO'
  419. CALL ECCTAB(KTIME,'MOT ',IVALI,XVALI,CHARI,LOGII,IRETI,
  420. $ 'ENTIER ',ITI3-ITI2,XVALR,CHARR,LOGIR,IRETR)
  421. CHARI='TOTAL '
  422. CALL ECCTAB(KTIME,'MOT ',IVALI,XVALI,CHARI,LOGII,IRETI,
  423. $ 'ENTIER ',ITI3-ITI1,XVALR,CHARR,LOGIR,IRETR)
  424. SEGDES KTIME
  425. CALL ECROBJ('TABLE ',KTIME)
  426. ENDIF
  427. IOUBE=IOUBL/10
  428. *! WRITE(IOIMP,*) 'IOUBE=',IOUBE
  429. IF (IOUBE.GE.1) THEN
  430. call ooohor(0)
  431. SEGACT MATRIK*MOD
  432. IF (IOUBE.EQ.2) THEN
  433. PMORS=KIDMAT(4)
  434. IF (PMORS.NE.0) THEN
  435. IF (IMPR.GT.2) THEN
  436. WRITE(IOIMP,*) 'Destruction du profil morse'
  437. ENDIF
  438. SEGSUP PMORS
  439. KIDMAT(4)=0
  440. ENDIF
  441. ENDIF
  442. IZA=KIDMAT(5)
  443. IF (IZA.NE.0) THEN
  444. IF (IMPR.GT.2) THEN
  445. WRITE(IOIMP,*) 'Destruction des valeurs morses'
  446. ENDIF
  447. SEGSUP IZA
  448. KIDMAT(5)=0
  449. ENDIF
  450. PMORS=KIDMAT(6)
  451. IF (PMORS.NE.0) THEN
  452. IF (IMPR.GT.2) THEN
  453. WRITE(IOIMP,*) 'Destruction du profil du precon'
  454. ENDIF
  455. SEGSUP PMORS
  456. KIDMAT(6)=0
  457. ENDIF
  458. IZA=KIDMAT(7)
  459. IF (IZA.NE.0) THEN
  460. IF (IMPR.GT.2) THEN
  461. WRITE(IOIMP,*) 'Destruction des valeurs du precon'
  462. ENDIF
  463. SEGSUP IZA
  464. KIDMAT(7)=0
  465. ENDIF
  466. SEGDES MATRIK
  467. ELSEIF (IOUBE.NE.0) THEN
  468. WRITE(IOIMP,*) 'IOUBL=',IOUBL, ' non prevu'
  469. GOTO 9999
  470. ENDIF
  471. *
  472. * On décondense si nécessaire
  473. *
  474. * write (6,*) ' resou - mchsol '
  475. * call ecchpo(mchsol,0)
  476. * call mucpri(mchsol,mrigid,iresi)
  477. * write (6,*) ' kres - iresi '
  478. * call ecchpo(iresi,0)
  479. * WRITE(IOIMP,*) 'Avant KRES7'
  480. IF (MRIGI0.NE.0.AND.LDEPE) THEN
  481. MSOLC=MCHSOL
  482. * CALL KRES7(MRIGID,IDEPE,KSMBR0,KSMBR1,
  483. * $ MSOLC,
  484. * $ MCHSOL)
  485. CALL KRES7(MSOLC,MRIGI0,KSMBR0,KSMBR1,1,
  486. $ MCHSOL)
  487. IF (IERR.NE.0) RETURN
  488. ENDIF
  489. CALL ECROBJ('CHPOINT ',MCHSOL)
  490. *STAT CALL SUMSTA(MSTAT,IMPR)
  491. *
  492. * Normal termination
  493. *
  494. RETURN
  495. *
  496. * Format handling
  497. *
  498. *
  499. * Error handling
  500. *
  501. 9999 CONTINUE
  502. WRITE(IOIMP,*) 'An error was detected in kres2.eso'
  503. * 153 2
  504. * Opération illicite dans ce contexte
  505. CALL ERREUR(153)
  506. RETURN
  507. *
  508. * End of KRES2
  509. *
  510. END
  511.  
  512.  
  513.  

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