      subroutine melmit(grr,grrz,grrzz,grs,grsz,grszz,gss,gssz,gsszz,
     *                  gzz,zint,npiz,npi,enu,e,tt,l,lz,deltat,elas,
     *                  girr,giss,girs)
c    ...............................................................
c    calcul de la matrice d'elasticite
c    ..............................................................
c
      double precision elas(9),grr,grrz,grrzz,grs,grsz,grszz,
     *                 gss,gssz,gsszz, delta, girr,giss,girs ,
     *                 e,enu,tt(npi), gzz(npi),deltat
      double precision zint(npiz)

            delta = (grr+grrz*zint(lz)+grrzz*zint(lz)**2) *
     *              (gss+gssz*zint(lz)+gsszz*zint(lz)**2) -
     *              (grs+grsz*zint(lz)+grszz*zint(lz)**2)**2
            girr =  (gss+gssz*zint(lz)+gsszz*zint(lz)**2) /delta
            giss =  (grr+grrz*zint(lz)+grrzz*zint(lz)**2) /delta
            girs =  -(grs+grsz*zint(lz)+grszz*zint(lz)**2) /delta

               

            deltat = sqrt(delta*gzz(l))

c
c           matrice d'elasticite 2 blocs diagonaux
c           --------------------------------------

c           c11 c12 c13                el(1)  el(2)  el(4)
c           c12 c22 c23                       el(3)  el(5) 
c           c13 c23 c33          =                   el(6)
c                      c44 c45                            el(7) el(8)
c                      c45 c55                                  el(9)


            elas(1) = e/(2*(1+enu)) * (2*girr*girr + 
     *                                 2*enu/(1-enu) * girr*girr)
            elas(2) = e/(2*(1+enu)) * (2*girs*girs +
     *                                 2*enu/(1-enu) * girr*giss)
            elas(4) = e/(2*(1+enu)) * (girr*girs + girs*girr +
     *                                 2*enu/(1-enu) * girr*girs)
            elas(3) = e/(2*(1+enu)) * (giss*giss + giss*giss +
     *                                 2*enu/(1-enu) * giss*giss)
            elas(5) = e/(2*(1+enu)) * (girs*giss + giss*girs +
     *                                 2*enu/(1-enu) * giss*girs)
            elas(6) = e/(2*(1+enu)) * (girr*giss + girs*girs +
     *                                 2*enu/(1-enu) * girs*girs)
            elas(7) = (4/tt(l)**2) * e/(2*(1+enu)) * girr 
            elas(8) = (4/tt(l)**2) * e/(2*(1+enu)) * girs 
            elas(9) = (4/tt(l)**2) * e/(2*(1+enu)) * giss 

       end
