Télécharger la plaquette Cast3M


fr en

Télécharger chole.eso

Retour à la liste

Numérotation des lignes : 

chole

  1. C CHOLE     SOURCE    PV090527  26/04/30    21:15:17     12529          
  2. *
  3.       subroutine chole(mmatrx,prec,istab,nbnnma,nligra,xmatri)
  4. *
  5. * octobre 2022 version unifie chole chomod: normal et condensation superelement
  6. *
  7.       IMPLICIT INTEGER(I-N)
  8.       IMPLICIT REAL*8(A-H,O-Z)
  9. C  TANT QUE OOOVAL(1,4) NE MARCHE PAS SUR CRAY
  10.       PARAMETER (LPCRAY=10000)
  11.       INTEGER OOOVAL,OOOLEN
  12.       dimension ittime(4)
  13.  
  14.       POINTEUR LILIGN.MILIGN
  15. C
  16. C  ****  MISE SOUS FORME A=Lt D L  DE LA MATRICE MMATRX
  17. C
  18.  
  19. -INC PPARAM
  20. -INC CCOPTIO
  21. -INC CCREEL
  22. -INC SMMATRI
  23. -INC SMRIGID
  24.  
  25. -INC CCASSIS
  26. -INC CCHOLE
  27.       SEGMENT KIVPO(IIMAX)
  28.       SEGMENT KIVLO(IIMAX)
  29.       segment immt(nblig)
  30.       segment ireser(nvstrm)
  31.       external chole3i
  32.       SAVE IPASV
  33.       DATA IPASV/0/
  34. *  ngmpet dit si on tient en memoire (false) ou si on deborde (true)
  35.       logical ngmpet
  36. C     character*8 zen
  37. C     equivalence (zen,izen)
  38.       logical lsgdes,pasfait,ngdyn
  39.       ireser=0
  40.       nbres=0
  41.       matric=xmatri
  42.       maitre=0
  43.       nelrig=1
  44.       nligrp=nligra
  45.       nligrd=nligra
  46.       if (xmatri.ne.0) then
  47.        rigrel=0
  48.        segini xmatri
  49.       endif
  50.       nbnnmc=nbnnma
  51.       matric=xmatri
  52.  
  53.  
  54.  
  55. 10000 continue
  56.       pasfait=.true.
  57.       lsgdes=.false.
  58. *  faire attention a respecter l'ordre des segdes par la suite
  59.       call ooomru(1)
  60.       if (istab.ne.0.and.istab.ne.1) call erreur(5)
  61.       condmax=0.d0
  62.       condmin=xgrand*xzprec
  63.       ngmpet=.false.
  64.       ngdyn=.true.
  65.       call timespv(ittime,oothrd)
  66.       kcour=(ittime(1)+ittime(2))/10
  67.       kcouri=kcour
  68.       kcour=0
  69.       nbopit=0
  70.  
  71.       iposm=0
  72. C     zen='CPU'//char(0)
  73. C     le=4
  74.       nvaor=0
  75.       nbthro=nbthrs
  76.       ithrd=0
  77.       if (nbthro.gt.1) then
  78.        ithrd=1
  79.        call threadii
  80.        call oooprl(1)
  81.       endif
  82.       nbthr=nbthro
  83.       do ith=1,nbthr
  84.        nbop(ith)=0
  85.       enddo
  86.       stmult=1d-5
  87.  
  88. C  nouvelle methode de gestion de l'espace memoire necessitee par la parallelisation
  89. C  memoire vive totale
  90.       MACTIT=OOOVAL(1,1)
  91. **    write(6,*) ' mactit igrand ',mactit,igrand
  92. C  un bloc de memoire fera au plus macti/2
  93.       call intpdo(inpdo)
  94.       nvstrm=mactit/10
  95.       MMATRI=MMATRX
  96.       SEGACT,MMATRI*MOD
  97.       PRCHLV=PREC
  98.       MILIGN=IILIGN
  99.  
  100.       SEGACT,MILIGN*MOD
  101.       INO=ILIGN(/1)
  102.       MDIAG=IDIAG
  103.       SEGACT,MDIAG*MOD
  104.       NBLIG=INO
  105.       segini immt
  106.       precc=prec
  107.       INC=DIAG(/1)
  108.       if (matric.eq.0) nbnnmc=inc+1
  109.       nvstrm=max(inc*inpdo,nvstrm)
  110. **      write(6,*) ' nvstrm ',nvstrm
  111.       INCC=INC
  112.       MIMIK=IIMIK
  113.       MINCPO=IINCPO
  114.       SEGACT,MINCPO,MIMIK
  115.       IPLUMI=IMIK(/2)*2 +4
  116.       IL2=0
  117.       IIMAX=IJMAX+IPLUMI
  118.       SEGINI KIVPO,KIVLO
  119.       INEG=0
  120.       NBLAG=0
  121.       NENSLX=0
  122.       NVSTOC=0
  123.       NVSTOR=0
  124.       diagmax=XPETIT/XZPREC
  125.       diagmin=xgrand
  126.       do i=1,diag(/1)
  127.        if (ittr(i).eq.0) diagmax=max(diagmax,abs(diag(i)))
  128.        if (ittr(i).eq.0.and.abs(diag(i)).gt.xpetit/xzprec)
  129.      >   diagmin=min(diagmin,abs(diag(i)))
  130.       enddo
  131.       if (diagmax.le.xpetit/xzprec) then
  132.          do i=1,diag(/1)
  133.            diagmax=max(diagmax,abs(diag(i)))
  134.            if (abs(diag(i)).gt.xpetit/xzprec)
  135.      >       diagmin=min(diagmin,abs(diag(i)))
  136.          enddo
  137.       endif
  138.       diagmin=min(diagmin,diagmax)
  139. ***    write (6,*) ' chole diagmin diagmax ',diagmin,diagmax,diag(/1)
  140.  
  141.  
  142. C
  143. C
  144. C
  145. C  **** DEBUT DE LA TRIANGULARISATION. ON PREND NOEUD A NOEUD,
  146. C  **** DECOMPACTAGE PUIS TRAVAIL SUR LES LIGNES DU NOEUDS
  147. C
  148. C  ****  LA LONGUEUR DE LA PLUS GRANDE LIGNE EST DONNEE PAR IMAX
  149. C
  150. C SP  indicateurs pour impression message "stabilisation RESO..."
  151.       isr=0
  152.       isrl=0
  153.  
  154.   1   CONTINUE
  155.       IVALMA=IJMAX+IPLUMI
  156.       IL1=IL2+1
  157.       IMIN=IL1
  158. *  reserver de la place ou mettre les lignes superieures dans le cas debordement
  159.       if (ngmpet) then
  160.         if (ireser.eq.0) segini ireser
  161.       endif
  162.       DO 2 I=IL1,INO
  163.       ngdyn=.true.
  164. **    imasq=0
  165.       LIGN=0  
  166.  
  167.       LLIGN=    ILIGN(I)
  168.       SEGACT /ERR=32/LLIGN
  169.       goto 31
  170.   32  continue
  171. **       write(6,*) ' segact llign erreur',i,il1,lsgdes
  172.       if (.not.lsgdes) then
  173. **         write(6,*) ' lsgdes 1 '
  174.         lsgdes=.true.
  175. *****   ngmpet=.true.
  176. **      write(6,*) 'desactivation-1 ',1,il1-1
  177.         do it=il1-1,1,-1
  178.          lign=ilign(it)
  179.          segdes lign
  180.         enddo
  181.       else
  182.         goto 3
  183.       endif
  184.       SEGACT /ERR=3/LLIGN
  185.   31  continue
  186.       NA=  IMMMM(/1)
  187. *     write (6,*) ' chole ligne noeud inconnues ',i,ipno(i),na
  188.       NBPAR=NA+1
  189.       NVALL=  NJMAX
  190.       LMASQ=masqa(nvall)+1
  191.       SEGINI /ERR=33/ LIGN
  192.       goto 34
  193.   33  continue
  194. **     write(6,*) ' segini lign erreur',il1
  195.       if (.not.lsgdes) then
  196. **        write(6,*) ' lsgdes 2 '
  197.         lsgdes=.true.
  198. *****   ngmpet=.true.
  199. **      write(6,*) 'desactivation-2 ',1,il1-1
  200.         do it=il1-1,1,-1
  201.          lign=ilign(it)
  202.          segdes lign
  203.         enddo
  204.       else
  205.         goto 3
  206.       endif
  207.       SEGINI /ERR=3/ LIGN
  208. **     write(6,*) 'deuxieme essai lign ok'
  209.   34  continue
  210. C  recuperer la longueur du segment
  211.       lglig=na*(nvall/na)**1.333333333333333333
  212.       NVSTOC=NVSTOC + NVALL
  213.       IVALMA=IVALMA + NVALL
  214.  
  215.  
  216.       nvaor = nvaor + xxva(/1)
  217. C
  218. C  ****  DECOMPACTAGE
  219. C
  220.       IPA=1
  221.       DO 121 JPA=1,NA
  222.             IVPO(JPA)=IPA
  223.         KPA =   IPPO(JPA+1)-  IPPO(JPA)
  224.         IPP =   IPPO(JPA)
  225.             IPPVV(JPA)=IPA-1
  226.         LPA  =   LDEB(JPA)
  227.         LPA1 = LPA-IPA
  228.  
  229.         DO 122 MPA=1,KPA
  230.           LL   =   LINC(MPA+IPP)
  231.           IPLA = LL -LPA1
  232.            xxv=xxva(mpa+ipp)
  233.            if (abs(xxv).gt.xpetit) then
  234.               VAL(IPLA)=  XXV
  235.               IMASQ(masqa(IPLA))=1
  236.               if (ipla-ipa+1.ge.1) IMASQ(masqa(IPLA-ipa+1))=1
  237.            endif
  238. 122     CONTINUE
  239.  
  240.         IPA=IPA+  IMMMM(JPA)-LPA + 1
  241. Cpv         IMMM(JPA)=   IPNO(LPA)
  242.         IMMM(JPA)=LPA
  243.         IF(IMIN  .GT.   IPNO(LPA )) IMIN  =    IPNO(LPA)
  244.  
  245.  
  246.  
  247. 121   CONTINUE
  248.  
  249. * indexation de imasq
  250.         ipln=lmasq/na
  251.         iplp=lmasq/na
  252.         do 123 ipl=lmasq/na,1,-1
  253.           if (imasq(ipl).gt.0) then
  254.             imasq(ipl)=masqi(iplp+1)
  255.             ipln=ipl-1
  256.           else
  257.             imasq(ipl)=-masqi(ipln+1)
  258.             iplp=ipl-1
  259.           endif
  260.  123    continue
  261. **      write (6,*) ' imasq ',lmasq/na
  262. **      write (6,*) (imasq(ipl),ipl=1,lmasq/na)
  263.  
  264.  
  265.  
  266.  
  267.  
  268.  
  269.  
  270.  
  271. C*** **** ****
  272.  
  273.  
  274.  
  275.  
  276.  
  277.  
  278.  
  279.  
  280.  
  281.  
  282.  
  283.       if (na.gt.0) then
  284.           IPREL=  IMMMM(1)
  285.           IDERL=  IMMMM(NA)
  286.           lcara(2,i)=iprel
  287.           lcara(3,i)=iderl
  288.       endif
  289.           IPPVV(NA+1)=IPA-1
  290.  
  291.  
  292.  
  293.  
  294.  
  295.       SEGSUP LLIGN
  296.          ILIGN(I)=LIGN
  297.  
  298.  
  299. C*    write (6,*) 'longueur ligne ',nvall
  300. C nb de ligne multiple du nb de threads
  301. C blocage ligne lecture-ecriture pour minimiser le cache
  302. C on note si on est au minimum de lignes
  303.       nbthro=min(nbthrs,lglig/1200+1)
  304.       if (i+1-il1.ge.nbthro.and.(.not.ngmpet)) then
  305.        nbthro=min(nbthrs,i+1-il1)
  306.        ngdyn=.true.
  307.         if(i+1-il1.eq.nbthrs) ngdyn=.false.
  308.          il2=i
  309. **       write(6,*) ' nbthro il1 il2 ',nbthro,il1,il2
  310.  
  311.  
  312.          GOTO 4
  313.       endif
  314.  
  315.  
  316.     2 CONTINUE
  317.       IL2=INO
  318.       GO TO 4
  319.     3 IL2=I-1
  320. **      write(6,*) 'desactivation-4 ',llign
  321.       segdes llign
  322.       if (               lign.ne.0) segsup lign
  323.     4 CONTINUE
  324.       nbthro=min(nbthrs,nbthro)
  325.       nbthr=nbthro
  326.       if(ireser.ne.0) segsup ireser
  327. C
  328.       IF(IL2.GE.IL1) GO TO 40
  329. C
  330. C ****  APPEL AUX ERREURS  MESSAGE PAS ASSEZ DE PLACE MEMOIRE
  331. C
  332. C      ITYP=48
  333.       CALL ERREUR(48)
  334.         call ooodmp(0)
  335.       if (ithrd.eq.1) then
  336.        call threadis
  337.        call oooprl(0)
  338.       endif
  339.       call ooomru(0)
  340.       RETURN
  341.    40 CONTINUE
  342.       IM=INC
  343.  
  344.  
  345.  
  346.       DO 352 IH=IL2 ,IL1,-1
  347.       LIGN=   ILIGN(IH  )
  348.       IL=INC
  349.       DO 354 JH=1,    IMMM(/1)
  350.       IM=MIN(IM,    IMMM(JH))
  351.       IL=MIN(IL,    IMMM(JH))
  352.  354  CONTINUE
  353.       IML=IL
  354.       lcara(1,ih)=iml
  355.       immt(ih)=ipno(im)
  356.   352 CONTINUE
  357. C  353 CONTINUE
  358.       LIGN=ILIGN(IL1)
  359.  
  360.  
  361.  
  362.  
  363.  
  364.  
  365.  
  366.  
  367.  
  368.  
  369.  
  370.  
  371.  
  372.  
  373.  
  374.  
  375.       IL11=    IPREL
  376.  
  377.  
  378. C
  379. C  ****  BOUCLE  *5*  TRAVAILLE SUR LE NOEUD I QUI EST EN LECTURE
  380. C
  381. C    >      il1,il2
  382.       lig1=   ilign(imin)
  383.       ipos=0
  384.       iper=imin
  385.       ider=imin-1
  386.       iderac=imin-1
  387.  
  388.       IL1i=il1
  389.       DO 5 I=IMIN ,IL2
  390.  
  391.       LIG1=   ILIGN(I)
  392.  
  393.       IF(I.LT.IL1)  GO TO 7
  394. C____________
  395. C
  396. C  ******* LE NOEUD I EST EN MEMOIRE IL EST TRIANGULE JUSQU'A
  397. C  ******* IPREL  IL FAUT CONTINUER TOUTE LES LIGNES PUIS CALCULER
  398. C  ******* LE TERME DIAGONAL
  399. C
  400.       LIGN=LIG1
  401.       DO 156 KHG=1,IMMM(/1)
  402.  
  403.       II=IPREL-1+KHG
  404.       IMMM(KHG)=0
  405.       NN=IPPVV(KHG+1)
  406.  
  407.       NNM1=IPPVV(KHG)
  408.  
  409.       N=NN-NNM1
  410.       DIAG(II)=VAL(NN)
  411.        if (n.eq.1) then
  412.         if(ii-nbnnmc.gt.0) then
  413.          re(ii-nbnnmc,ii-nbnnmc,1)=val(nn)
  414.          goto 41
  415.         else
  416.          goto 8
  417.         endif
  418.        endif
  419.  
  420.  
  421.       NMI=N-II
  422. *  attention subtile difference noeud maitre et non maitre
  423.       if (ii-nbnnmc.le.0) then
  424.        IDEP=MAX(IL11,2-NMI)
  425.       else
  426.        IDEP=MAX(IL11,1-NMI)
  427.       endif
  428.       KIDEP=IDEP+NMI
  429.       KI1=N-1
  430.       KQ=-NMI
  431.       VAL(NN)=VAL(NN)+
  432.      > CHOLE1(ILIGN,IPREL,imasq(1),VAL(1+IPPVV(KHG)),DIAG(1-NMI),
  433.      > IPNO(1-NMI),IPPVV(1),KHG,IVPO(1),KIDEP,KI1,KQ,1+IPPVV(KHG),
  434.      > PREC,nbop(1))
  435.        imasq(masqa(nn))=1
  436.        imasq(masqa(n))=1
  437.        if(ii-nbnnmc.gt.0) then
  438.          re(ii-nbnnmc,ii-nbnnmc,1)=val(nn)
  439.          goto 41
  440.        endif
  441.    8  CONTINUE
  442.       diagref=max(abs(diag(ii)),diagmin)
  443.       diagcmp=diagref*5d-12
  444.       IF(   ITTR(II).EQ.0.AND.
  445.      & ABS(     VAL(NN)).GT.diagcmp)  GO TO 12
  446.       IF(   ITTR(II).NE.0.AND.
  447.      & ABS(     VAL(NN)).GT.diagcmp)  GO TO 12
  448. **      write(6,*) ' ittr val diagcmp ',ittr(ii),val(nn),diagcmp
  449. C  il faut mettre une valeur plus grande sur les LX car on a un probleme de conditionnement
  450. C  sur le calcul des reactions en cas de 2 relations presque identique
  451. C
  452. C  **** ON VIENT DE DETECTER UN MODE D'ENSEMBLE
  453. C  **** ON AJOUTE A LA STRUCTURE UN RESSORT EGAL A CELUI QUI EXISTAIT
  454. C  **** AU PREALABLE SUR CETTE INCONNUE.
  455. C
  456. *     write (6,*) ' chole mode d ensemble ittr ligne ',
  457. *    >      ittr(ii),ii,diag(ii),val(nn),diagref
  458. C  on garde le signe car il fau un moins sur les ML
  459.        vmaxi=diagref
  460.        do ipv=1+ippvv(khg),nn
  461.         vmaxi=max(vmaxi,abs(val(ipv)))
  462.        enddo
  463.       if(   ittr(ii).NE.0) then
  464.        VAL(NN)=val(nn)-1.30D0*diagref
  465.        NENSLX=NENSLX+1
  466.       else
  467.        val(nn)=vmaxi
  468.       endif
  469.       NENS=NENS+1
  470.            IMMM(KHG)=NENS
  471.    12 CONTINUE
  472.  
  473.  
  474. *  stabilisation
  475.       IF (ISTAB.NE.0) THEN
  476. *  elimination des Nan
  477.        if (.not.(abs(val(nn)).lt.xgrand*xzprec).and.pasfait) then
  478.          pasfait=.false.
  479.          val(nn)=xgrand*xzprec
  480.          write (6,*) 'Nan dans chole ligne',ii
  481.        endif
  482. *
  483.        diagcmp=abs(diagmin)*1d-5+xpetit/xzprec
  484.        if (val(nn).lt.-diagmax*1d-3.and.ittr(ii).eq.0) then
  485.         val(nn)=abs(val(nn))
  486. ***     val(nn)=diagmax*1d-3
  487.         if (immm(khg).eq.0) NENS=NENS+1
  488.         IMMM(KHG)=NENS
  489.        elseif (val(nn).le.diagcmp.and.ittr(ii).eq.0 ) then
  490.         if (isr.eq.0.or.iimpi.gt.0)
  491.      &     write (ioimp,*) ' stabilisation RESO ',ii,val(nn),diag(ii)
  492.         isr=isr+1
  493.         val(nn)=max(diagcmp,-val(nn))
  494.         if (immm(khg).eq.0) NENS=NENS+1
  495.         IMMM(KHG)=NENS
  496.        else
  497.         if (ittr(ii).eq.0) diagmin=min(diagmin,abs(val(nn)))
  498.        endif
  499.        if (val(nn).ge.abs(diag(ii))*stmult.and.ittr(ii).ne.0) then
  500.         if (isrl.eq.0.or.iimpi.gt.0)
  501.      &     write (ioimp,*) ' stabilisation RESO lagrange ',ii,val(nn)
  502.         isrl=isrl+1
  503.         val(nn)=-abs(diag(ii))*stmult
  504.         if (immm(khg).eq.0) then
  505.            NENS=NENS+1
  506.            NENSLX=NENSLX+1
  507.         endif
  508.         IMMM(KHG)=NENS
  509.        endif
  510.       ENDIF
  511.  
  512.  
  513.  
  514.  
  515.  
  516.       DIAG(II)=     VAL(NN)
  517.       IF(abs(DIAG(II)).gt.xpetit)  GO TO 41
  518.       DIAG(II)=diagmax
  519.       if (ittr(ii).ne.0) diag(ii)=-diagmax
  520.       VAL(NN)=DIAG(II)
  521.       GO TO 41
  522. C
  523. C  ****  ENVOI ERREUR MATRICE SINGUIERE
  524. C
  525. C      ITYP=49
  526.       INTERR(1)=I
  527.       CALL ERREUR(49)
  528.       if (ithrd.eq.1) then
  529.         call threadis
  530.        call oooprl(0)
  531.       endif
  532.       call ooomru(0)
  533.       RETURN
  534. C
  535. C  ****  ON COMPTE LE NOMBRE DE TERMES DIAGONAUX NEGATIFS
  536. C        ET LE NOMBRE DE MULTIPLICATEUR DE LAGRANGE
  537. C
  538.    41 IF(DIAG(II).LT.0.D0)  INEG=INEG+1
  539.       IF(   ITTR(II).NE.0) NBLAG=NBLAG+1
  540.  
  541.       if (ii.le.nbnnmc) condmin=min(condmin,abs(diag(ii)))
  542.       condmax=max(condmax,abs(diag(ii)))
  543.       if (ii.le.nbnnmc) diag(ii)=1.d0/diag(ii)
  544.   156 CONTINUE
  545.       NA=IMMM(/1)
  546. C
  547. C  RECOMPACTAGE DE LIGN (DEJA ENTIEREMENT TRAITEE)
  548. C
  549. *       write(6,*) 'chole ligne lpl',i,na,ippvv(2)-ippvv(1)
  550.       if (na.gt.0)
  551.      >CALL COMPAC(     VAL(1),NBPAR,KIVPO(1),KIVLO(1),
  552.      #    NVALL,     IPPVV(1),IZROSF,NA,PREC,imasq(1),iprel,iderl)
  553.       lmasq=0
  554. C on recree lign car la compacter en place emiette la memoire
  555.        lig1=lign
  556.        segini /err=1431/ lig1
  557.        goto 1432
  558.  1431  continue
  559.        write(6,*) ' segini echoue'
  560.  1432  continue
  561.       nbres=max(nbres,nvall)
  562. *  deplacement fait ici maintenant, avec unrolling
  563.       do 300 nbp=1,nbpar-1
  564.        kdif =kivpo(nbp)-kivlo(nbp)
  565.        do     iv=kivlo(nbp),kivlo(nbp+1)-4,4
  566.           lig1.val(iv)=val(iv+kdif )
  567.           lig1.val(iv+1)=val(iv+1+kdif )
  568.           lig1.val(iv+2)=val(iv+2+kdif )
  569.           lig1.val(iv+3)=val(iv+3+kdif )
  570.        enddo
  571.        do     iv1=iv,kivlo(nbp+1)-1
  572.           lig1.val(iv1)=val(iv1+kdif )
  573.        enddo
  574.  300  continue
  575. *      do it=1,nvall
  576. *       lig1.val(it)=val(it)
  577. *      enddo
  578.        do it=1,na
  579.         lig1.immm(it)=immm(it)
  580.         lig1.ippvv(it)=ippvv(it)
  581.        enddo
  582.         lig1.ippvv(na+1)=ippvv(na+1)
  583.         lig1.iml=iml
  584.         lig1.iprel=iprel
  585.         lig1.iderl=iderl
  586.         lig1.iml=iml
  587.         lcara(2,i)=lig1.iprel
  588.         lcara(3,i)=lig1.iderl
  589.         lcara(1,i)=lig1.iml
  590.        if (lign.ne.lig1) then
  591.         segsup lign
  592.         lign=lig1
  593.        else
  594.         segadj lign
  595.        endif
  596.           ilign(i)=lign
  597.       NVSTOR=NVSTOR+NVALL
  598.       nvstrm=max(nvstrm,nvall*inpdo)
  599.       DO 143 LHG=1,NBPAR
  600.             IVPO(2*LHG-1)=KIVPO(LHG)
  601.             IVPO(2*LHG)  =KIVLO(LHG)
  602.  143  CONTINUE
  603. C   Si la ligne est petite, il n'y a rien a gagner a paralleliser
  604.  
  605.       if (i.gt.1) then
  606.        lig1=ilign(i-1)
  607. **      if(lsgdes) write(6,*) 'desactivation-5 ',i
  608.        if (lsgdes)  segdes lig1
  609.        iderac=min(iderac,i-2)
  610.       endif
  611.  
  612. C
  613. C  ****  ON TRIANGULARISE LES AUTRES LIGNES
  614. C
  615.       il1=il1+1
  616.       if (il1.gt.il2) goto 5
  617.       LIG1=ILIGN(I)
  618.       lign=ilign(il1)
  619.       IL11=IPREL
  620.       goto 7
  621. C  72  continue
  622.   71  continue
  623. *  passage en superlent
  624.       if (ider.lt.il1-1.and..not.ngmpet) then
  625. **      write(6,*) 'ngmpet true ',iper,ider,il1,il2
  626.         ngmpet=.true.
  627. **      do ipv=1,iper-1
  628. **        lig2=ilign(ipv)
  629. **       call oooeta(lig2,ieta,imod)
  630. **         if (ieta.eq.1) write(6,*) 'ligne active ',ipv,iper
  631. **      enddo
  632.       endif
  633.       if (iper.gt.ider) then
  634.         call erreur(48)
  635.         write(6,*) ' 1 chole iper ider i il1 il2 ',
  636.      >         iper,ider,i,il1,il2
  637.         do ipv=1,ino
  638.          lig2=ilign(ipv)
  639.          call oooeta(lig2,ieta,imod)
  640.          if (ipv.lt.il1.or.ipv.gt.il2) then
  641.            if (ieta.eq.1) write(6,*) 'ligne active ',ipv
  642.          endif
  643.         enddo
  644.         call ooodmp(0)
  645.       if (ithrd.eq.1) then
  646.        call threadis
  647.        call oooprl(0)
  648.       endif
  649.       call ooomru(0)
  650.         return
  651.       endif
  652.  
  653. C      soit parce qu'on a fini, soit parce qu'on manque de memoire
  654. C      il faut executer les lignes activees puis les desactiver
  655. C      lancer les chole3 et attendre qu'ils soient finis
  656.       if (ipos.ne.0) then
  657. **      write (6,*) ' lancement thread ',iper,ider,il1,il2
  658.         if (iper.gt.ider) then
  659.           write (6,*) ' erreur interne chole '
  660.           call erreur(5)
  661.         endif
  662.         nbthr=min(nbthr,il2-il1+1)
  663. **      write(6,*) 'nbthr ',nbthr,il2-il1+1
  664.          LILIGN=MILIGN
  665. *  blocage pour rester dans le cache secondaire
  666.         ipers=iper
  667.         iders=ider
  668.         ipas=1500
  669.         if(nbthr.eq.1) ipas=igrand
  670. *       do 400 iper=ipers,iders,ipas
  671.  401    continue
  672.         ider=min(iders,iper+ipas-1)
  673.           do ith=1,nbthr-1
  674.             call threadid(ith,chole3i)
  675.           enddo
  676.           call chole3i(nbthr)
  677.           do ith=nbthr-1,1,-1
  678.             call threadif(ith)
  679.           enddo
  680.         iper=iper+ipas
  681.         ipas=ipas/2
  682.         ipas=max(ipas,750)
  683.         if (iper.le.iders) goto 401
  684.  400    continue
  685.         iper=ipers
  686.         ider=iders
  687.       endif
  688.  
  689. C  test ctrlC
  690.       if (ierr.ne.0) goto 9999
  691.       iposm=max(iposm,ipos)
  692.       ipos=0
  693.       iderf=ider-1
  694.       if (ider.ne.il1-1) iderf=ider
  695.       if (lsgdes) then
  696. **      write(6,*) 'desactivation 7 iper iderf',iper,iderf
  697.       do il=iderf,iper,-1
  698.         lign=ilign(il)
  699.         segdes lign
  700. C       write (6,*) ' desactivation ligne ',il
  701.       enddo
  702.       endif
  703.        iderac=min(iderac,iper-1)
  704.       iper=ider+1
  705. C     write (6,*) ' iper ider il1 ',iper,ider,il1
  706.       if (iper.ne.il1) goto 7
  707.  
  708.       goto 5
  709.     7 CONTINUE
  710.  
  711.       if (lsgdes) SEGACT/err=71/LIG1
  712.       ipos=ipos+1
  713.       ider=i
  714.       if (i.gt.iderac) iderac=i
  715.       if (i.eq.il1-1)   goto 71
  716.     5 CONTINUE
  717. **    write (6,*) ' il1 il2 apres 5 ',il1,il2
  718.       if (lsgdes) then
  719. **      write(6,*) 'desactivation 8 il1 il2',il1,il2
  720.       DO I=min(il1,il2),max(il1,il2)
  721.         LIGN=   ILIGN(I)
  722.         if(lign.ne.0) SEGDES,LIGN
  723.       enddo
  724.       endif
  725.       nbopt=0
  726.       do ith=1,nbthro
  727.        nbopt=nbopt+nbop(ith)
  728.        nbop(ith)=0
  729.       enddo
  730.       nbopin=nbopt
  731.       nbopit=nbopit+nbopin
  732.       call timespv(ittime,oothrd)
  733.       kcour=(ittime(1)+ittime(2))/10
  734.  
  735.  
  736.  
  737.  
  738.  
  739.       iderac=min(iderac,il1-1)
  740.       if (ireser.ne.0) segsup ireser
  741.       IF(IL2.LT.INO)  GO TO 1
  742. C  ON MET A JOUR LE NOMBRE DE TERMES DIAGONAUX NEGATIF
  743. C   ON ENLEVE LE NOMBRE DE MULTIPLICATEUR DE LAGRANGE
  744. C       INEG=INEG-NBLAG
  745. C     on ne compte pas 2 fois les multiplicateurs qui vont etre
  746. C     elimines lors de la resolution car mode d'ensemble
  747.       INEG=INEG-(NBLAG-NENSLX)
  748.       if (iimpi.ne.0.and.NENSLX.gt.0) WRITE(IOIMP,4820) NENSLX
  749.  4820 FORMAT(I12,' MODES D ENSEMBLE PORTANT SUR DES MULTIPLICATEURS',
  750.      1' DE LAGRANGE DETECTES')
  751.  
  752.       IF (IIMPI.EQ.1) WRITE(IOIMP,4821) NVSTOC
  753.  4821 FORMAT( ' NOMBRE DE VALEURS  DANS LE PROFIL',I12)
  754.       IF (IIMPI.EQ.1) WRITE(IOIMP,4822) NVSTOR
  755.  4822 FORMAT( ' NOMBRE DE VALEURS STOCKEES DANS LE PROFIL',I12)
  756.       IF (IIMPI.EQ.1) WRITE(IOIMP,4825) Nbopit/1000000
  757.  4825 FORMAT( ' NOMBRE DE GIGA OPERATIONS FMA',I40)
  758.       IF (IIMPI.EQ.1) WRITE(IOIMP,4823) NVaor
  759.  4823 FORMAT( ' NOMBRE DE VALEURS initiales',I9)
  760. C     IF (IIMPI.EQ.1) WRITE(IOIMP,4824) nbopit/1d6/(kcour-kcouri)
  761. C 4824 FORMAT( ' Performance en Gigaflops ',F8.1)
  762.       INTERR(1)=NVSTOR
  763.       reaerr(1)=nvstor/inc**(4.D0/3.D0)
  764.       reaerr(2)=2.D0*nbopit/1.D6/max(1,(kcour-kcouri))
  765.       reaerr(3)=condmax/condmin
  766.       IF (IPASV.EQ.0.or.reaerr(3).gt.1.D30) CALL ERREUR(-278)
  767.       IPASV=1
  768.       call ooomru(0)
  769.       if(lsgdes) then
  770.         do ipv=1,ino
  771.           lign=ilign(ipv)
  772.           segdes lign
  773.         enddo
  774.       endif
  775.       SEGDES,MINCPO
  776.       SEGDES,MIMIK
  777.       SEGDES,MMATRI
  778.       SEGDES,MILIGN
  779.       MMATRX=MMATRI
  780.       SEGSUP KIVPO,KIVLO
  781.       segsup immt
  782. C     write (6,*) ' chole ipos max ',iposm
  783.  9999 continue
  784.       if (ithrd.eq.1) then
  785.         call threadis
  786.        call oooprl(0)
  787.       endif
  788.       if(iimpi.ne.0)  write (6,*) 
  789.      >      ' chole condmin condmax ',condmin,condmax,diag(/1)
  790.       SEGDES,MDIAG
  791. **    write(6,*) 'kseq kpar',xkseq,xkpar
  792.  
  793.  
  794.       RETURN
  795.       END
  796.  
  797.  
  798.  
  799.  
  800.  
  801.  
  802.  
  803.  
  804.  
  805.  

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