o1 = 2.00000000000000d0*phi(i,j,k)
        o2 = cos(o1)
        o3 = cos(phi(i,j,k))
        o4 = o3**2
        o5 = c*o4
        o6 = 1.00000000000000d0 + o5
        o7 = -eta0
        o8 = eta(i,j,k) + o7
        o9 = o8**2
        o10 = sigma**2
        o11 = 1/o10
        o12 = -(o11*o9)
        o13 = exp(o12)
        o14 = eta(i,j,k) + eta0
        o15 = o14**2
        o16 = -(o11*o15)
        o17 = exp(o16)
        o18 = o13 + o17
        o19 = sin(q(i,j,k))
        o20 = o19**n
        o21 = 2.00000000000000d0*amp*o18*o20*o6
        o22 = exp(o21)
        o23 = -1.00000000000000d0 + o22
        o24 = -eta(i,j,k)
        o25 = exp(o24)
        o26 = sin(phi(i,j,k))
        o27 = o26**2
        o28 = 1/o19
        o29 = sin(o1)
        o30 = cos(q(i,j,k))
        o31 = o30**2
        o32 = -1.00000000000000d0 + n
        o33 = o19**o32
        o34 = amp*n*o18*o22*o3*o31*o33*o6
        o35 = 2.00000000000000d0*amp*c*o18*o22*o27*o3*o33
        o36 = 4.0000000000000d0*eta(i,j,k)*eta0*o11
        o37 = exp(o36)
        o38 = eta(i,j,k)*o37
        o39 = -(eta0*o37)
        o40 = eta(i,j,k) + eta0 + o38 + o39
        o41 = o16 + o21
        o42 = exp(o41)
        o43 = 1.00000000000000d0 + n
        o44 = o19**o43
        o45 = -2.00000000000000d0*amp*o11*o3*o40*o42*o44*o6
        o46 = -2.00000000000000d0*amp*c*o18*o22*o26*o33*o4
        o47 = amp*n*o18*o22*o26*o31*o33*o6
        o48 = -2.00000000000000d0*amp*o11*o26*o40*o42*o44*o6
        o49 = c*o2
        o50 = 2.00000000000000d0 + c + o49
        o51 = -1.00000000000000d0 + o37
        o52 = -2.00000000000000d0*eta0*o51
        o53 = 1.00000000000000d0 + o37
        o54 = 2.00000000000000d0*eta(i,j,k)*o53
        o55 = n*o10*o53
        o56 = o52 + o54 + o55
        o57 = o16 + o21 + o24
        o58 = exp(o57)
        o59 = o11*o9
        o60 = exp(o59)
        o61 = o11*o15
        o62 = exp(o61)
        o63 = o60 + o62
        o64 = o19**2
        gxx(i,j,k) = 5.0000000000000d-1*(1.00000000000000d0 + o22 + o2*o
     &  23)
c        dxgxx(i,j,k) = 5.0000000000000d-1*o25*(1.00000000000000d0*o22*o2
c     &  6*o28*o29 - 2.00000000000000d0*o27*o28*o3 + 2.00000000000000d0*a
c     &  mp*c*o18*o2*o22*o27*o3*o33 + o34 + o35 + o45 + amp*n*o18*o2*o22*
c     &  o3*o31*o33*o6 - 2.00000000000000d0*amp*o11*o2*o3*o40*o42*o44*o6)
c        dygxx(i,j,k) = 5.0000000000000d-1*o25*(1.00000000000000d0*o28*o2
c     &  9*o3 - 2.00000000000000d0*o22*o26*o28*o4 - 2.00000000000000d0*am
c     &  p*c*o18*o2*o22*o26*o33*o4 + o46 + o47 + o48 + amp*n*o18*o2*o22*o
c     &  26*o31*o33*o6 - 2.00000000000000d0*amp*o11*o2*o26*o40*o42*o44*o6
c     &  )
c        dzgxx(i,j,k) = -5.0000000000000d-1*amp*o11*o20*o30*o4*o50*o56*o5
c     &  8
        gxy(i,j,k) = 5.0000000000000d-1*o23*o29
c        dxgxy(i,j,k) = 5.0000000000000d-1*o25*(1.00000000000000d0*o2*o26
c     &  *o28 - o2*o22*o26*o28 + 2.00000000000000d0*amp*c*o18*o22*o27*o29
c     &  *o3*o33 + amp*n*o18*o22*o29*o3*o31*o33*o6 - 2.00000000000000d0*a
c     &  mp*o11*o29*o3*o40*o42*o44*o6)
c        dygxy(i,j,k) = 5.0000000000000d-1*o25*(-(o2*o28*o3) + 1.00000000
c     &  000000d0*o2*o22*o28*o3 - 2.00000000000000d0*amp*c*o18*o22*o26*o2
c     &  9*o33*o4 + amp*n*o18*o22*o26*o29*o31*o33*o6 - 2.00000000000000d0
c     &  *amp*o11*o26*o29*o40*o42*o44*o6)
c        dzgxy(i,j,k) = -5.0000000000000d-1*amp*o11*o20*o29*o30*o56*o58*o
c     &  6
        gxz(i,j,k) = 0
c        dxgxz(i,j,k) = 0
c        dygxz(i,j,k) = 0
c        dzgxz(i,j,k) = 0
        gyy(i,j,k) = 5.0000000000000d-1*(1.00000000000000d0 + o22 - o2*o
     &  23)
c        dxgyy(i,j,k) = 5.0000000000000d-1*o25*(1.00000000000000d0*o26*o2
c     &  8*o29 - 2.00000000000000d0*o22*o27*o28*o3 - 2.00000000000000d0*a
c     &  mp*c*o18*o2*o22*o27*o3*o33 + o34 + o35 + o45 - amp*n*o18*o2*o22*
c     &  o3*o31*o33*o6 + 2.00000000000000d0*amp*o11*o2*o3*o40*o42*o44*o6)
c        dygyy(i,j,k) = 5.0000000000000d-1*o25*(1.00000000000000d0*o22*o2
c     &  8*o29*o3 - 2.00000000000000d0*o26*o28*o4 + 2.00000000000000d0*am
c     &  p*c*o18*o2*o22*o26*o33*o4 + o46 + o47 + o48 - amp*n*o18*o2*o22*o
c     &  26*o31*o33*o6 + 2.00000000000000d0*amp*o11*o2*o26*o40*o42*o44*o6
c     &  )
c        dzgyy(i,j,k) = -5.0000000000000d-1*amp*o11*o20*o27*o30*o50*o56*o
c     &  58
        gyz(i,j,k) = 0
c        dxgyz(i,j,k) = 0
c        dygyz(i,j,k) = 0
c        dzgyz(i,j,k) = 0
        gzz(i,j,k) = o22
c        dxgzz(i,j,k) = amp*o11*o3*o33*(n*o10*o31*o6*o63 + 2.000000000000
c     &  00d0*(c*o10*o27*o63 - o40*o6*o60*o64))*exp(o12 + o16 + o21 + o24
c     &  )
c        dygzz(i,j,k) = amp*o26*o33*(-2.00000000000000d0*c*o18*o4 + n*o18
c     &  *o31*o6 - 2.00000000000000d0*o11*o17*o40*o6*o64)*exp(o21 + o24)
c        dzgzz(i,j,k) = -(amp*o11*o20*o30*o56*o58*o6)