C LDMTS     SOURCE    PV090527  26/04/30    21:15:49     12529          
      SUBROUTINE LDMTS(MMATRX,PREC,ISTAB,NBNNMA,NLIGRA,XMATRI)
*
      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8(A-H,O-Z)
C  TANT QUE OOOVAL(1,4) NE MARCHE PAS SUR CRAY
      PARAMETER (LPCRAY=10000)
      INTEGER OOOVAL,OOOLEN
      dimension ittime(4)
C
C  ****  MISE SOUS FORME A=LDM DE LA MATRICE MMATRX
C
-INC PPARAM
-INC CCOPTIO
-INC CCREEL
-INC SMMATRI
-INC SMRIGID
-INC CCASSIS
-INC CCHOLE
-INC SMCOORD
C
      POINTEUR LILIGN.MILIGN
      POINTEUR LLIGL.LLIGN,LLIGM.LLIGN
      POINTEUR LIGL.LIGN,LIGM.LIGN,LIG5.LIGN,LIG3.LIGN
      SEGMENT KIVPO(IIMAX)
      SEGMENT KIVLO(IIMAX)
      segment immt(nblig)
      segment xreser              
       integer ireser(ivstrm)
       real yreser(nvstrm)
      endsegment
      SEGMENT ILR
        integer ILIGR(NNOE)
      ENDSEGMENT
      POINTEUR ILIGNS.ILR
      external chole4i
      external shole3i
      SAVE IPASV
      DATA IPASV/0/
*  ngmpet dit si on tient en memoire (false) ou si on deborde (true)
      logical ngmpet
      logical lsgdes,pasfait
C
      SEGDES,MCOORD
      maitre=0                                                                              
      xreser=0
      pasfait=.true.
      lsgdes=.false.
*  faire attention a respecter l'ordre des segdes par la suite
      call ooomru(1)
      if (istab.ne.0.and.istab.ne.1) call erreur(5)
      condmax=0.d0
      condmin=xgrand*xzprec
      ngmpet=.false.
      nvallp=nvall
      call timespv(ittime,oothrd)
      kcour=(ittime(1)+ittime(2))/10
      kcouri=kcour
      kcour=0
      nbopit=0
      nvaor=0
      nbthro=nbthrs
      ithrd=0
      if (nbthro.gt.1) then
       ithrd=1
       call threadii
       call oooprl(1)
      endif
      nbthr=nbthro
      do ith=1,nbthr
       nbop(ith)=0
      enddo
      stmult=1d-5

      nvstrm=0             
      ivstrm=oooval(1,4)*2
C     write(6,*) 'nvstrm initial',nvstrm
      MMATRI=MMATRX
      SEGACT,MMATRI*MOD
      PRCHLV=PREC
      MILIGN=IILIGN
      LILIGN=IILIGS
      SEGACT,MILIGN*MOD,LILIGN*MOD
C
      INO=MILIGN.ILIGN(/1)
      NNOE=INO
      SEGINI,ILIGNS
C
      MDIAG=IDIAG
      SEGACT,MDIAG*MOD
      NBLIG=INO
      segini immt
      precc=prec
      INC=DIAG(/1)
C
C     Cas du super-element
      nbnnmc=inc+1
      if (xmatri.ne.0) then
        nelrig=1
        nligrp=nligra
        nligrd=nligra
        rigrel=0
        segini xmatri
        nbnnmc=nbnnma
      endif
      matric=xmatri
      INCC=INC
      MIMIK=IIMIK
      MINCPO=IINCPO
      SEGACT,MINCPO,MIMIK
      IPLUMI=IMIK(/2)*2 +4
      IL2=0
      IIMAX=IJMAX+IPLUMI
      SEGINI KIVPO,KIVLO
      INEG=0
      NBLAG=0
      NENSLX=0
      NVSTOC=0
      NVSTOR=0
      diagmax=XPETIT/XZPREC
      diagmin=xgrand
      do i=1,diag(/1)
       if (ittr(i).eq.0) diagmax=max(diagmax,abs(diag(i)))
       if (ittr(i).eq.0.and.abs(diag(i)).gt.xpetit/xzprec)
     >   diagmin=min(diagmin,abs(diag(i)))
      enddo
      if (diagmax.le.xpetit/xzprec) then
         do i=1,diag(/1)
           diagmax=max(diagmax,abs(diag(i)))
           if (abs(diag(i)).gt.xpetit/xzprec)
     >       diagmin=min(diagmin,abs(diag(i)))
         enddo
      endif
      diagmin=min(diagmin,diagmax)
***    write (6,*) ' ldmts diagmin diagmax ',diagmin,diagmax,diag(/1)
C
C  **** DEBUT DE LA TRIANGULARISATION. ON PREND NOEUD A NOEUD,
C  **** DECOMPACTAGE PUIS TRAVAIL SUR LES LIGNES DU NOEUDS
C
C  ****  LA LONGUEUR DE LA PLUS GRANDE LIGNE EST DONNEE PAR IMAX
C
C SP  indicateurs pour impression message "stabilisation RESO..."
      isr=0
      isrl=0

  1   CONTINUE
      IL1=IL2+1
      IMIN=IL1
