C NTAPE6    SOURCE    CHAT      05/01/13    02:02:22     5004
      SUBROUTINE NTAPE6(MCP,MCQ,IVMINU,IVMINL,IVMAXU,IVMAXL,IVLAMB,
     * M,N,NVD,IVFP,IVFQ,MVDU,MVDL,IVB,IVD,IVN,II,KK,IVDR,IDVD,
     * NDR,TERMIN,IVLL,IVUL,IPBASE)
      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8(A-H,O-Z)
      LOGICAL TERMIN
-INC TMXMAT
-INC SMLENTI

-INC PPARAM
-INC CCOPTIO
-INC SMLREEL
      POINTEUR MLREE4.MLREEL,MLREE5.MLREEL,MLREE6.MLREEL,MLREE7.MLREEL
      IVXL1=0
      IVXU1=0
      IVU1=0
      IVN1=0
      IVD1=0
      IVGM1=0
      IVGE1=0
      IVLAM1=0
      N11 = N + 1
      II=-1
      KK=-1
      MLREEL=IVDR
      JG=N11
      SEGINI MLREE1,MLREE2
      IVPZ=MLREE2
      IVQZ=MLREE1
      MXMAT=MCQ
      CALL MATVE1(XMAT,PROG,M,N11,MLREE1.PROG,1)
      MXMAT=MCP
      CALL MATVE1(XMAT,PROG,M,N11,MLREE2.PROG,1)
      ALPMAX=1.E25
      MLREE1=IVLAMB
      IPB=0
      DO 1 I=1,M
      IF(PROG(I).LT.-1.D-20)   THEN
          IF(ALPMAX+(MLREE1.PROG(I)/PROG(I)).GT.0.D0) THEN
            ALPMAX= -( MLREE1.PROG(I)/PROG(I))
            IPB=I
          ENDIF
      ENDIF
    1 CONTINUE
      IF(IIMPI.EQ.1799)
     *WRITE(IOIMP,FMT='('' ALPHAMAX  IPB = '',E12.5,2X,I3)')ALPMAX,IPB
      NDIS=0
      MLENTI=IDVD
      DO 7 I=1,NVD
         NDIS=NDIS+LECT(I)-1
    7 CONTINUE
      JG=(NDIS-NVD)+2*(N11-NVD)+2
      SEGINI MLREE6,MLENT1,MLENT2
      MLREE6.PROG(1)=0.D0
      MLENT1.LECT(1)=-3
      MLENT2.LECT(1)=0
*
*  CALCUL DES ALPHA CRITIQUES
*  ON COMMENCE PAR LES VARIABLES CONTINUES
*
      IP=1
      MLREEL=IVD
      MLREE1=IVN
      MLREE2=IVMINU
      MLREE3=IVMINL
      MLREE4=IVQZ
      MLREE5=IVPZ
      DO 2 I=NVD+1,N11
         CN=PROG(I)*(MLREE2.PROG(I)**2)
         CN=CN-(MLREE1.PROG(I)*(MLREE3.PROG(I)**2))
         CD=(MLREE3.PROG(I) ** 2) * MLREE5.PROG(I)
         CD=CD-((MLREE2.PROG(I)**2)*MLREE4.PROG(I))
         IP=IP+1
         IF(CD.EQ.0.D0) CD=1.D-35
         MLREE6.PROG(IP)=CN/CD
         MLENT1.LECT(IP)=I-NVD
         MLENT2.LECT(IP)=-1
    2 CONTINUE
      MLREE2=IVMAXU
      MLREE3=IVMAXL
      DO 3 I=NVD+1,N11
         CN=PROG(I)*(MLREE2.PROG(I)**2)
         CN=CN-(MLREE1.PROG(I)*(MLREE3.PROG(I)**2))
         CD=(MLREE3.PROG(I) ** 2) * MLREE5.PROG(I)
         CD=CD-((MLREE2.PROG(I)**2)*MLREE4.PROG(I))
         IP=IP+1
         IF(CD.EQ.0.D0) CD=1.D-35
         MLREE6.PROG(IP)=CN/CD
         MLENT1.LECT(IP)=I-NVD
         MLENT2.LECT(IP)=-1
    3 CONTINUE
