C ELNLIN    SOURCE    AM        11/05/26    21:15:22     6982           
       subroutine ELNLIN (da, uel, force, XMKel,
     &             delta_prec, Ftr, delta,deltamax)

        IMPLICIT INTEGER(I-N)
        IMPLICIT REAL*8(A-H,O-Z)
*
        DIMENSION da(20), XMKel(5,5), uel(5), force(5), Ftr(5)
        DIMENSION Res(3), AA(3,3),AA1(3,3),DVM(3)
*
        a=da(8)
        b=da(9)
        c=da(10)
        d=da(11)
        e=da(12)
        f=da(13)
        a9=da(15)
*
        A1=2.5
        epsilon2=1E-12
C       Calcul de la limite élastique en décollement.
        uelu=sqrt(uel(3)**2+uel(5)**2)
        forceu=sqrt(force(3)**2+force(5)**2)
        XMKelu=1/sqrt(2.)*sqrt(XMKel(3,3)**2+XMKel(5,5)**2)
*
        if (a9.eq.1) then
        q_1=4.
        q_2=1.
        else
        if (a9.eq.2) then
        q_1=4.
        q_2=1.
        else
        q_1=6.
        q_2=2.
        endif
        endif
*
        fe=forceu**2.-(force(1)/q_1*EXP(-A1*force(1)))**2.
*
        if (fe.le.0) then
        call MATVE1 (XMKel, uel, 5,5, Ftr,2)
        delta=delta_prec
        else
*
        deltamax=1
*
                ff1=(force(1)*(EXP(-A1*force(1))
     & +q_2*(delta_prec*deltamax)))/q_1-forceu
                gg1=-(((delta_prec/1.)/(1-(delta_prec/1.))
     & +LOG(1.-(delta_prec/1.)))
     & *force(1)*EXP(-A1*force(1)))
     & /(2*q_1*XMKelu)+force(1)/XMKel(1,1)-uel(1)
        Xhh1=((delta_prec/1.)**2.*force(1)*EXP(-A1*force(1)))/q_1
     & +(1.-(delta_prec/1.))*(forceu-XMKelu*uelu)
        Res(1)=-ff1
        Res(2)=-gg1
        Res(3)=-Xhh1
        XNorm_Res=SQRT(Res(1)**2.+Res(2)**2.+Res(3)**2.)
                do while (XNorm_Res.gt.epsilon2)
                dff1dV=(EXP(-A1*force(1))+q_2*(delta_prec*deltamax))/q_1
     & -(A1*force(1)*EXP(-A1*force(1)))/q_1
                dff1dM=-1.
                dff1ddelta=1*deltamax*q_2*force(1)/q_1
                dgg1dV=(((delta_prec/1)/(1.-(delta_prec/1))
     & +LOG(1.-(delta_prec/1)))*A1*force(1)*EXP(-A1*force(1)))
     & /(2*q_1*XMKelu)
     & -(((delta_prec/1)/(1.-(delta_prec/1))
     & +LOG(1.-(delta_prec/1)))
     & *EXP(-A1*force(1)))/(2*q_1*XMKelu)+1./XMKel(1,1)
                dgg1dM=0
                dgg1ddelta=-1/1*(((delta_prec/1)*force(1)
     & *EXP(-A1*force(1)))/(2*q_1*(1.-(delta_prec/1))**2.*XMKelu))
                dhh1dV=((delta_prec/1)**2.*EXP(-A1*force(1)))/q_1
     & -((delta_prec/1)**2.*A1*force(1)*EXP(-A1*force(1)))/q_1
                dhh1dM=1.-(delta_prec/1)
                dhh1ddelta=1/1*(2*((delta_prec/1)*force(1)
     & *EXP(-A1*force(1)))/q_1-forceu+XMKelu*uelu)
                AA(1,1)=dff1dV
                AA(1,2)=dff1dM
                AA(1,3)=dff1ddelta
                AA(2,1)=dgg1dV
                AA(2,2)=dgg1dM
                AA(2,3)=dgg1ddelta
                AA(3,1)=dhh1dV
                AA(3,2)=dhh1dM
                AA(3,3)=dhh1ddelta
                call HINV33(AA,AA1)
                call MATVE1(AA1,Res,3,3,DVM,2)
                force(1)=force(1)+DVM(1)
                forceu=forceu+DVM(2)
                delta_prec=delta_prec+DVM(3)
*
                ff1=(force(1)*(EXP(-A1*force(1))
     & +q_2*(delta_prec*deltamax)))/q_1-forceu
                gg1=-(((delta_prec/1.)/(1.-(delta_prec/1.))
     & +LOG(1.-(delta_prec/1.)))
     & *force(1)*EXP(-A1*force(1)))
     & /(2*q_1*XMKelu)+force(1)/XMKel(1,1)-uel(1)
        Xhh1=((delta_prec/1.)**2.*force(1)*EXP(-A1*force(1)))/q_1
     & +(1.-(delta_prec/1.))*(forceu-XMKelu*uelu)
        Res(1)=-ff1
        Res(2)=-gg1
        Res(3)=-Xhh1
        XNorm_Res=SQRT(Res(1)**2.+Res(2)**2.+Res(3)**2.)
                enddo
        force(2)=XMKel(2,2)*uel(2)
        force(4)=XMKel(2,2)*uel(4)
        force(3)=uel(3)/uelu*forceu
        force(5)=uel(5)/uelu*forceu
        uelu1=sqrt(uel(3)**2+uel(5)**2)
       DO I=1,5
        Ftr(I)=force(I)
       ENDDO
        delta=delta_prec
        endif
        return
        end