C
C     Reserver de la place ou mettre les lignes superieures dans le cas debordement
      if (ngmpet) then
        if (xreser.eq.0) then
CC      write(6,*)'segini xreser',nvstrm
        segini xreser
        endif
      endif
C
      DO 2 I=IL1,INO
C
      LIGL=0
      LIGM=0
      LGLIG=0
C
      LLIGM=LILIGN.ILIGN(I)
      LLIGL=MILIGN.ILIGN(I)
      SEGACT /ERR=32/ LLIGM
      SEGACT /ERR=32/ LLIGL
      GOTO 31
  32  CONTINUE
      IF (LSGDES) GOTO 3
**    write(6,*) 'desactivation-1 ',1,il1-1
      lsgdes=.true.
      do it=il1-1,1,-1
        lig1=ILIGNS.ILIGR(it)
        segdes lig1
        lig1=LILIGN.ilign(it)
        segdes lig1
        lig1=MILIGN.ilign(it)
        segdes lig1
      enddo
      SEGACT /ERR=3/LLIGM
      SEGACT /ERR=3/LLIGL
  31  CONTINUE
** essai pv
      NVMAX=MAX(LLIGM.NJMAX,LLIGL.NJMAX)
      if (ngmpet) then
        if (nvmax.gt.2*njmaxp.and.i-il1+1.gt.nbthrs/2) then
**        write(6,*) 'fin de bloc il1  i nvmax',il1,i,nvmax
          goto 3
        endif
      endif
      NJMAXP=NVMAX
** fin essai
C
      NA1  =LLIGM.IMMMM(/1)
      NVAL1=LLIGM.IMMMM(NA1)-LLIGM.LDEB(1)
      LDEB1=LLIGM.LDEB(1)
      NA2  =LLIGL.IMMMM(/1)
      NVAL2=LLIGL.IMMMM(NA2)-LLIGL.LDEB(1)
      LDEB2=LLIGL.LDEB(1)
C
C     Il peut arriver que les segments ligne et colonne 
C     d'un meme noeud ne demarrent pas de la meme colonne/ligne.
C     Decalage pour que ce ne soit plus le cas
      LDEBM=MIN(LDEB1,LDEB2)
      IF      (LDEB1-LDEBM.NE.0) THEN
        SEGACT,LLIGM*MOD
        LLIGM.NJMAX=LLIGL.NJMAX
        DO IMV=1,NA1
          LLIGM.LDEB(IMV)=LDEBM
        ENDDO
      ELSE IF (LDEB2-LDEBM.NE.0) THEN
        SEGACT,LLIGL*MOD
        LLIGL.NJMAX=LLIGM.NJMAX
        DO IMV=1,NA2
          LLIGL.LDEB(IMV)=LDEBM
        ENDDO
      ENDIF
C
      NVALL=0
      NVALT=MAX(NVAL1,NVAL2)
      if (ngmpet.and.nvalt.lt.nvallp/3) then
        write(6,*) 'passage en lent'
        ngmpet=.false.
      endif
      NVALLP=NVALT
      if (ngmpet) then
C       nvall=(NVMAX*nvstor)/nvstoc
        NVALL=NVMAX
        if (nvall.gt.NVSTRM) then
C       write(*,*)'redimensionnement de xreser',nvstrm,(4*NVALL),il1,il2
          nvstrm=4*nvall
          segsup,xreser
          segini,xreser
        endif
      endif
C
C     Le masque est uniquement present dans LIGM
      NA=NA1
      NBPAR=NA+1
      LMASQ=MASQA(NVALT+MASDIM)
      SEGINI /ERR=33/ LIGM
C
      NA=NA2
      NBPAR=NA+1
      LMASQ=0
      SEGINI /ERR=33/ LIGL
C
      GOTO 34
  33  CONTINUE
      IF (LSGDES) GOTO 3
**    write(6,*) 'desactivation-2 ',1,il1-1
      lsgdes=.true.
      do it=il1-1,1,-1
        lig1=ILIGNS.ILIGR(it)
        segdes lig1
        lig1=LILIGN.ilign(it)
        segdes lig1
        lig1=MILIGN.ilign(it)
        segdes lig1
      enddo
      IF (LIGM.EQ.0) THEN
        NA=NA1
        NBPAR=NA+1
        LMASQ=MASQA(NVAL1+MASDIM)
        SEGINI /ERR=3/ LIGM
      ENDIF
      IF (LIGL.EQ.0) THEN
        NA=NA2
        NBPAR=NA+1
        LMASQ=0
        SEGINI /ERR=3/ LIGL
      ENDIF
  34  CONTINUE
C
C     Travail sur la partie superieure
      MILIGN=IILIGS
      LLIGN =LLIGM
      LIGN  =LIGM
      ILIGNS.ILIGR(I)=LIGM
C
  35  CONTINUE
C
C     Recuperer la longueur du segment
      NA    =LIGN.IMMM(/1)
      LGLIG =LGLIG  + (NA*(NJMAX/NA)**1.333333333333333333)
      NVSTOC=NVSTOC + LLIGN.NJMAX
      NVAOR =NVAOR  + LLIGN.XXVA(/1)