*
* POUR LES VARIABLES DISCRETES
*
      IF(NVD.NE.0) THEN
       IF(IIMPI.EQ.1799)
     * WRITE(IOIMP,FMT='('' IL Y A'',I4,'' VARIABLES DISCRETES'')')NVD
        MXMAT=MVDU
        MXMA1=MVDL
        DO 4 I=1,NVD
          DO 4 J = 3,LECT(I)
            CN=PROG(I)*XMAT(I,J)*XMAT(I,J-1)
            CN=CN-(MLREE1.PROG(I)*MXMA1.XMAT(I,J)*MXMA1.XMAT(I,J-1))
            CD=MLREE4.PROG(I)* MXMA1.XMAT(I,J)*MXMA1.XMAT(I,J-1)
            CD=CD-(MLREE5.PROG(I)*XMAT(I,J)*XMAT(I,J-1))
            IP=IP+1
            IF(CD.EQ.0.D0) CD=1.D-35
            MLREE6.PROG(IP)= CN/CD
            IF(ABS(MLREE6.PROG(IP)).LT.1D-10)MLREE6.PROG(IP)=0.D0
            MLENT1.LECT(IP)=I
            MLENT2.LECT(IP)=J
    4   CONTINUE
      ENDIF
      MAC2=MLENT1
      MAC3=MLENT2
      SEGSUP MLREE4,MLREE5
*
* ON ENTRE ALPHA MAX DANS LE TABLEAU DES ALPHA
*
      IP=IP+1
      MLREE6.PROG(IP)=ALPMAX
*
      MLREEL=MLREE6
*
*  TRI DES ALPHA
*
      JG=PROG(/1)
      SEGINI MLREE1
      SEGINI MLENT1,MLENT2
      MORDRE=MLENT1
      DO 5 I=1,JG
         MLENT1.LECT(I)=I
   5  CONTINUE
      CALL TRIFLO(PROG,MLREE1.PROG,MLENT1.LECT,MLENT2.LECT,JG)
      IF(IIMPI.EQ.1799) WRITE(IOIMP,FMT='('' LISTE DES ALPHA APRES TRI
     * '',/,(1X,5E12.5))')(PROG(I),I=1,IP)
      SEGSUP MLREE1,MLENT2
*
*   ELIMINATION DES ALPHA NON ADMISSIBLES
*
      IDL=0
      IFL=IP+2
      DO 6 I=1,IP
         IF(PROG(I).LT.1.D-20) IDL=I
         IF(PROG(I).GT.ALPMAX)THEN
             IF(IFL.EQ.IP+2) IFL=I-1
             MLENT1.LECT(I)=-2
         ENDIF
   6  CONTINUE
      IFL=MIN(IFL,IP)
      IF(IIMPI.EQ.1799) WRITE(IOIMP,FMT='('' LISTE DES ALPHA APRES TRI
     * ET ELIMINATION '',/,(1X,5E12.5))')(PROG(I),I=IDL,IFL)
*
*   RECHERCHE DU PAS ALPHA OPTIMUM
*
      MLREE1=IVLAMB
      MLREE2=IVDR
      SEGINI,MLREE3=MLREE1
      IVLAM1=MLREE3
*
*  CALCUL DU PREMIER POINT
*
      I=IDL
      INPINF=IDL
      QQZ=PROG(I)
      DO 8 J=1,MLREE1.PROG(/1)
         MLREE3.PROG(J)=QQZ*MLREE2.PROG(J)+MLREE1.PROG(J)
    8 CONTINUE
      IF(IIMPI.EQ.1799) WRITE(IOIMP,FMT='('' VALEUR DU LAMBDA POUR ALPHA
     * =0 '',/,(1X,5E12.5))')(MLREE3.PROG(I),I=1,M)
      CALL NTAPE1(MCP,MCQ,IVFP,IVFQ,IVLAM1,NVD,M,N,MVDU,MVDL,IVMINU,
     * IVMINL,IVMAXU,IVMAXL,IVU1,IVN1,IVD1,IVUL,IVLL,IVXU1,IVXL1)
      CALL NTAPE2(MCP,MCQ,IVXU1,IVXL1,IVB,N,M,IVGE1,
     *IVGM1,IVLAM1,IPBASE)
      XLP1=0.D0
      MLREE4=IVGE1
      DO 9 J=1,MLREE1.PROG(/1)
      XLP1=MLREE2.PROG(J)*MLREE4.PROG(J) +XLP1
    9 CONTINUE
      IF(IIMPI.EQ.1799) WRITE(IOIMP,FMT='(''  1ER ALPHA , PROD SCAL :''
     *  ,1X,E12.5,'' , '',E12.5)')QQZ,XLP1
         MLREE4=IVXU1
         MLREE5=IVU1
         MLREE6=IVN1
         SEGSUP MLREE4,MLREE5,MLREE6
         MLREE4=IVD1
         MLREE5=IVGM1
         MLREE6=IVGE1
         MLREE7=IVXL1
         SEGSUP MLREE4,MLREE5,MLREE6,MLREE7
      IVXL1=0
      IVXU1=0
      IVU1=0
      IVN1=0
      IVD1=0
      IVGM1=0
      IVGE1=0
