c234567===================================================== SUB DMAT ==
      subroutine dmat(ym,pr,itdp,d)
C     d= matrice [3x3], sigma= d*epsilon
      dimension d(3,3)
      if(itdp.EQ.0) goto 100
      if (itdp.eq.1) goto 200
100   continue
      ctp=ym/(1-pr**2)
      g=ym/((1+pr)*2)
      d(1,1)=1
      d(1,2)=pr
      d(2,1)=pr
      d(2,2)=1
      do 10,I10=1,2
      do 20,I20=1,2
      d (I10,I20)=d(I10,I20)*ctp
20    continue
10    continue
      d(1,3)=0
      d(2,3)=0
      d(3,1)=0
      d(3,2)=0
      d(3,3)=g
      return
200   continue
      cdp =ym*(1-pr)/((1+pr)*(1-2*pr))
      cpr=pr/(1-pr)
      g=ym/((1+pr)*2)
      d(1,1)=1
      d(1,2)=cpr
      d(2,1)=cpr
      d(2,2)=1
      do 30,I30=1,2
      do 40,I40=1,2
      d(I30,I40)=d(I30,I40)*cdp
40    continue
30    continue
      d(1,3)=0
      d(2,3)=0
      d(3,1)=0
      d(3,2)=0
      d(3,3)=g
      return
      end