C
C  ****  DECOMPACTAGE
C
      IPA=1
      DO 121 JPA=1,NA
        KPA =LLIGN.IPPO(JPA+1)-LLIGN.IPPO(JPA)
        IPP =LLIGN.IPPO(JPA)
        LPA =LLIGN.LDEB(JPA)
        LPA1=LPA-IPA
        LIGN.IMMM(JPA) =LPA
        LIGN.IVPO(JPA) =IPA
        LIGN.IPPVV(JPA)=IPA-1
        DO 122 MPA=1,KPA
          IPLA=LLIGN.LINC(MPA+IPP)-LPA1
          ICC =IPLA-IPA+1
          ICC1=ICC-MASQD(ICC)+1
          MICC=MASQA(ICC)
          IMSQ=LIGM.IMASQ(MICC)
          IF (IMSQ.LE.0) THEN
            MSQH=icc1
            MSQB=icc1
          ELSE
            MSQH=MAX(MASQH(IMSQ),icc1)
            MSQB=MIN(MASQB(IMSQ),icc1)
          ENDIF
          LIGM.IMASQ(MICC)=MASQV(MSQB,MSQH)
122     CONTINUE
        IPA=IPA+LLIGN.IMMMM(JPA)-LPA + 1
        IF(IMIN.GT.MILIGN.IPNO(LPA)) IMIN=MILIGN.IPNO(LPA)
121   CONTINUE
      LIGN.IPPVV(NA+1)=IPA-1
C
      IF (NA.GT.0) THEN
        LIGN.IPREL=LLIGN.IMMMM(1)
        LIGN.IDERL=LLIGN.IMMMM(NA)
        MILIGN.LCARA(2,I)=LIGN.IPREL
        MILIGN.LCARA(3,I)=LIGN.IDERL
      ENDIF
C
C     Travail sur la partie inferieure
      IF (MILIGN.EQ.IILIGS) THEN
        MILIGN=IILIGN
        LLIGN =LLIGL
        LIGN  =LIGL
        GOTO 35
      ENDIF
C
C     Indexation de imasq
      LMASQ=LIGM.IMASQ(/1)
      ipln=lmasq
      idec=1
      do 123 ipl=lmasq,1,-1
        if (LIGM.imasq(ipl).gt.0) then
          ipln=ipl-1
          idec=MASQB(LIGM.IMASQ(IPL))
        else
          LIGM.IMASQ(IPL)=-(IPLN*MASDIM+IDEC)
        endif
 123  continue
CCC   write (6,*) ' imasq ',lmasq
CCC   write (6,*) (LIGM.imasq(ipl),ipl=1,lmasq)
C
C*    write (6,*) 'longueur ligne ',nvall
C nb de ligne multiple du nb de threads
C blocage ligne lecture-ecriture pour minimiser le cache
C on note si on est au minimum de lignes
      nbthro=min(nbthrs,lglig/800+1)
      nbthro2=min(nbthrs,lglig/10000+1)
      IF (.not.ngmpet) THEN
        if (i+1-il1.ge.nbthro2) then
          nbthro2=min(nbthrs,i+1-il1)
        endif
        if (i+1-il1.ge.nbthro) then
          nbthro=min(nbthrs,i+1-il1)
          il2=i
**        write(6,*) ' nbthro il1 il2 ',nbthro,il1,il2
          GOTO 4
        endif
      ENDIF

    2 CONTINUE
      IL2=INO
      GOTO 4
    3 CONTINUE
      IL2=I-1
**    write(6,*) 'desactivation-4 ',llign
      SEGDES LLIGM,LLIGL
      IF (LIGL.NE.0) SEGSUP,LIGL
      IF (LIGM.NE.0) SEGSUP,LIGM
    4 CONTINUE
      nbthro=min(nbthrs,nbthro)
      nbthro2=min(nbthrs,nbthro2)
      nbthr=nbthro
      if(xreser.ne.0) then
CC      write(6,*)'segsup xreser 1',il1,il2
        segsup xreser
        xreser=0
      endif
C
      IF(IL2.GE.IL1) GOTO 40
C
C ****  MESSAGE PAS ASSEZ DE PLACE MEMOIRE
C
      CALL ERREUR(48)
**    call ooodmp(0)
      if (ithrd.eq.1) then
       call threadis
       call oooprl(0)
      endif
      call ooomru(0)
      RETURN
   40 CONTINUE
C
C     Travail sur la partie superieure
      MILIGN=IILIGS
C
   39 CONTINUE
C
      IM=INC
      DO 352 IH=IL2,IL1,-1
        LIGN=ILIGNS.ILIGR(IH)
        IL=INC
        DO 354 JH=1,LIGN.IMMM(/1)
          IM=MIN(IM,LIGN.IMMM(JH))
          IL=MIN(IL,LIGN.IMMM(JH))
 354    CONTINUE
        LIGN.IML=IL
        MILIGN.LCARA(1,IH)=LIGN.IML
        IMMT(IH)=MILIGN.IPNO(IM)
  352 CONTINUE
C
C     Travail sur la partie inferieure
      IF (MILIGN.EQ.IILIGS) THEN
        MILIGN=IILIGN
        GOTO 39
      ENDIF