*
* CALCUL SUR LE DERNIER POINT
*
      INPSUP=IFL
      I=IFL
      QQZ=PROG(I)
      DO 10 J=1,MLREE1.PROG(/1)
         MLREE3.PROG(J)=QQZ*MLREE2.PROG(J)+MLREE1.PROG(J)
   10 CONTINUE
      CALL NTAPE1(MCP,MCQ,IVFP,IVFQ,IVLAM1,NVD,M,N,MVDU,MVDL,IVMINU,
     * IVMINL,IVMAXU,IVMAXL,IVU1,IVN1,IVD1,IVUL,IVLL,IVXU1,IVXL1)
      CALL NTAPE2(MCP,MCQ,IVXU1,IVXL1,IVB,N,M,IVGE1,
     *IVGM1,IVLAM1,IPBASE)
      XLP2=0.D0
      MLREE4=IVGE1
      DO 11 J=1,MLREE1.PROG(/1)
         XLP2=MLREE2.PROG(J)*MLREE4.PROG(J) +XLP2
   11 CONTINUE
      IF(IIMPI.EQ.1799) WRITE(IOIMP,FMT='('' DERNIER ALPHA , PROD SCAL
     * :'',1X,E12.5,'' , '',E12.5)')QQZ,XLP2
*
* RECHERCHE  PAR DICHOTOMIE SUR LES INDICES pour trouver l'interval
* de alpha interessant
*
      IF(XLP1 * XLP2.GE.0.D0) THEN
          ALPHA=ALPMAX
          IF (ALPMAX.LT.1E25) THEN
             KK = -3
             II=IPB
          ENDIF
          GO TO 20
      ENDIF
   12 CONTINUE
      IK=(INPINF+INPSUP)/2
      IF(IK.EQ.INPINF) GO TO 15
         MLREE4=IVXU1
         MLREE5=IVU1
         MLREE6=IVN1
         SEGSUP MLREE4,MLREE5,MLREE6
         MLREE4=IVD1
         MLREE5=IVGM1
         MLREE6=IVGE1
         MLREE7=IVXL1
         SEGSUP MLREE4,MLREE5,MLREE6,MLREE7
      IVXL1=0
      IVXU1=0
      IVU1=0
      IVN1=0
      IVD1=0
      IVGM1=0
      IVGE1=0
      I=IK
      QQZ=PROG(I)
      DO 13 J=1,MLREE1.PROG(/1)
         MLREE3.PROG(J)=QQZ*MLREE2.PROG(J)+MLREE1.PROG(J)
   13 CONTINUE
      CALL NTAPE1(MCP,MCQ,IVFP,IVFQ,IVLAM1,NVD,M,N,MVDU,MVDL,IVMINU,
     * IVMINL,IVMAXU,IVMAXL,IVU1,IVN1,IVD1,IVUL,IVLL,IVXU1,IVXL1)
      CALL NTAPE2(MCP,MCQ,IVXU1,IVXL1,IVB,N,M,IVGE1,
     *IVGM1,IVLAM1,IPBASE)
      XLP3=0.D0
      MLREE4=IVGE1
      DO 14 J=1,MLREE1.PROG(/1)
         XLP3=MLREE2.PROG(J)*MLREE4.PROG(J) +XLP3
   14 CONTINUE
      IF(XLP1*XLP3.LE.0.D0) THEN
         XLP2=XLP3
         INPSUP=I
      ELSE
         XLP1=XLP3
         INPINF=I
      ENDIF
      GO TO 12
   15 CONTINUE
