Télécharger kres2.eso

Retour à la liste

Numérotation des lignes :

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

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