C
      IMAX=IL1-1
      IF (LSGDES) THEN
        DO IT=IMIN,IL1-1
          LIG1=ILIGNS.ILIGR(IT)
          SEGACT /ERR=333/ LIG1
        ENDDO
  333   CONTINUE
        IMAX=IT-1
      ENDIF
C
      nbthrp=nbthr
      nbthr=min(nbthrp,il2-il1+1)
      nbthr=nbthro2 
      IPER=IMIN
      IDER=IMAX
      IL1P=IL1
      DO 13 I=IL1,IL2
C
        MILIGN=IILIGS
  336   CONTINUE
C
C       Travail sur la partie inferieure
        LILIGN=IILIGN
        DO ITH=1,NBTHR-1
          CALL THREADID(ITH,SHOLE3I)
        ENDDO
        CALL SHOLE3I(NBTHR)
        DO ITH=NBTHR-1,1,-1
          CALL THREADIF(ITH)
        ENDDO
C
C       Travail sur la partie superieure
        LILIGN=IILIGS
        DO ITH=1,NBTHR-1
          CALL THREADID(ITH,SHOLE3I)
        ENDDO
        CALL SHOLE3I(NBTHR)
        DO ITH=NBTHR-1,1,-1
          CALL THREADIF(ITH)
        ENDDO
C
        IF (I.EQ.IL1P.AND.IDER.NE.IL1P-1) THEN
          DO IT=IDER,IPER,-1
            LIG1=ILIGNS.ILIGR(IT)
            SEGDES LIG1
          ENDDO
          IPER=IDER+1
          IDER=IL1-1
          DO IT=IPER,IDER
            LIG1=ILIGNS.ILIGR(IT)
            SEGACT /ERR=335/ LIG1
          ENDDO
  335     CONTINUE
          IDER=IT-1
          IF (IDER.LT.IPER) THEN
            CALL ERREUR(48)
            RETURN
          ENDIF
          GOTO 336
        ENDIF
C
        NBTHR=1
C
C       Travail sur la partie superieure
        MILIGN=IILIGS
        LLIGN=MILIGN.ILIGN(I)
        LIGN =ILIGNS.ILIGR(I)
        NA=LIGN.IDERL-LIGN.IPREL+1
        CALL SOMPAC(IPPVV(1),IMASQ(1),NA,KIVPO(1),KIVLO(1),NBPAR,IZROSF)
C
        NVALL=KIVLO(NBPAR)-1
        LMASQ=MASQA(LLIGN.NJMAX)+1
C
  444   CONTINUE
C
C       LIGN est deja dimensionne a NJMAX
        IF (ngmpet) THEN
          SEGADJ,LIGN
          LIG5=LIGN
          GOTO 17
        ENDIF
C
        SEGINI /ERR=109/ LIG5
        GOTO 108
 109    CONTINUE
        IF (LSGDES) GOTO 14
C       write(6,*) 'desactivation-3 ',1,il1p-1
        lsgdes=.true.
        do it=il1p-1,1,-1
          lig1=ILIGNS.ILIGR(it)
          segdes lig1
          lig1=MILIGN.ilign(it)
          segdes lig1
        enddo
        SEGINI /ERR=14/ LIG5
 108    CONTINUE
        LIG5.IML  =LIGN.IML
        LIG5.IMM  =NBPAR
        LIG5.IPREL=LIGN.IPREL
        LIG5.IDERL=LIGN.IDERL
        DO 111 LHE=1,NA
          LIG5.IMMM(LHE)=LIGN.IMMM(LHE)
 111    CONTINUE
        DO 112 LHF=1,NA+1
          LIG5.IPPVV(LHF)=LIGN.IPPVV(LHF)
 112    CONTINUE
C
 17     CONTINUE
C
        DO 113 LHG=1,NBPAR
          LIG5.IVPO(2*LHG-1)=KIVPO(LHG)
          LIG5.IVPO(2*LHG)  =KIVLO(LHG)
 113    CONTINUE
C
        DO IT=NBPAR-1,1,-1
          IVI=LIG5.IVPO(2*IT)
          IVF=LIG5.IVPO(2*(IT+1))-1
          ICI=LIG5.IVPO(2*IT-1)
          DO IC=MASQA(IVF-IVI+ICI),MASQA(ICI),-1
            LIG5.IMASQ(IC)=IT
          ENDDO
        ENDDO
C
        IPA=1
        DO 119 JPA=1,NA
          IPP=LLIGN.IPPO(JPA)
          KPA=LLIGN.IPPO(JPA+1)-IPP
          LPA=LLIGN.LDEB(JPA)
          LPA1=LPA-IPA
          DO 120 MPA=1,KPA
            XXV=LLIGN.XXVA(MPA+IPP)
            IF (ABS(XXV).LE.XPETIT) GOTO 120
            IPLA=LLIGN.LINC(MPA+IPP)-LPA1
            IMSQ=LIG5.IMASQ(MASQA(IPLA))
            IF (IMSQ.GT.0) THEN
114           CONTINUE
              IF (IPLA.GE.LIG5.IVPO(2*(IMSQ+1)-1)) THEN
                IMSQ=IMSQ+1
                GOTO 114
              ENDIF
              IDBC=LIG5.IVPO(2*IMSQ-1)
              IDBV=LIG5.IVPO(2*IMSQ  )
              IPOS=IDBV+(IPLA-IDBC)
              LIG5.VAL(IPOS)=XXV
            ENDIF