*
* RECHERCHE PAR DICHOTOMIE pour determiner la valeur de alpha donnant le
* min de l(lambda + alpha * grad)
*
      IPP=0
      ALPE=PROG(INPINF)
      ALGR=PROG(INPSUP)
      ALPHA=(ALGR + ALPE)/2.D0
   60 CONTINUE
      IPP=IPP + 1
      IF(IPP.GT.100) THEN
        CALL ERREUR ( 603)
        RETURN
       ENDIF
      ALPHA0=ALPHA
         MLREE4=IVXU1
         MLREE5=IVU1
         MLREE6=IVN1
         SEGSUP MLREE4,MLREE5,MLREE6
         MLREE4=IVD1
         MLREE5=IVGM1
         MLREE6=IVGE1
         MLREE7=IVXL1
         SEGSUP MLREE4,MLREE5,MLREE6,MLREE7
      IVXL1=0
      IVXU1=0
      IVU1=0
      IVN1=0
      IVD1=0
      IVGM1=0
      IVGE1=0
*      IF(IIMPI.EQ.1799) WRITE(IOIMP,FMT='('' ALPHA '',E12.5)')ALPHA
      QQZ=ALPHA
      DO 17 J=1,MLREE1.PROG(/1)
         MLREE3.PROG(J)=QQZ*MLREE2.PROG(J)+MLREE1.PROG(J)
  17  CONTINUE
      CALL NTAPE1(MCP,MCQ,IVFP,IVFQ,IVLAM1,NVD,M,N,MVDU,MVDL,IVMINU,
     *    IVMINL,IVMAXU,IVMAXL,IVU1,IVN1,IVD1,IVUL,IVLL,IVXU1,IVXL1)
      CALL NTAPE2(MCP,MCQ,IVXU1,IVXL1,IVB,N,M,IVGE1,
     *IVGM1,IVLAM1,IPBASE)
      XLP3=0.D0
      MLREE4=IVGE1
      DO 50 J=1,MLREE1.PROG(/1)
         XLP3=MLREE2.PROG(J)*MLREE4.PROG(J) +XLP3
   50 CONTINUE
      IF(XLP1*XLP3.LE.0.D0) THEN
         XLP2=XLP3
         ALGR=ALPHA
      ELSE
         XLP1=XLP3
         ALPE=ALPHA
      ENDIF
      ALPHA =(ALGR + ALPE)/2
      EPSI=1.D-10
      IF(ALPHA.LT.1.D-12.AND.NDR.GE.1)THEN
         TERMIN=.TRUE.
         GO TO 20
      ENDIF
      IF((ABS(ALPHA-ALPHA0)/ALPHA).GE.EPSI)   GO TO 60
   20 CONTINUE
         MLENT1=MAC2
         MLENT2=MAC3
         MLENT3=MORDRE
       IF(IIMPI.EQ.1799)WRITE(IOIMP,FMT='('' ALPHA '',E12.5)')ALPHA
      IF (KK.NE.-3)THEN
*
****************** TEST SI ARRET SUR PLAN DE DISCONTINUITE ************
*
         ALPHA1=PROG(INPSUP)
         IF((ABS(ALPHA1-ALPHA)/ALPHA).LE.EPSI)THEN
              KK=MLENT2.LECT(MLENT3.LECT(INPSUP))
              II=MLENT1.LECT(MLENT3.LECT(INPSUP))
              ALPHA=ALPHA1
         ENDIF
         ALPHA2=PROG(INPINF)
         IF(((ABS(ALPHA2-ALPHA)/ALPHA).LE.EPSI.AND.ALPHA2.NE.0.D0)
     *     .OR.ALPHA.LT.1.D-12)THEN
              KK=MLENT2.LECT(MLENT3.LECT(INPINF))
              II=MLENT1.LECT(MLENT3.LECT(INPINF))
              ALPHA=ALPHA2
         ENDIF
      ENDIF
      DO 21 J=1,MLREE1.PROG(/1)
         MLREE3.PROG(J)=ALPHA * MLREE2.PROG(J)+MLREE1.PROG(J)
   21 CONTINUE
      IF (KK.EQ.-3)THEN
         MLREE3.PROG(II)=0.
      ENDIF
      ALPHA0=ALPHA
         MLREE4=IVXU1
         MLREE5=IVU1
         MLREE6=IVN1
         SEGSUP MLREE4,MLREE5,MLREE6
         MLREE4=IVD1
         MLREE5=IVGM1
         MLREE6=IVGE1
         MLREE7=IVXL1
         SEGSUP MLREE4,MLREE5,MLREE6,MLREE7
      IVLAMB=IVLAM1
      SEGSUP MLENT1,MLENT2,MLREEL,MLENT3
      RETURN
      END