Télécharger la plaquette Cast3M


fr en

Télécharger kres8.eso

Retour à la liste

Numérotation des lignes : 

kres8

  1. C KRES8     SOURCE    GOUNAND   25/03/24    21:15:06     12216          
  2.       SUBROUTINE KRES8(MRIGID,KSMBR,INORMU,
  3.      $     KTYPI,ITER,RESID,ICALRS,IRSTRT,LBCG,BRTOL,IDDOT,IMVEC,
  4.      $     IORINC,
  5.      $     KPREC,RXMILU,RXILUP,XLFIL,XDTOL,XSPIV,
  6.      $     KTIME,LTIME,
  7.      $     MCHSOL,LRES,LNMV,ICVG,IMPR)
  8.       IMPLICIT REAL*8 (A-H,O-Z)
  9.       IMPLICIT INTEGER (I-N)
  10. C***********************************************************************
  11. C NOM         : KRES8
  12. C DESCRIPTION : - Assemblage par RESOU
  13. C               - Conversion au format Morse de la matrice
  14. C               - Conversion du second membre en MVECTD
  15. C               - Construction du préconditionneur
  16. C               - Appel des solveurs itératifs
  17. C               - Conversion du résultat en CHPOINT
  18. C
  19. C
  20. C LANGAGE     : ESOPE
  21. C AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  22. C               mél : gounand@semt2.smts.cea.fr
  23. C***********************************************************************
  24. C VERSION    : v1, 04/08/2011, version initiale
  25. C HISTORIQUE : v1, 04/08/2011, création
  26. C HISTORIQUE :
  27. C HISTORIQUE :
  28. C***********************************************************************
  29.  
  30. -INC PPARAM
  31. -INC CCOPTIO
  32. C On inclue SMCOORD car MCHSOL doit avoir la configuration courante
  33. -INC SMCOORD
  34. -INC SMCHPOI
  35.       POINTEUR MCHSOL.MCHPOI
  36. -INC SMRIGID
  37. -INC SMVECTD
  38.       POINTEUR ISMBR.MVECTD
  39.       POINTEUR INCX.MVECTD
  40.       POINTEUR IR.MVECTD
  41. -INC SMMATRI
  42.       SEGMENT PMORS
  43.       INTEGER IA (NTT+1)
  44.       INTEGER JA (NJA)
  45.       ENDSEGMENT
  46.       POINTEUR PMS1.PMORS,PMS2.PMORS
  47.       POINTEUR KMORS.PMORS
  48. C Segment de stokage
  49.       SEGMENT IZA
  50.       REAL*8  A(NBVA)
  51.       ENDSEGMENT
  52.       POINTEUR IZA1.IZA,IZA2.IZA,IZAU.IZA,IZAL.IZA,ISA.IZA
  53.       POINTEUR KIZA.IZA
  54.  
  55. -INC SMLENTI
  56.       POINTEUR KTYP.MLENTI
  57.       POINTEUR KNOD.MLENTI
  58.       POINTEUR KRINC.MLENTI
  59. -INC SMLMOTS
  60.       POINTEUR IORINC.MLMOTS
  61.       POINTEUR IORINU.MLMOTS
  62. -INC SMTABLE
  63.       POINTEUR KTIME.MTABLE
  64.       DIMENSION ITTIME(4)
  65.       CHARACTER*16 CHARI
  66.       CHARACTER*1 CCOMP
  67.       LOGICAL LTIME,LOGII
  68.       REAL*8 GNRM2
  69. C     ..
  70. C     .. External subroutines and functions..
  71. *inutile      EXTERNAL GAXPY,GCOPY,GDOT,GNRM2
  72.  
  73.       IVALI=0
  74.       XVALI=0.D0
  75.       LOGII=.FALSE.
  76.       IRETI=0
  77.       XVALR=0.D0
  78. *inutile      IOBRE=0
  79.       IRETR=0
  80. C
  81. C Executable statements
  82. C
  83.       IF (LTIME) THEN
  84.          CALL CRTABL(KTIME)
  85.          call timespv(ittime,oothrd)
  86.          ITI1=(ITTIME(1)+ITTIME(2))/10
  87.       ELSE
  88.          KTIME=0
  89.       ENDIF
  90. C
  91. C  CAS PARTICULIER : Si la matrice est vide (toutes les inconnues
  92. C      éliminées, par exemple)
  93. C
  94.       SEGACT MRIGID
  95.       IF (IRIGEL(/2).EQ.0) THEN
  96.          NSOUPO=0
  97.          NAT=0
  98.          SEGINI MCHSOL
  99.          SEGDES MCHSOL
  100.          ICVG=0
  101.          LNMV=0
  102.          LRES=0
  103.          IF (LTIME) THEN
  104.             call timespv(ittime,oothrd)
  105.             ITI2=(ITTIME(1)+ITTIME(2))/10
  106.             CHARI='MATVIDE'
  107.             CALL ECCTAB(KTIME,'MOT     ',IVALI,XVALI,CHARI,LOGII,IRETI,
  108.      $           'ENTIER  ',ITI2-ITI1,XVALR,CHARR,LOGIR,IRETR)
  109.             SEGDES KTIME
  110.          ENDIF
  111.          SEGDES MRIGID
  112.          RETURN
  113.       ENDIF
  114. C
  115. C               - Assemblage par RESOU
  116. C
  117. C  old INORMU=1 : Normalisation des mutiplicateurs de Lagrange
  118. *     INORMU est transmis à la subroutine
  119. *     Le problème est que si MRIGID est deja assemblée, INORMU n'est pas
  120. *     pris en compte... mais où le stocker ??
  121.       CALL KRES9(MRIGID,INORMU)
  122.       IF (IERR.NE.0) RETURN
  123.       IF (LTIME) THEN
  124.          call timespv(ittime,oothrd)
  125.          ITI2=(ITTIME(1)+ITTIME(2))/10
  126.       ENDIF
  127. C
  128. C               - Conversion au format Morse de la matrice
  129. C
  130.       CALL KRES10(MRIGID,KMORS,KIZA)
  131.       IF (IERR.NE.0) RETURN
  132.       IF (LTIME) THEN
  133.          call timespv(ittime,oothrd)
  134.          ITI3=(ITTIME(1)+ITTIME(2))/10
  135.       ENDIF
  136. C
  137. C On donne des infos sur la matrice
  138. C
  139. *      SEGACT MRIGID
  140. *     ICHOLX=ICHOLE
  141. **    INFDDL.ESO est dans ~/triou/p1nc
  142. **      CALL INFDDL(ICHOLX)
  143. C      WRITE(IOIMP,*) 'IMPR=',IMPR
  144.       CALL INFMAT(KMORS,KIZA,IMPR,IRET)
  145. C      IF (IRET.NE.0) GOTO 9999
  146. C      WRITE(IOIMP,*) 'Apres KRES10'
  147. C      WRITE(IOIMP,*) 'KMORS=',KMORS
  148. C      WRITE(IOIMP,*) 'KIZA=',KIZA
  149.  
  150. C
  151. C               - Conversion du second membre en MVECTD
  152. C                 et initialisation du résultat
  153. C
  154.       SEGACT MRIGID
  155.       ICHOLX=ICHOLE
  156.       ISECO=KSMBR
  157. C On ne vérifie pas que le second membre doit être dans le dual
  158.       NOID=1
  159.       CALL CHVNS(ICHOLX,ISECO,ISMBR,NOID)
  160.       IF (IERR.NE.0) RETURN
  161.       IF (LTIME) THEN
  162.          call timespv(ittime,oothrd)
  163.          ITI4=(ITTIME(1)+ITTIME(2))/10
  164.       ENDIF
  165.  
  166. C      SEGACT ISMBR
  167. C      WRITE(IOIMP,*) 'Second membre sous forme vecteur'
  168. C      INC=ISMBR.VECTBB(/1)
  169. C      WRITE(IOIMP,*) '     ISMBR, INC=',INC
  170. C      WRITE(IOIMP,2022) (ISMBR.VECTBB(II),II=1,ISMBR.VECTBB(/1))
  171. C
  172. C     Gestion normalisation Lagrange (repris de MONDES)
  173. C
  174. *      IF (INORMU.EQ.1) THEN
  175.          SEGACT ISMBR*MOD
  176.          MMATRI=ICHOLE
  177.          SEGACT MMATRI
  178.          IF(IDNORD.GT.0) THEN
  179.             MDNO1=IDNORD
  180.          ELSE
  181.             MDNO1=IDNORM
  182.          ENDIF
  183.          SEGACT MDNO1
  184.          INC=MDNO1.DNOR(/1)
  185.          DO 45 I=1,INC
  186.             ISMBR.VECTBB(I)=ISMBR.VECTBB(I)*MDNO1.DNOR(I)
  187.  45      CONTINUE
  188.          SEGDES MDNO1
  189.          SEGDES MMATRI
  190.          SEGDES ISMBR
  191. *      ENDIF
  192. C
  193. C               - Construction du préconditionneur (repris sur kres5)
  194. C               - Appel des solveurs itératifs
  195. C
  196. C Si solveur multigrille, il faut un segment permettant de distinguer
  197. C les inconnues (cf. kres5.eso -> inctyp.eso)
  198.       IF (KTYPI.EQ.7.OR.KTYPI.EQ.8.OR.KTYPI.EQ.10.OR.KTYPI.EQ.11) THEN
  199.          MMATRI=ICHOLE
  200.          SEGACT MMATRI
  201.          MINCPO=IINCPO
  202.          SEGACT MINCPO
  203.          NCOMP=INCPO(/1)
  204.          NNOE=INCPO(/2)
  205.          SEGACT ISMBR
  206.          INC=ISMBR.VECTBB(/1)
  207.          SEGDES ISMBR
  208.          JG=INC
  209.          SEGINI KTYP
  210.          SEGINI KNOD
  211.          MIMIK=IIMIK
  212. *         write(ioimp,*) 'coucou mimik'
  213.          SEGACT,MIMIK
  214. *         SEGPRT,MIMIK
  215.          JG=NCOMP
  216.          SEGINI KRINC
  217. *         WRITE(IOIMP,*) 'NCOMP,IORINC= ',NCOMP,IORINC
  218. *
  219.          IF (IORINC.NE.0) THEN
  220.             SEGACT IORINC
  221.             SEGINI,IORINU=IORINC
  222.             JGN=IORINC.MOTS(/1)
  223.             JGM=IORINC.MOTS(/2)
  224. *            write(ioimp,*) 'JGN,JGM=',JGN,JGM
  225.             CALL CUNIQ(IORINC.MOTS,JGN,JGM,IORINU.MOTS,JGMU,IMPR,IRET)
  226.             IF (IRET.NE.0) GOTO 9999
  227.             IF (JGM.NE.JGMU) THEN
  228.                WRITE(IOIMP,*) 'IORINC ne doit pas avoir de doublons'
  229.                GOTO 9999
  230.             ENDIF
  231.             SEGSUP IORINU
  232.             IF (JGM.NE.NCOMP) THEN
  233.                WRITE(IOIMP,*)
  234.      $      'IORINC doit referencer toutes les inconnues de la matrice'
  235.                GOTO 9999
  236.             ENDIF
  237.             CALL CREPER(JGN,NCOMP,JGM,MIMIK.IMIK,IORINC.MOTS,KRINC.LECT
  238.      $           ,IMPR,IRET)
  239.             IF (IRET.NE.0) GOTO 9999
  240.          ELSE
  241.             DO ICOMP=1,NCOMP
  242.                KRINC.LECT(ICOMP)=ICOMP
  243.             ENDDO
  244.          ENDIF
  245.          IF (IMPR.GT.2) THEN
  246.             WRITE(IOIMP,*) 'NCOMP= ',NCOMP
  247.             WRITE(IOIMP,*) 'MIMIK.IMIK(1..',NCOMP,')= '
  248.             WRITE(IOIMP,*)(MIMIK.IMIK(II),II=1,NCOMP)
  249.             IF (IORINC.NE.0) THEN
  250.                WRITE(IOIMP,*) 'IORINC.MOTS(1..',NCOMP,')= '
  251.                WRITE(IOIMP,*)(IORINC.MOTS(II),II=1,NCOMP)
  252.             ENDIF
  253.             WRITE(IOIMP,*) 'KRINC.LECT(1..',NCOMP,')= '
  254.             WRITE(IOIMP,*)(KRINC.LECT(II),II=1,NCOMP)
  255.          ENDIF
  256. *
  257. *         write(ioimp,*) 'KPREC=',KPREC
  258.          DO ICOMP=1,NCOMP
  259.             DO INOE=1,NNOE
  260.                IG=INCPO(ICOMP,INOE)
  261.                IF (IG.GT.0) THEN
  262.                   KTYP.LECT(IG)=KRINC.LECT(ICOMP)
  263.                   KNOD.LECT(IG)=INOE
  264.                ENDIF
  265.             ENDDO
  266.          ENDDO
  267.          SEGSUP KRINC
  268.          SEGDES KTYP
  269.          SEGDES KNOD
  270.          SEGDES MINCPO
  271.          SEGDES MMATRI
  272.       ELSE
  273.          KTYP=0
  274.       ENDIF
  275. C
  276. C Warning KMORS, KIZA et KTYP sont détruits dans KRES11 et KRES12
  277. C si inodet=0
  278.       INODET=1
  279. C         CALL ECMORS(KMORS,KIZA,4)
  280. C      SEGACT ISMBR
  281. C      WRITE(IOIMP,*) 'Second membre sous forme vecteur'
  282. C      INC=ISMBR.VECTBB(/1)
  283. C      WRITE(IOIMP,*) '     ISMBR, INC=',INC
  284. C      WRITE(IOIMP,2022) (ISMBR.VECTBB(II),II=1,ISMBR.VECTBB(/1))
  285. C Solveur Direct
  286.       IF (KTYPI.EQ.1) THEN
  287.          SEGINI,INCX=ISMBR
  288.          CALL KRES12(KMORS,KIZA,INCX,
  289. C         CALL KRES12(KMORS,KIZA,ISMBR,
  290.      $        KTIME,LTIME,
  291.      $        INCX,LRES,LNMV,ICVG,IMPR,INODET)
  292.       ELSE
  293. C Solveur Itératif
  294.          CALL KRES11(KMORS,KIZA,KTYP,KNOD,ISMBR,
  295.      $        KTYPI,ITER,RESID,ICALRS,IRSTRT,LBCG,BRTOL,IDDOT,IMVEC,
  296.      $        KPREC,RXMILU,RXILUP,XLFIL,XDTOL,XSPIV,
  297.      $        KTIME,LTIME,
  298.      $        INCX,LRES,LNMV,ICVG,IMPR,INODET)
  299. C      WRITE(IOIMP,*) 'Apres KRES11'
  300.       ENDIF
  301.       IF(IERR.NE.0) RETURN
  302. C      SEGACT INCX
  303. C      WRITE(IOIMP,*) 'Inconnue sous forme vecteur'
  304. C      INC=INCX.VECTBB(/1)
  305. C      WRITE(IOIMP,*) '     INCX, INC=',INC
  306. C      WRITE(IOIMP,2022) (INCX.VECTBB(II),II=1,INCX.VECTBB(/1))
  307. C      IF(IERR.NE.0) RETURN
  308. C                    r(0)=b
  309. C      SEGINI,IR=ISMBR
  310. C      SEGACT INCX
  311. C      SEGACT KMORS
  312. C      SEGACT KIZA
  313. CC                    r(0)=b-Ax
  314. C      CALL GMOMV(IMVEC,'N',-1.D0,KMORS,KIZA,INCX,1.D0,IR)
  315. C      RNRM2 = GNRM2(IR)
  316. C      WRITE(IOIMP,*) '||R||=',RNRM2
  317. C
  318. C     Gestion normalisation Lagrange (repris de MONDES)
  319. C     + égalité multiplicateurs
  320. *      IF (INORMU.EQ.1) THEN
  321.          SEGACT INCX*MOD
  322.          MMATRI=ICHOLE
  323.          SEGACT MMATRI
  324.          MDNOR=IDNORM
  325.          SEGACT MDNOR
  326.          INC=DNOR(/1)
  327.          DO 35 I=1,INC
  328.             INCX.VECTBB(I)=INCX.VECTBB(I)*DNOR(I)
  329.  35      CONTINUE
  330.          SEGDES MDNOR
  331.          MILIGN=IILIGN
  332.          SEGACT,MILIGN
  333.          DO 36 I = 1, INC
  334.             if (ITTR(I).ne.0) then
  335. *            write (6,*) ' dans mondes ',i,ittr(i)
  336.                if (incx.vectbb(i).eq.0.d0
  337.      $              .or.incx.vectbb(ittr(i)).eq.0.d0) then
  338. *      write (6,*) ' mondes vectbbs ',vectbb(i+k),vectbb(ittr(i)+k)
  339.                   incx.vectbb(i)=0.d0
  340.                   incx.vectbb(ittr(i))=0.d0
  341.                   goto 36
  342.                endif
  343.                incx.vectbb(i)=(incx.vectbb(i)+incx.vectbb(ittr(i)))/2
  344.                incx.vectbb(ittr(i))=incx.vectbb(i)
  345.             endif
  346.  36      CONTINUE
  347.          SEGDES MILIGN
  348.          SEGDES MMATRI
  349.          SEGDES INCX
  350. *      ENDIF
  351. C
  352. C
  353. C      SEGACT INCX
  354. C      WRITE(IOIMP,*) 'Inconnue sous forme vecteur'
  355. C      INC=INCX.VECTBB(/1)
  356. C      WRITE(IOIMP,*) '     INCX, INC=',INC
  357. C      WRITE(IOIMP,2022) (INCX.VECTBB(II),II=1,INCX.VECTBB(/1))
  358. C      IF(IERR.NE.0) RETURN
  359.  
  360.       IF (LTIME) THEN
  361.          call timespv(ittime,oothrd)
  362.          ITI5=(ITTIME(1)+ITTIME(2))/10
  363.       ENDIF
  364. C
  365. C               - Conversion du résultat en CHPOINT
  366. C
  367.       CALL VCH1(ICHOLX,INCX,MCHSOL,MRIGID)
  368. C      WRITE(IOIMP,*) 'Apres VCH1'
  369.       IF(IERR.NE.0) RETURN
  370.       IF (LTIME) THEN
  371.          call timespv(ittime,oothrd)
  372.          ITI6=(ITTIME(1)+ITTIME(2))/10
  373.          CHARI='RESO ASS+RENU'
  374.          CALL ECCTAB(KTIME,'MOT     ',IVALI,XVALI,CHARI,LOGII,IRETI,
  375.      $                     'ENTIER  ',ITI2-ITI1,XVALR,CHARR,LOGIR,IRETR)
  376.          CHARI='CONVMORS'
  377.          CALL ECCTAB(KTIME,'MOT     ',IVALI,XVALI,CHARI,LOGII,IRETI,
  378.      $                     'ENTIER  ',ITI3-ITI2,XVALR,CHARR,LOGIR,IRETR)
  379. C         CHARI='CONVSMB '
  380. C         CALL ECCTAB(KTIME,'MOT     ',IVALI,XVALI,CHARI,LOGII,IRETI,
  381. C     $                     'ENTIER  ',ITI4-ITI3,XVALR,CHARR,LOGIR,IRETR)
  382.          IF (KTYPI.EQ.1) THEN
  383.             CHARI='KRES FAC+RESO'
  384.          ELSE
  385.             CHARI='KRES PRE+RESO'
  386.          ENDIF
  387.          CALL ECCTAB(KTIME,'MOT     ',IVALI,XVALI,CHARI,LOGII,IRETI,
  388.      $                     'ENTIER  ',ITI5-ITI4,XVALR,CHARR,LOGIR,IRETR)
  389. C         CHARI='CONVINC'
  390. C         CALL ECCTAB(KTIME,'MOT     ',IVALI,XVALI,CHARI,LOGII,IRETI,
  391. C     $                     'ENTIER  ',ITI6-ITI5,XVALR,CHARR,LOGIR,IRETR)
  392.          CHARI='TOTAL   '
  393.          CALL ECCTAB(KTIME,'MOT     ',IVALI,XVALI,CHARI,LOGII,IRETI,
  394.      $                     'ENTIER  ',ITI6-ITI1,XVALR,CHARR,LOGIR,IRETR)
  395.          SEGDES KTIME
  396.       ENDIF
  397. C     Le solveur direct surcharge le second membre
  398.       IF (ISMBR.NE.INCX) SEGSUP ISMBR
  399.       SEGSUP INCX
  400.       SEGDES MRIGID
  401. C
  402. C Normal termination
  403. C
  404.       RETURN
  405. C
  406. C Error Handling
  407. C
  408.  9999 CONTINUE
  409.       MOTERR(1:8)='KRES8   '
  410.       CALL ERREUR(1127)
  411.       RETURN
  412. C
  413. C Format handling
  414. C
  415.  2022 FORMAT(10(1X,1PG12.5))
  416. C
  417. C End of subroutine KRES8
  418. C
  419.       END
  420.  
  421.  

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