120       CONTINUE
          IPA=IPA+LLIGN.IMMMM(JPA)-LPA+1
119     CONTINUE
C
        ITR=NBPAR
        DO 125 IPL=LMASQ,1,-1
          IF (LIG5.IMASQ(IPL).GT.0) THEN
            ITR=LIG5.IMASQ(IPL)
          ELSE
            LIG5.IMASQ(IPL)=-ITR
          ENDIF
 125    CONTINUE
C
        SEGSUP,LLIGN
        MILIGN.ILIGN(I)=LIG5
        ILIGNS.ILIGR(I)=0
C
C       Travail sur la partie inferieure
        IF (MILIGN.EQ.IILIGS) THEN
          MILIGN=IILIGN
          LLIGN=MILIGN.ILIGN(I)
          GOTO 444
        ENDIF
        IF (.not.ngmpet) SEGSUP,LIGN
C
        IDER=I
        IPER=IDER
        IL1=I+1
  13  CONTINUE
      GOTO 15
  14  CONTINUE
      LIGN=ILIGNS.ILIGR(I)
      SEGSUP,LIGN
      IL2=I-1
  15  CONTINUE
      IL1=IL1P
C
      if (lsgdes) then
C       write(6,*) 'desactivation-4 ',1,il1-1
        do it=il1-1,imin,-1
          lig1=ILIGNS.ILIGR(it)
          if (lig1.ne.0) segdes lig1
        enddo
      endif
      nbthr=nbthrp
C
C  ****  BOUCLE  *5*  TRAVAILLE SUR LE NOEUD I QUI EST EN LECTURE
C
      ipos=0
      iper=imin
      ider=imin-1
      DO 5 I=IMIN,IL2

      LIGM=LILIGN.ILIGN(I)
      LIGL=MILIGN.ILIGN(I)
      IF(I.LT.IL1) GOTO 7
C
C  ******* LE NOEUD I EST EN MEMOIRE IL EST TRIANGULE JUSQU'A
C  ******* IPREL  IL FAUT CONTINUER TOUTE LES LIGNES PUIS CALCULER
C  ******* LE TERME DIAGONAL
C
      LIGN=LIGM
      NBCOL=LIGN.IVPO(2*LIGN.IPPVV(2)-1)-1
      DO 156 KHG=1,LIGN.IMMM(/1)
C
C       Travail sur la partie superieure
        LIGN=LIGM
        II=LIGN.IPREL-1+KHG
        LIGN.IMMM(KHG)=0
C
        MM=LIGN.IVPO(2*LIGN.IPPVV(KHG))
        NN=LIGN.IVPO(2*LIGN.IPPVV(KHG+1))
        N=NN-MM
        NN=NN-1
C
        DIAG(II)=LIGN.VAL(NN)
        IF (N.EQ.1) THEN
          IF (II-NBNNMC.GT.0) THEN
            RE(II-NBNNMC,II-NBNNMC,1)=LIGN.VAL(NN)
            GOTO 41
          ENDIF
          GOTO 8
        ENDIF
C
        LIGN.VAL(NN)=LIGN.VAL(NN)+
     &               CHOLE5(LIGN,LIGL,KHG,DIAG(1),NBCOL,0,NBOP(1))
   8    CONTINUE
C
C       Travail sur la partie inferieure + terme diagonal
        LIGN=LIGL
        II=LIGN.IPREL-1+KHG
        LIGN.IMMM(KHG)=0
C
        MM=LIGN.IVPO(2*LIGN.IPPVV(KHG))
        NN=LIGN.IVPO(2*LIGN.IPPVV(KHG+1))
        N=NN-MM
        NN=NN-1
C
        LIGN.VAL(NN)=LIGN.VAL(NN)+
     &               CHOLE5(LIGN,LIGM,KHG,DIAG(1),NBCOL,1,NBOP(1))
        IF (II-NBNNMC.GT.0) THEN
          RE(II-NBNNMC,II-NBNNMC,1)=LIGN.VAL(NN)
          GOTO 41
        ENDIF
C
C       Terme diagonal ok?
        diagref=max(abs(diag(ii)),diagmin)
        diagcmp=diagref*5d-12
        IF (ABS(LIGN.VAL(NN)).GT.diagcmp) GOTO 12
C
C  il faut mettre une valeur plus grande sur les LX car on a un probleme de conditionnement
C  sur le calcul des reactions en cas de 2 relations presque identique
C
C  **** ON VIENT DE DETECTER UN MODE D'ENSEMBLE
C  **** ON AJOUTE A LA STRUCTURE UN RESSORT EGAL A CELUI QUI EXISTAIT
C  **** AU PREALABLE SUR CETTE INCONNUE.
C
*     write (6,*) ' ldmts mode d ensemble ittr ligne ',
*    >      ittr(ii),ii,diag(ii),val(nn),diagref
C  on garde le signe car il fau un moins sur les ML
        vmaxi=diagref
        do ipv=MM,NN
          vmaxi=max(vmaxi,abs(lign.val(ipv)))
        enddo
        if(ittr(ii).NE.0) then
          lign.VAL(NN)=lign.val(nn)-1.30D0*diagref
          NENSLX=NENSLX+1
        else
          lign.val(nn)=vmaxi
        endif
        NENS=NENS+1
        LIGL.IMMM(KHG)=NENS
        LIGM.IMMM(KHG)=NENS
   12   CONTINUE
C
C       Stabilisation
        IF (ISTAB.NE.0) THEN
C
C         Elimination des Nan
          if(.not.(abs(LIGN.val(nn)).lt.xgrand*xzprec).and.pasfait) then
            pasfait=.false.
            LIGN.val(nn)=xgrand*xzprec
            write (6,*) 'Nan dans ldmts ligne',ii
          endif
C
          diagcmp=abs(diagmin)*1d-5+xpetit/xzprec
          IF (ITTR(II).EQ.0) THEN
            IF (LIGN.VAL(NN).LT.-DIAGMAX*1D-3) THEN
              LIGN.VAL(NN)=ABS(LIGN.VAL(NN))
***           val(nn)=diagmax*1d-3
              IF (LIGN.IMMM(KHG).EQ.0) NENS=NENS+1
              LIGN.IMMM(KHG)=NENS
            ELSEIF (LIGN.VAL(NN).LE.DIAGCMP) THEN
              IF (ISR.EQ.0.OR.IIMPI.GT.0)
     &        write (ioimp,*) 'stabilisation RESO',ii,val(nn),diag(ii)
              isr=isr+1
              LIGN.VAL(NN)=MAX(DIAGCMP,-LIGN.VAL(NN))
              IF (LIGN.IMMM(KHG).EQ.0) NENS=NENS+1
              LIGN.IMMM(KHG)=NENS
            ELSE
              DIAGMIN=MIN(DIAGMIN,ABS(LIGN.VAL(NN)))
            ENDIF
          ELSE IF (ITTR(II).NE.0) THEN
            IF (LIGN.VAL(NN).GE.ABS(DIAG(II))*STMULT) THEN
              IF (ISRL.EQ.0.OR.IIMPI.GT.0)
     &        write (ioimp,*) 'stabilisation RESO lagrange',ii,val(nn)
              isrl=isrl+1
              LIGN.VAL(NN)=-ABS(DIAG(II))*STMULT
              IF (LIGN.IMMM(KHG).EQ.0) THEN
                NENS=NENS+1
                NENSLX=NENSLX+1
              ENDIF
              LIGN.IMMM(KHG)=NENS
            ENDIF
          ENDIF
        ENDIF
C
        DIAG(II)=LIGN.VAL(NN)
        IF(ABS(DIAG(II)).LE.XPETIT) THEN
          DIAG(II)=diagmax
          IF(ITTR(II).NE.0) DIAG(II)=-DIAGMAX
          LIGL.VAL(NN)=DIAG(II)
        ENDIF
        LIGM.VAL(NN)=DIAG(II)
        GOTO 41
C
C  ****  ERREUR MATRICE SINGUIERE
C
        INTERR(1)=I
        CALL ERREUR(49)
        if (ithrd.eq.1) then
          call threadis
         call oooprl(0)
        endif
        call ooomru(0)
        RETURN
C
C  ****  ON COMPTE LE NOMBRE DE TERMES DIAGONAUX NEGATIFS
C        ET LE NOMBRE DE MULTIPLICATEUR DE LAGRANGE
C
   41   CONTINUE
        IF(DIAG(II).LT.0.D0) INEG=INEG+1
        IF(ITTR(II).NE.0)   NBLAG=NBLAG+1
C
        CONDMAX=MAX(CONDMAX,ABS(DIAG(II)))
        IF (II.LE.NBNNMC) THEN
          CONDMIN=MIN(CONDMIN,ABS(DIAG(II)))
          DIAG(II)=1.D0/DIAG(II)
        ENDIF
C
  156 CONTINUE
C
C     Suppression du masque
      LMASQ=0
C
      LIGN=LIGL
      NVALL=LIGN.VAL(/1)
      NBPAR=LIGN.IVPO(/1)/2
      NA   =LIGN.IMMM(/1)
      SEGADJ,LIGN
C
      LIGN=LIGM
      NVALL=LIGN.VAL(/1)
      NBPAR=LIGN.IVPO(/1)/2
      NA   =LIGN.IMMM(/1)
      SEGADJ,LIGN
C
      NVSTOR=NVSTOR+NVALL
      NVSTRM=MAX(NVSTRM,NVALL*4)
C
      NVALL=0
      SEGINI,LIG3
      lig3.iml=iml
      lig3.imm=imm
      lig3.iprel=iprel
      lig3.iderl=iderl
      do ipv=1,immm(/1)
       lig3.immm(ipv)=immm(ipv)
      enddo
      do ipv=1,ippvv(/1)
       lig3.ippvv(ipv)=ippvv(ipv)
      enddo
      do ipv=1,ivpo(/1)
       lig3.ivpo(ipv)=ivpo(ipv)
      enddo
      ILIGNS.ILIGR(I)=LIG3
      IF (lsgdes) SEGDES,LIG3
C
      IF (I.EQ.IDBGL) THEN
        LIG1=LILIGN.ILIGN(I)
        WRITE(*,*) 'LIGNE FINALE SUP',I,LIG1
        segprt,lig1
        LIG1=MILIGN.ILIGN(I)
        WRITE(*,*) 'LIGNE FINALE INF',I,LIG1
        segprt,lig1
      ENDIF
C
C  ****  ON TRIANGULARISE LES AUTRES LIGNES
C
      IL1=IL1+1
      IF (IL1.GT.IL2) GOTO 5
      LIGM=LILIGN.ILIGN(I)
      LIGL=MILIGN.ILIGN(I)
      lign=MILIGN.ilign(il1)
      GOTO 7
C
  71  CONTINUE
C
      IF (IPER.GT.IDER) THEN
        call erreur(48)
        write(6,*) ' 1 ldmts iper ider i il1 il2 ',
     >         iper,ider,i,il1,il2
        do ipv=1,ino
          lig1=MILIGN.ilign(ipv)
          call oooeta(lig1,ieta,imod)
          if (ipv.lt.il1.or.ipv.gt.il2) then
**         if (ieta.eq.1) write(6,*) 'ligne active ',ipv
          endif
        enddo
**      call ooodmp(0)
        if (ithrd.eq.1) then
          call threadis
          call oooprl(0)
        endif
        call ooomru(0)
        RETURN
      ENDIF
C
C     Passage en superlent
      if (ider.lt.il1-1.and..not.ngmpet) then
        write(6,*) 'passage en superlent',ider,i,il1-1
        ngmpet=.true.
      endif
C
C     soit parce qu'on a fini, soit parce qu'on manque de memoire
C     il faut executer les lignes activees puis les desactiver
C     lancer les chole4 et attendre qu'ils soient finis
      if (xreser.ne.0) then
CC      write(6,*)'segsup xreser 2'              
        segsup xreser
        xreser=0
      endif

      if (ipos.ne.0) then
**      write (6,*) ' lancement thread ',iper,ider,il1,il2
        nbthr=min(nbthr,il2-il1+1)
        ipers=iper
        iders=ider
        ipert=iper
        ifois=1
*       write(6,*) 'il2 kidepb ider',lcara(1,il2),kidepb,ider
**      write(6,*) 'i il1 il2',i,il1,il2
 401    continue
        ipas=max(1200/ifois,80)
        if (ifois.ge.6 ) ipas= 80
        if (nbthr.eq.1) ipas=igrand
*402    continue
        liper=lcara(2,iper)
        iper2=il1-1
**      call ooosur(MILIGN.ilign(il2))
        do il=il1,il2
          lign=MILIGN.ilign(il)
          kidepb=lcara(1,il)-1
          iperC=liper-kidepb
          iperC=max(1,iperC)
          imsq=LIGN.imasq(masqa(iperC))
          if (imsq.lt.0) then
CC          iperC=max(LIGN.ivpo(2*(-imsq)-1),iperC+1)
            iperC=LIGN.ivpo(2*(-imsq)-1)
          endif
          iper1=ipno(iperC+kidepb)
          iper2=min(iper1,iper2)
        enddo
        iper=iper2
        ipert=iper
        ider=min(iders,iper+ipas-1)
        lider=lcara(3,ider)
        ider2=il1-1
        isaug=igrand
        if (il2.lt.il1) then
          write(6,*) 'il1 > il2',il1,il2
          call erreur(48)
        endif
        do il=il1,il2
          lign=MILIGN.ilign(il)
          kidepb=lcara(1,il)-1
          iderC=lider-kidepb
          iderC=max(1,iderC)
          imsq=LIGN.imasq(masqa(iderC))
          isaut=0
          if (imsq.lt.0) then
            iderC=LIGN.ivpo(2*(-imsq)-1)-1
**          isaut=0
**          do mpv=masqa(iderC),masqa(iper-kidepb),-1
**           if (LIGN.imasq(mpv).lt.0) isaut=isaut+masdim
**          enddo
            isaut=ipas/2
**          isaut=-1
          endif
          ider1=ipno(iderC+kidepb)
          ider2=min(ider1,ider2)
          isaug=min(isaug,isaut)
        enddo
        ider=min(ider2+isaug,iders)
        ifois=ifois+1
C
C       Travail sur la partie superieure
        MILIGN=IILIGS
        LILIGN=IILIGN
        iths=nbthr
        do ith=1,nbthr
          if (ith.ne.iths) call threadid(ith,chole4i)
        enddo
        call chole4i(iths)
        do ith=nbthr,1,-1
          if (ith.ne.iths) call threadif(ith)
        enddo
C
C       Travail sur la partie inferieure
        MILIGN=IILIGN
        LILIGN=IILIGS
        do ith=1,nbthr
          if (ith.ne.iths) call threadid(ith,chole4i)
        enddo
        call chole4i(iths)
        do ith=nbthr,1,-1
          if (ith.ne.iths) call threadif(ith)
        enddo
C
        iper=ider+1
        if (iper.le.iders) goto 401
        iper=ipers
        ider=iders
      ENDIF

C  test ctrlC
      if (ierr.ne.0) goto 9999
      ipos=0
      IF (LSGDES) THEN
C       write(6,*) 'desactivation 5 iper ider',iper,ider
        do ipv=ider,iper,-1
          lig1=MILIGN.ilign(ipv)
          segdes lig1
          lig1=LILIGN.ilign(ipv)
          segdes lig1
        enddo
        IF (I.NE.IL1-1) THEN
          LIG1=MILIGN.ILIGN(I)
          SEGACT,LIG1*MOD
          LIG1=LILIGN.ILIGN(I)
          SEGACT,LIG1*MOD
        ENDIF
      ENDIF
      iper=ider+1
      GOTO 5
C
    7 CONTINUE
      IF (LSGDES.AND.(I.LT.IL1)) THEN
        IF (.not.ngmpet) THEN
          iperI=lcara(2,i)
          do il=il1,il2
            lign=MILIGN.ilign(il)
            kidepb=lcara(1,il)-1
            iperC=iperI-kidepb
            iperC=max(1,iperC)
            if (LIGN.imasq(masqa(iperC)).gt.0) goto 432
          enddo
          GOTO 5
        ENDIF
 432    CONTINUE
        SEGACT/ERR=71/LIGM
        SEGACT/ERR=71/LIGL
      ENDIF
      ipos=ipos+1
      ider=i
      IF (I.EQ.IL1-1) GOTO 71
C
    5 CONTINUE
C
**    write (6,*) ' il1 il2 apres 5 ',il1,il2
      if (lsgdes) then
C       write(6,*) 'desactivation 8 il1 il2',il1,il1p,il2
        DO I=IL2,IL1P,-1
          LIGN=MILIGN.ILIGN(I)
          SEGDES,LIGN
          LIGN=LILIGN.ILIGN(I)
          SEGDES,LIGN
        ENDDO
      ENDIF
      nbopt=0
      do ith=1,nbthro
       nbopt=nbopt+nbop(ith)
       nbop(ith)=0
      enddo
      nbopin=nbopt
      nbopit=nbopit+nbopin
      call timespv(ittime,oothrd)
      kcour=(ittime(1)+ittime(2))/10
      if (xreser.ne.0) then 
CC      write(6,*)'segsup xreser 3'              
        segsup xreser
        xreser=0
      endif
      IF(IL2.LT.INO) GOTO 1
C  ON MET A JOUR LE NOMBRE DE TERMES DIAGONAUX NEGATIF
C   ON ENLEVE LE NOMBRE DE MULTIPLICATEUR DE LAGRANGE
C       INEG=INEG-NBLAG
C     on ne compte pas 2 fois les multiplicateurs qui vont etre
C     elimines lors de la resolution car mode d'ensemble
      INEG=INEG-(NBLAG-NENSLX)
      if (iimpi.ne.0.and.NENSLX.gt.0) WRITE(IOIMP,4820) NENSLX
 4820 FORMAT(I12,' MODES D ENSEMBLE PORTANT SUR DES MULTIPLICATEURS',
     1' DE LAGRANGE DETECTES')

      IF (IIMPI.EQ.1) WRITE(IOIMP,4821) NVSTOC
 4821 FORMAT( ' NOMBRE DE VALEURS  DANS LE PROFIL',I12)
      IF (IIMPI.EQ.1) WRITE(IOIMP,4822) NVSTOR*2
 4822 FORMAT( ' NOMBRE DE VALEURS STOCKEES DANS LE PROFIL',I12)
      IF (IIMPI.EQ.1) WRITE(IOIMP,4825) Nbopit/1000000
 4825 FORMAT( ' NOMBRE DE GIGA OPERATIONS FMA',I40)
      IF (IIMPI.EQ.1) WRITE(IOIMP,4823) NVaor
 4823 FORMAT( ' NOMBRE DE VALEURS initiales',I9)
C     IF (IIMPI.EQ.1) WRITE(IOIMP,4824) nbopit/1d6/(kcour-kcouri)
C 4824 FORMAT( ' Performance en Gigaflops ',F8.1)
      INTERR(1)=NVSTOR*2
      REAERR(1)=NVSTOR/INC**(4.D0/3.D0)
      REAERR(2)=2.D0*NBOPIT/1.D6/MAX(1,(KCOUR-KCOURI))
      REAERR(3)=CONDMAX/CONDMIN
      IF (IPASV.EQ.0.or.reaerr(3).gt.1.D30) CALL ERREUR(-278)
      IPASV=1
      call ooomru(0)
      if(lsgdes) then
        do ipv=1,ino
          lign=MILIGN.ilign(ipv)
          segdes lign
          lign=LILIGN.ilign(ipv)
          segdes lign
        enddo
      endif
      SEGDES,MINCPO
      SEGDES,MIMIK
      SEGDES,MMATRI
      SEGDES,MILIGN,LILIGN
      MMATRX=MMATRI
      SEGSUP KIVPO,KIVLO
      do ipv=1,nnoe
        lig3=ILIGNS.ILIGR(ipv)
        segsup lig3
      enddo
      SEGSUP,ILIGNS
      segsup immt
C     write (6,*) ' ldmts ipos max ',iposm
 9999 continue
      if (ithrd.eq.1) then
        call threadis
       call oooprl(0)
      endif
      if(iimpi.ne.0)  write (6,*) 
     >      ' ldmts condmin condmax ',condmin,condmax,diag(/1)
      SEGDES,MDIAG
      RETURN
      END