TABLE OF CONTENTS


ABINIT/contistr01 [ Functions ]

[ Top ] [ Functions ]

NAME

 contistr01

FUNCTION

 Carries out specialized metric tensor operations needed for contraction
 of the 2nd strain derivative of the l=0,1,2,3 nonlocal Kleinman-Bylander
 pseudopotential operation.  Derivatives are wrt a pair of cartesian
 strain components.
 Full advantage is taken of the full permutational symmetry of these
 tensors.

COPYRIGHT

 Copyright (C) 1998-2018 ABINIT group (DRH)
 This file is distributed under the terms of the
 GNU General Public License, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .
 For the initials of contributors, see ~abinit/doc/developers/contributors.txt.

INPUTS

  istr=1,...6 specifies cartesian strain component 11,22,33,32,31,21
  rank=angular momentum
  gm(3,3)=metric tensor (array is symmetric but stored as 3x3)
  gprimd(3,3)=reciprocal space dimensional primitive translations
  aa(2,*)=unique elements of complex right-hand tensor
  bb(2,*)=unique elements of complex left-hand tensor

OUTPUT

  eisnl(3)=contraction for nonlocal internal strain derivative energy

NOTES

 All tensors are stored in a compressed storage mode defined below;
 input and output conform to this scheme.
 When tensor elements occur repeatedly due to symmetry, the
 WEIGHT IS INCLUDED in the output tensor element to simplify later
 contractions with other tensors of the same rank and form, i.e. the
 next contraction is then simply a dot product over the unique elements.

PARENTS

      nonlop_pl

CHILDREN

SOURCE

  47 #if defined HAVE_CONFIG_H
  48 #include "config.h"
  49 #endif
  50 
  51 #include "abi_common.h"
  52 
  53 
  54 module m_contistr01
  55 contains
  56 !!**
  57 
  58 subroutine contistr01(istr,rank,gm,gprimd,eisnl,aa,bb)
  59 
  60  use defs_basis
  61  use m_abicore
  62 
  63 !This section has been created automatically by the script Abilint (TD).
  64 !Do not modify the following lines by hand.
  65 #undef ABI_FUNC
  66 #define ABI_FUNC 'contistr01'
  67 !End of the abilint section
  68 
  69  implicit none
  70 
  71 !Arguments ------------------------------------
  72 !scalars
  73  integer,intent(in) :: istr,rank
  74 !arrays
  75  real(dp),intent(in) :: aa(2,((rank+1)*(rank+2))/2),bb(2,((rank+2)*(rank+3))/2)
  76  real(dp),intent(in) :: gm(3,3),gprimd(3,3)
  77  real(dp),intent(out) :: eisnl(3)
  78 
  79 !Local variables-------------------------------
  80 !scalars
  81  integer :: ii,jj,ka,kb
  82 !arrays
  83  integer,parameter :: idx(12)=(/1,1,2,2,3,3,3,2,3,1,2,1/)
  84  real(dp) :: dgm(3,3),tmp(2,3)
  85  real(dp),allocatable :: cm(:,:,:)
  86 
  87 ! *************************************************************************
  88 
  89  ABI_ALLOCATE(cm,(3,((rank+1)*(rank+2))/2,((rank+2)*(rank+3))/2))
  90 
  91  ka=idx(2*istr-1);kb=idx(2*istr)
  92 
  93  do ii = 1,3
  94    dgm(:,ii)=-(gprimd(ka,:)*gprimd(kb,ii)+gprimd(kb,:)*gprimd(ka,ii))
  95  end do
  96 
  97  cm(:,:,:)=0.d0
  98 !
  99 !The code below was written by a Mathematica program and formatted by
 100 !a combination of editing scripts.  It is not intended to be read
 101 !by human beings, and certainly not to be modified by one.  Conceivably
 102 !it could be shortened somewhat by identifying common subexpressions.
 103 !
 104  if(rank==1)then
 105    cm(1,1,1)=dgm(1,1)
 106    cm(1,2,1)=dgm(1,2)
 107    cm(1,3,1)=dgm(1,3)
 108    cm(1,1,5)=dgm(1,3)
 109    cm(1,2,5)=dgm(2,3)
 110    cm(1,3,5)=dgm(3,3)
 111    cm(1,1,6)=dgm(1,2)
 112    cm(1,2,6)=dgm(2,2)
 113    cm(1,3,6)=dgm(2,3)
 114    cm(2,1,2)=dgm(1,2)
 115    cm(2,2,2)=dgm(2,2)
 116    cm(2,3,2)=dgm(2,3)
 117    cm(2,1,4)=dgm(1,3)
 118    cm(2,2,4)=dgm(2,3)
 119    cm(2,3,4)=dgm(3,3)
 120    cm(2,1,6)=dgm(1,1)
 121    cm(2,2,6)=dgm(1,2)
 122    cm(2,3,6)=dgm(1,3)
 123    cm(3,1,3)=dgm(1,3)
 124    cm(3,2,3)=dgm(2,3)
 125    cm(3,3,3)=dgm(3,3)
 126    cm(3,1,4)=dgm(1,2)
 127    cm(3,2,4)=dgm(2,2)
 128    cm(3,3,4)=dgm(2,3)
 129    cm(3,1,5)=dgm(1,1)
 130    cm(3,2,5)=dgm(1,2)
 131    cm(3,3,5)=dgm(1,3)
 132  elseif(rank==2)then
 133    cm(1,1,1)=2*gm(1,1)*dgm(1,1)
 134    cm(1,2,1)=-0.5d0*gm(2,2)*dgm(1,1)+3*gm(1,2)*dgm(1,2)-0.5d0*gm(1,1)&
 135 &   *dgm(2,2)
 136    cm(1,3,1)=-0.5d0*gm(3,3)*dgm(1,1)+3*gm(1,3)*dgm(1,3)-0.5d0*gm(1,1)&
 137 &   *dgm(3,3)
 138    cm(1,4,1)=-gm(2,3)*dgm(1,1)+3*gm(1,3)*dgm(1,2)+3*gm(1,2)*dgm(1,3)&
 139 &   -gm(1,1)*dgm(2,3)
 140    cm(1,5,1)=2*(gm(1,3)*dgm(1,1)+gm(1,1)*dgm(1,3))
 141    cm(1,6,1)=2*(gm(1,2)*dgm(1,1)+gm(1,1)*dgm(1,2))
 142    cm(1,1,2)=-0.5d0*gm(2,2)*dgm(1,1)+3*gm(1,2)*dgm(1,2)-0.5d0*gm(1,1)&
 143 &   *dgm(2,2)
 144    cm(1,2,2)=2*gm(2,2)*dgm(2,2)
 145    cm(1,3,2)=-0.5d0*gm(3,3)*dgm(2,2)+3*gm(2,3)*dgm(2,3)-0.5d0*gm(2,2)&
 146 &   *dgm(3,3)
 147    cm(1,4,2)=2*(gm(2,3)*dgm(2,2)+gm(2,2)*dgm(2,3))
 148    cm(1,5,2)=3*gm(2,3)*dgm(1,2)-gm(2,2)*dgm(1,3)-gm(1,3)*dgm(2,2)&
 149 &   +3*gm(1,2)*dgm(2,3)
 150    cm(1,6,2)=2*(gm(2,2)*dgm(1,2)+gm(1,2)*dgm(2,2))
 151    cm(1,1,3)=-0.5d0*gm(3,3)*dgm(1,1)+3*gm(1,3)*dgm(1,3)-0.5d0*gm(1,1)&
 152 &   *dgm(3,3)
 153    cm(1,2,3)=-0.5d0*gm(3,3)*dgm(2,2)+3*gm(2,3)*dgm(2,3)-0.5d0*gm(2,2)&
 154 &   *dgm(3,3)
 155    cm(1,3,3)=2*gm(3,3)*dgm(3,3)
 156    cm(1,4,3)=2*(gm(3,3)*dgm(2,3)+gm(2,3)*dgm(3,3))
 157    cm(1,5,3)=2*(gm(3,3)*dgm(1,3)+gm(1,3)*dgm(3,3))
 158    cm(1,6,3)=-gm(3,3)*dgm(1,2)+3*gm(2,3)*dgm(1,3)+3*gm(1,3)*dgm(2,3)&
 159 &   -gm(1,2)*dgm(3,3)
 160    cm(1,1,4)=-gm(2,3)*dgm(1,1)+3*gm(1,3)*dgm(1,2)+3*gm(1,2)*dgm(1,3)&
 161 &   -gm(1,1)*dgm(2,3)
 162    cm(1,2,4)=2*(gm(2,3)*dgm(2,2)+gm(2,2)*dgm(2,3))
 163    cm(1,3,4)=2*(gm(3,3)*dgm(2,3)+gm(2,3)*dgm(3,3))
 164    cm(1,4,4)=3*gm(3,3)*dgm(2,2)+2*gm(2,3)*dgm(2,3)+3*gm(2,2)*dgm(3,3)
 165    cm(1,5,4)=3*gm(3,3)*dgm(1,2)+gm(2,3)*dgm(1,3)+gm(1,3)*dgm(2,3)&
 166 &   +3*gm(1,2)*dgm(3,3)
 167    cm(1,6,4)=gm(2,3)*dgm(1,2)+3*gm(2,2)*dgm(1,3)+3*gm(1,3)*dgm(2,2)&
 168 &   +gm(1,2)*dgm(2,3)
 169    cm(1,1,5)=2*(gm(1,3)*dgm(1,1)+gm(1,1)*dgm(1,3))
 170    cm(1,2,5)=3*gm(2,3)*dgm(1,2)-gm(2,2)*dgm(1,3)-gm(1,3)*dgm(2,2)&
 171 &   +3*gm(1,2)*dgm(2,3)
 172    cm(1,3,5)=2*(gm(3,3)*dgm(1,3)+gm(1,3)*dgm(3,3))
 173    cm(1,4,5)=3*gm(3,3)*dgm(1,2)+gm(2,3)*dgm(1,3)+gm(1,3)*dgm(2,3)&
 174 &   +3*gm(1,2)*dgm(3,3)
 175    cm(1,5,5)=3*gm(3,3)*dgm(1,1)+2*gm(1,3)*dgm(1,3)+3*gm(1,1)*dgm(3,3)
 176    cm(1,6,5)=3*gm(2,3)*dgm(1,1)+gm(1,3)*dgm(1,2)+gm(1,2)*dgm(1,3)&
 177 &   +3*gm(1,1)*dgm(2,3)
 178    cm(1,1,6)=2*(gm(1,2)*dgm(1,1)+gm(1,1)*dgm(1,2))
 179    cm(1,2,6)=2*(gm(2,2)*dgm(1,2)+gm(1,2)*dgm(2,2))
 180    cm(1,3,6)=-gm(3,3)*dgm(1,2)+3*gm(2,3)*dgm(1,3)+3*gm(1,3)*dgm(2,3)&
 181 &   -gm(1,2)*dgm(3,3)
 182    cm(1,4,6)=gm(2,3)*dgm(1,2)+3*gm(2,2)*dgm(1,3)+3*gm(1,3)*dgm(2,2)&
 183 &   +gm(1,2)*dgm(2,3)
 184    cm(1,5,6)=3*gm(2,3)*dgm(1,1)+gm(1,3)*dgm(1,2)+gm(1,2)*dgm(1,3)&
 185 &   +3*gm(1,1)*dgm(2,3)
 186    cm(1,6,6)=3*gm(2,2)*dgm(1,1)+2*gm(1,2)*dgm(1,2)+3*gm(1,1)*dgm(2,2)
 187    cm(2,1,2)=2*(gm(1,2)*dgm(1,1)+gm(1,1)*dgm(1,2))
 188    cm(2,2,2)=2*(gm(2,2)*dgm(1,2)+gm(1,2)*dgm(2,2))
 189    cm(2,3,2)=-gm(3,3)*dgm(1,2)+3*gm(2,3)*dgm(1,3)+3*gm(1,3)*dgm(2,3)&
 190 &   -gm(1,2)*dgm(3,3)
 191    cm(2,4,2)=gm(2,3)*dgm(1,2)+3*gm(2,2)*dgm(1,3)+3*gm(1,3)*dgm(2,2)&
 192 &   +gm(1,2)*dgm(2,3)
 193    cm(2,5,2)=3*gm(2,3)*dgm(1,1)+gm(1,3)*dgm(1,2)+gm(1,2)*dgm(1,3)&
 194 &   +3*gm(1,1)*dgm(2,3)
 195    cm(2,6,2)=3*gm(2,2)*dgm(1,1)+2*gm(1,2)*dgm(1,2)+3*gm(1,1)*dgm(2,2)
 196    cm(2,1,4)=2*(gm(1,3)*dgm(1,1)+gm(1,1)*dgm(1,3))
 197    cm(2,2,4)=3*gm(2,3)*dgm(1,2)-gm(2,2)*dgm(1,3)-gm(1,3)*dgm(2,2)&
 198 &   +3*gm(1,2)*dgm(2,3)
 199    cm(2,3,4)=2*(gm(3,3)*dgm(1,3)+gm(1,3)*dgm(3,3))
 200    cm(2,4,4)=3*gm(3,3)*dgm(1,2)+gm(2,3)*dgm(1,3)+gm(1,3)*dgm(2,3)&
 201 &   +3*gm(1,2)*dgm(3,3)
 202    cm(2,5,4)=3*gm(3,3)*dgm(1,1)+2*gm(1,3)*dgm(1,3)+3*gm(1,1)*dgm(3,3)
 203    cm(2,6,4)=3*gm(2,3)*dgm(1,1)+gm(1,3)*dgm(1,2)+gm(1,2)*dgm(1,3)&
 204 &   +3*gm(1,1)*dgm(2,3)
 205    cm(2,1,6)=2*gm(1,1)*dgm(1,1)
 206    cm(2,2,6)=-0.5d0*gm(2,2)*dgm(1,1)+3*gm(1,2)*dgm(1,2)-0.5d0*gm(1,1)&
 207 &   *dgm(2,2)
 208    cm(2,3,6)=-0.5d0*gm(3,3)*dgm(1,1)+3*gm(1,3)*dgm(1,3)-0.5d0*gm(1,1)&
 209 &   *dgm(3,3)
 210    cm(2,4,6)=-gm(2,3)*dgm(1,1)+3*gm(1,3)*dgm(1,2)+3*gm(1,2)*dgm(1,3)&
 211 &   -gm(1,1)*dgm(2,3)
 212    cm(2,5,6)=2*(gm(1,3)*dgm(1,1)+gm(1,1)*dgm(1,3))
 213    cm(2,6,6)=2*(gm(1,2)*dgm(1,1)+gm(1,1)*dgm(1,2))
 214    cm(2,1,7)=-0.5d0*gm(2,2)*dgm(1,1)+3*gm(1,2)*dgm(1,2)-0.5d0*gm(1,1)&
 215 &   *dgm(2,2)
 216    cm(2,2,7)=2*gm(2,2)*dgm(2,2)
 217    cm(2,3,7)=-0.5d0*gm(3,3)*dgm(2,2)+3*gm(2,3)*dgm(2,3)-0.5d0*gm(2,2)&
 218 &   *dgm(3,3)
 219    cm(2,4,7)=2*(gm(2,3)*dgm(2,2)+gm(2,2)*dgm(2,3))
 220    cm(2,5,7)=3*gm(2,3)*dgm(1,2)-gm(2,2)*dgm(1,3)-gm(1,3)*dgm(2,2)&
 221 &   +3*gm(1,2)*dgm(2,3)
 222    cm(2,6,7)=2*(gm(2,2)*dgm(1,2)+gm(1,2)*dgm(2,2))
 223    cm(2,1,8)=-0.5d0*gm(3,3)*dgm(1,1)+3*gm(1,3)*dgm(1,3)-0.5d0*gm(1,1)&
 224 &   *dgm(3,3)
 225    cm(2,2,8)=-0.5d0*gm(3,3)*dgm(2,2)+3*gm(2,3)*dgm(2,3)-0.5d0*gm(2,2)&
 226 &   *dgm(3,3)
 227    cm(2,3,8)=2*gm(3,3)*dgm(3,3)
 228    cm(2,4,8)=2*(gm(3,3)*dgm(2,3)+gm(2,3)*dgm(3,3))
 229    cm(2,5,8)=2*(gm(3,3)*dgm(1,3)+gm(1,3)*dgm(3,3))
 230    cm(2,6,8)=-gm(3,3)*dgm(1,2)+3*gm(2,3)*dgm(1,3)+3*gm(1,3)*dgm(2,3)&
 231 &   -gm(1,2)*dgm(3,3)
 232    cm(2,1,9)=-gm(2,3)*dgm(1,1)+3*gm(1,3)*dgm(1,2)+3*gm(1,2)*dgm(1,3)&
 233 &   -gm(1,1)*dgm(2,3)
 234    cm(2,2,9)=2*(gm(2,3)*dgm(2,2)+gm(2,2)*dgm(2,3))
 235    cm(2,3,9)=2*(gm(3,3)*dgm(2,3)+gm(2,3)*dgm(3,3))
 236    cm(2,4,9)=3*gm(3,3)*dgm(2,2)+2*gm(2,3)*dgm(2,3)+3*gm(2,2)*dgm(3,3)
 237    cm(2,5,9)=3*gm(3,3)*dgm(1,2)+gm(2,3)*dgm(1,3)+gm(1,3)*dgm(2,3)&
 238 &   +3*gm(1,2)*dgm(3,3)
 239    cm(2,6,9)=gm(2,3)*dgm(1,2)+3*gm(2,2)*dgm(1,3)+3*gm(1,3)*dgm(2,2)&
 240 &   +gm(1,2)*dgm(2,3)
 241    cm(3,1,3)=2*(gm(1,3)*dgm(1,1)+gm(1,1)*dgm(1,3))
 242    cm(3,2,3)=3*gm(2,3)*dgm(1,2)-gm(2,2)*dgm(1,3)-gm(1,3)*dgm(2,2)&
 243 &   +3*gm(1,2)*dgm(2,3)
 244    cm(3,3,3)=2*(gm(3,3)*dgm(1,3)+gm(1,3)*dgm(3,3))
 245    cm(3,4,3)=3*gm(3,3)*dgm(1,2)+gm(2,3)*dgm(1,3)+gm(1,3)*dgm(2,3)&
 246 &   +3*gm(1,2)*dgm(3,3)
 247    cm(3,5,3)=3*gm(3,3)*dgm(1,1)+2*gm(1,3)*dgm(1,3)+3*gm(1,1)*dgm(3,3)
 248    cm(3,6,3)=3*gm(2,3)*dgm(1,1)+gm(1,3)*dgm(1,2)+gm(1,2)*dgm(1,3)&
 249 &   +3*gm(1,1)*dgm(2,3)
 250    cm(3,1,4)=2*(gm(1,2)*dgm(1,1)+gm(1,1)*dgm(1,2))
 251    cm(3,2,4)=2*(gm(2,2)*dgm(1,2)+gm(1,2)*dgm(2,2))
 252    cm(3,3,4)=-gm(3,3)*dgm(1,2)+3*gm(2,3)*dgm(1,3)+3*gm(1,3)*dgm(2,3)&
 253 &   -gm(1,2)*dgm(3,3)
 254    cm(3,4,4)=gm(2,3)*dgm(1,2)+3*gm(2,2)*dgm(1,3)+3*gm(1,3)*dgm(2,2)&
 255 &   +gm(1,2)*dgm(2,3)
 256    cm(3,5,4)=3*gm(2,3)*dgm(1,1)+gm(1,3)*dgm(1,2)+gm(1,2)*dgm(1,3)&
 257 &   +3*gm(1,1)*dgm(2,3)
 258    cm(3,6,4)=3*gm(2,2)*dgm(1,1)+2*gm(1,2)*dgm(1,2)+3*gm(1,1)*dgm(2,2)
 259    cm(3,1,5)=2*gm(1,1)*dgm(1,1)
 260    cm(3,2,5)=-0.5d0*gm(2,2)*dgm(1,1)+3*gm(1,2)*dgm(1,2)-0.5d0*gm(1,1)&
 261 &   *dgm(2,2)
 262    cm(3,3,5)=-0.5d0*gm(3,3)*dgm(1,1)+3*gm(1,3)*dgm(1,3)-0.5d0*gm(1,1)&
 263 &   *dgm(3,3)
 264    cm(3,4,5)=-gm(2,3)*dgm(1,1)+3*gm(1,3)*dgm(1,2)+3*gm(1,2)*dgm(1,3)&
 265 &   -gm(1,1)*dgm(2,3)
 266    cm(3,5,5)=2*(gm(1,3)*dgm(1,1)+gm(1,1)*dgm(1,3))
 267    cm(3,6,5)=2*(gm(1,2)*dgm(1,1)+gm(1,1)*dgm(1,2))
 268    cm(3,1,8)=-gm(2,3)*dgm(1,1)+3*gm(1,3)*dgm(1,2)+3*gm(1,2)*dgm(1,3)&
 269 &   -gm(1,1)*dgm(2,3)
 270    cm(3,2,8)=2*(gm(2,3)*dgm(2,2)+gm(2,2)*dgm(2,3))
 271    cm(3,3,8)=2*(gm(3,3)*dgm(2,3)+gm(2,3)*dgm(3,3))
 272    cm(3,4,8)=3*gm(3,3)*dgm(2,2)+2*gm(2,3)*dgm(2,3)+3*gm(2,2)*dgm(3,3)
 273    cm(3,5,8)=3*gm(3,3)*dgm(1,2)+gm(2,3)*dgm(1,3)+gm(1,3)*dgm(2,3)&
 274 &   +3*gm(1,2)*dgm(3,3)
 275    cm(3,6,8)=gm(2,3)*dgm(1,2)+3*gm(2,2)*dgm(1,3)+3*gm(1,3)*dgm(2,2)&
 276 &   +gm(1,2)*dgm(2,3)
 277    cm(3,1,9)=-0.5d0*gm(2,2)*dgm(1,1)+3*gm(1,2)*dgm(1,2)-0.5d0*gm(1,1)&
 278 &   *dgm(2,2)
 279    cm(3,2,9)=2*gm(2,2)*dgm(2,2)
 280    cm(3,3,9)=-0.5d0*gm(3,3)*dgm(2,2)+3*gm(2,3)*dgm(2,3)-0.5d0*gm(2,2)&
 281 &   *dgm(3,3)
 282    cm(3,4,9)=2*(gm(2,3)*dgm(2,2)+gm(2,2)*dgm(2,3))
 283    cm(3,5,9)=3*gm(2,3)*dgm(1,2)-gm(2,2)*dgm(1,3)-gm(1,3)*dgm(2,2)&
 284 &   +3*gm(1,2)*dgm(2,3)
 285    cm(3,6,9)=2*(gm(2,2)*dgm(1,2)+gm(1,2)*dgm(2,2))
 286    cm(3,1,10)=-0.5d0*gm(3,3)*dgm(1,1)+3*gm(1,3)*dgm(1,3)-0.5d0*gm(1,1)&
 287 &   *dgm(3,3)
 288    cm(3,2,10)=-0.5d0*gm(3,3)*dgm(2,2)+3*gm(2,3)*dgm(2,3)-0.5d0*gm(2,2)&
 289 &   *dgm(3,3)
 290    cm(3,3,10)=2*gm(3,3)*dgm(3,3)
 291    cm(3,4,10)=2*(gm(3,3)*dgm(2,3)+gm(2,3)*dgm(3,3))
 292    cm(3,5,10)=2*(gm(3,3)*dgm(1,3)+gm(1,3)*dgm(3,3))
 293    cm(3,6,10)=-gm(3,3)*dgm(1,2)+3*gm(2,3)*dgm(1,3)+3*gm(1,3)*dgm(2,3)&
 294 &   -gm(1,2)*dgm(3,3)
 295  elseif(rank==3)then
 296    cm(1,1,1)=3*gm(1,1)**2*dgm(1,1)
 297    cm(1,2,1)=4.5d0*gm(1,2)**2*dgm(1,1)+9*gm(1,1)*gm(1,2)*dgm(1,2)&
 298 &   +gm(1,1)*(-3*gm(2,2)*dgm(1,1)-1.5d0*gm(1,1)*dgm(2,2))
 299    cm(1,3,1)=4.5d0*gm(1,3)**2*dgm(1,1)+9*gm(1,1)*gm(1,3)*dgm(1,3)&
 300 &   +gm(1,1)*(-3*gm(3,3)*dgm(1,1)-1.5d0*gm(1,1)*dgm(3,3))
 301    cm(1,4,1)=9*gm(1,2)*(gm(1,3)*dgm(1,1)+gm(1,1)*dgm(1,3))+gm(1,1)&
 302 &   *(-6*gm(2,3)*dgm(1,1)+9*gm(1,3)*dgm(1,2)-3*gm(1,1)*dgm(2,3))
 303    cm(1,5,1)=gm(1,1)*(6*gm(1,3)*dgm(1,1)+3*gm(1,1)*dgm(1,3))
 304    cm(1,6,1)=gm(1,1)*(6*gm(1,2)*dgm(1,1)+3*gm(1,1)*dgm(1,2))
 305    cm(1,7,1)=7.5d0*gm(1,2)**2*dgm(1,2)-1.5d0*gm(1,1)*gm(2,2)*dgm(1,2)&
 306 &   -1.5d0*gm(1,2)*(gm(2,2)*dgm(1,1)+gm(1,1)*dgm(2,2))
 307    cm(1,8,1)=7.5d0*gm(1,3)**2*dgm(1,2)+gm(1,1)*(-1.5d0*gm(3,3)*dgm(1,2)&
 308 &   -3*gm(2,3)*dgm(1,3))+gm(1,3)*(-3*gm(2,3)*dgm(1,1)+15*gm(1,2)&
 309 &   *dgm(1,3)-3*gm(1,1)*dgm(2,3))-1.5d0*gm(1,2)*(gm(3,3)*dgm(1,1)&
 310 &   +gm(1,1)*dgm(3,3))
 311    cm(1,9,1)=7.5d0*gm(1,2)**2*dgm(1,3)+gm(1,1)*(-3*gm(2,3)*dgm(1,2)&
 312 &   -1.5d0*gm(2,2)*dgm(1,3))+gm(1,3)*(-1.5d0*gm(2,2)*dgm(1,1)+15*gm(1,2)&
 313 &   *dgm(1,2)-1.5d0*gm(1,1)*dgm(2,2))-3*gm(1,2)*(gm(2,3)*dgm(1,1)&
 314 &   +gm(1,1)*dgm(2,3))
 315    cm(1,10,1)=7.5d0*gm(1,3)**2*dgm(1,3)-1.5d0*gm(1,1)*gm(3,3)*dgm(1,3)&
 316 &   -1.5d0*gm(1,3)*(gm(3,3)*dgm(1,1)+gm(1,1)*dgm(3,3))
 317    cm(1,1,2)=4.5d0*gm(1,2)**2*dgm(1,1)+9*gm(1,1)*gm(1,2)*dgm(1,2)&
 318 &   +gm(1,1)*(-3*gm(2,2)*dgm(1,1)-1.5d0*gm(1,1)*dgm(2,2))
 319    cm(1,2,2)=6*gm(2,2)**2*dgm(1,1)+3*gm(1,2)**2*dgm(2,2)+gm(2,2)&
 320 &   *(6*gm(1,2)*dgm(1,2)+12*gm(1,1)*dgm(2,2))
 321    cm(1,3,2)=7.5d0*gm(2,3)**2*dgm(1,1)-6*gm(1,2)*gm(3,3)*dgm(1,2)&
 322 &   -3*gm(1,3)**2*dgm(2,2)-1.5d0*gm(1,1)*gm(3,3)*dgm(2,2)+9*gm(1,2)&
 323 &   *gm(1,3)*dgm(2,3)+gm(2,3)*(9*gm(1,3)*dgm(1,2)+9*gm(1,2)*dgm(1,3)&
 324 &   +15*gm(1,1)*dgm(2,3))-3*gm(1,2)**2*dgm(3,3)+gm(2,2)*(-1.5d0*gm(3,3)&
 325 &   *dgm(1,1)-6*gm(1,3)*dgm(1,3)-1.5d0*gm(1,1)*dgm(3,3))
 326    cm(1,4,2)=12*gm(1,1)*gm(2,3)*dgm(2,2)+gm(1,2)*(6*gm(2,3)*dgm(1,2)&
 327 &   +3*gm(1,3)*dgm(2,2))+3*gm(1,2)**2*dgm(2,3)+gm(2,2)*(12*gm(2,3)&
 328 &   *dgm(1,1)+3*gm(1,3)*dgm(1,2)+3*gm(1,2)*dgm(1,3)+12*gm(1,1)*dgm(2,3))
 329    cm(1,5,2)=1.5d0*gm(1,2)**2*dgm(1,3)+gm(1,1)*(12*gm(2,3)*dgm(1,2)&
 330 &   -4.5d0*gm(2,2)*dgm(1,3))+gm(1,3)*(-4.5d0*gm(2,2)*dgm(1,1)+3*gm(1,2)&
 331 &   *dgm(1,2)-4.5d0*gm(1,1)*dgm(2,2))+12*gm(1,2)*(gm(2,3)*dgm(1,1)&
 332 &   +gm(1,1)*dgm(2,3))
 333    cm(1,6,2)=4.5d0*gm(1,2)**2*dgm(1,2)+7.5d0*gm(1,1)*gm(2,2)*dgm(1,2)&
 334 &   +7.5d0*gm(1,2)*(gm(2,2)*dgm(1,1)+gm(1,1)*dgm(2,2))
 335    cm(1,7,2)=gm(2,2)*(3*gm(2,2)*dgm(1,2)+6*gm(1,2)*dgm(2,2))
 336    cm(1,8,2)=1.5d0*gm(2,3)**2*dgm(1,2)-4.5d0*gm(1,2)*gm(3,3)*dgm(2,2)&
 337 &   +gm(2,3)*(12*gm(2,2)*dgm(1,3)+12*gm(1,3)*dgm(2,2)+3*gm(1,2)*dgm(2,3))&
 338 &   +gm(2,2)*(-4.5d0*gm(3,3)*dgm(1,2)+12*gm(1,3)*dgm(2,3)-4.5d0*gm(1,2)&
 339 &   *dgm(3,3))
 340    cm(1,9,2)=6*gm(2,2)**2*dgm(1,3)+3*gm(1,2)*gm(2,3)*dgm(2,2)+gm(2,2)&
 341 &   *(3*gm(2,3)*dgm(1,2)+12*gm(1,3)*dgm(2,2)+3*gm(1,2)*dgm(2,3))
 342    cm(1,10,2)=7.5d0*gm(2,3)**2*dgm(1,3)+gm(3,3)*(-1.5d0*gm(1,3)*dgm(2,2)&
 343 &   -3*gm(1,2)*dgm(2,3))+gm(2,3)*(-3*gm(3,3)*dgm(1,2)+15*gm(1,3)&
 344 &   *dgm(2,3)-3*gm(1,2)*dgm(3,3))-1.5d0*gm(2,2)*(gm(3,3)*dgm(1,3)&
 345 &   +gm(1,3)*dgm(3,3))
 346    cm(1,1,3)=4.5d0*gm(1,3)**2*dgm(1,1)+9*gm(1,1)*gm(1,3)*dgm(1,3)&
 347 &   +gm(1,1)*(-3*gm(3,3)*dgm(1,1)-1.5d0*gm(1,1)*dgm(3,3))
 348    cm(1,2,3)=7.5d0*gm(2,3)**2*dgm(1,1)-6*gm(1,2)*gm(3,3)*dgm(1,2)&
 349 &   -3*gm(1,3)**2*dgm(2,2)-1.5d0*gm(1,1)*gm(3,3)*dgm(2,2)+9*gm(1,2)&
 350 &   *gm(1,3)*dgm(2,3)+gm(2,3)*(9*gm(1,3)*dgm(1,2)+9*gm(1,2)*dgm(1,3)&
 351 &   +15*gm(1,1)*dgm(2,3))-3*gm(1,2)**2*dgm(3,3)+gm(2,2)*(-1.5d0*gm(3,3)&
 352 &   *dgm(1,1)-6*gm(1,3)*dgm(1,3)-1.5d0*gm(1,1)*dgm(3,3))
 353    cm(1,3,3)=6*gm(3,3)**2*dgm(1,1)+3*gm(1,3)**2*dgm(3,3)+gm(3,3)&
 354 &   *(6*gm(1,3)*dgm(1,3)+12*gm(1,1)*dgm(3,3))
 355    cm(1,4,3)=3*gm(1,3)**2*dgm(2,3)+gm(3,3)*(3*gm(1,2)*dgm(1,3)+12*gm(1,1)&
 356 &   *dgm(2,3))+gm(2,3)*(12*gm(3,3)*dgm(1,1)+6*gm(1,3)*dgm(1,3)+12*gm(1,1)&
 357 &   *dgm(3,3))+3*gm(1,3)*(gm(3,3)*dgm(1,2)+gm(1,2)*dgm(3,3))
 358    cm(1,5,3)=4.5d0*gm(1,3)**2*dgm(1,3)+7.5d0*gm(1,1)*gm(3,3)*dgm(1,3)&
 359 &   +7.5d0*gm(1,3)*(gm(3,3)*dgm(1,1)+gm(1,1)*dgm(3,3))
 360    cm(1,6,3)=1.5d0*gm(1,3)**2*dgm(1,2)+gm(1,1)*(-4.5d0*gm(3,3)*dgm(1,2)&
 361 &   +12*gm(2,3)*dgm(1,3))+gm(1,3)*(12*gm(2,3)*dgm(1,1)+3*gm(1,2)&
 362 &   *dgm(1,3)+12*gm(1,1)*dgm(2,3))-4.5d0*gm(1,2)*(gm(3,3)*dgm(1,1)&
 363 &   +gm(1,1)*dgm(3,3))
 364    cm(1,7,3)=7.5d0*gm(2,3)**2*dgm(1,2)-1.5d0*gm(1,2)*gm(3,3)*dgm(2,2)&
 365 &   +gm(2,3)*(-3*gm(2,2)*dgm(1,3)-3*gm(1,3)*dgm(2,2)+15*gm(1,2)*dgm(2,3))&
 366 &   +gm(2,2)*(-1.5d0*gm(3,3)*dgm(1,2)-3*gm(1,3)*dgm(2,3)-1.5d0*gm(1,2)&
 367 &   *dgm(3,3))
 368    cm(1,8,3)=6*gm(3,3)**2*dgm(1,2)+3*gm(1,3)*gm(2,3)*dgm(3,3)+gm(3,3)&
 369 &   *(3*gm(2,3)*dgm(1,3)+3*gm(1,3)*dgm(2,3)+12*gm(1,2)*dgm(3,3))
 370    cm(1,9,3)=1.5d0*gm(2,3)**2*dgm(1,3)+gm(3,3)*(-4.5d0*gm(1,3)*dgm(2,2)&
 371 &   +12*gm(1,2)*dgm(2,3))+gm(2,3)*(12*gm(3,3)*dgm(1,2)+3*gm(1,3)&
 372 &   *dgm(2,3)+12*gm(1,2)*dgm(3,3))-4.5d0*gm(2,2)*(gm(3,3)*dgm(1,3)&
 373 &   +gm(1,3)*dgm(3,3))
 374    cm(1,10,3)=gm(3,3)*(3*gm(3,3)*dgm(1,3)+6*gm(1,3)*dgm(3,3))
 375    cm(1,1,4)=9*gm(1,2)*(gm(1,3)*dgm(1,1)+gm(1,1)*dgm(1,3))+gm(1,1)&
 376 &   *(-6*gm(2,3)*dgm(1,1)+9*gm(1,3)*dgm(1,2)-3*gm(1,1)*dgm(2,3))
 377    cm(1,2,4)=12*gm(1,1)*gm(2,3)*dgm(2,2)+gm(1,2)*(6*gm(2,3)*dgm(1,2)&
 378 &   +3*gm(1,3)*dgm(2,2))+3*gm(1,2)**2*dgm(2,3)+gm(2,2)*(12*gm(2,3)&
 379 &   *dgm(1,1)+3*gm(1,3)*dgm(1,2)+3*gm(1,2)*dgm(1,3)+12*gm(1,1)*dgm(2,3))
 380    cm(1,3,4)=3*gm(1,3)**2*dgm(2,3)+gm(3,3)*(3*gm(1,2)*dgm(1,3)+12*gm(1,1)&
 381 &   *dgm(2,3))+gm(2,3)*(12*gm(3,3)*dgm(1,1)+6*gm(1,3)*dgm(1,3)+12*gm(1,1)&
 382 &   *dgm(3,3))+3*gm(1,3)*(gm(3,3)*dgm(1,2)+gm(1,2)*dgm(3,3))
 383    cm(1,4,4)=9*gm(2,3)**2*dgm(1,1)+18*gm(1,2)*gm(3,3)*dgm(1,2)+9*gm(1,3)&
 384 &   **2*dgm(2,2)+15*gm(1,1)*gm(3,3)*dgm(2,2)-6*gm(1,2)*gm(1,3)*dgm(2,3)&
 385 &   +gm(2,3)*(-6*gm(1,3)*dgm(1,2)-6*gm(1,2)*dgm(1,3)+18*gm(1,1)*dgm(2,3))&
 386 &   +9*gm(1,2)**2*dgm(3,3)+gm(2,2)*(15*gm(3,3)*dgm(1,1)+18*gm(1,3)&
 387 &   *dgm(1,3)+15*gm(1,1)*dgm(3,3))
 388    cm(1,5,4)=3*gm(1,3)**2*dgm(1,2)+gm(1,1)*(12*gm(3,3)*dgm(1,2)+3*gm(2,3)&
 389 &   *dgm(1,3))+gm(1,3)*(3*gm(2,3)*dgm(1,1)+6*gm(1,2)*dgm(1,3)+3*gm(1,1)&
 390 &   *dgm(2,3))+12*gm(1,2)*(gm(3,3)*dgm(1,1)+gm(1,1)*dgm(3,3))
 391    cm(1,6,4)=3*gm(1,2)**2*dgm(1,3)+gm(1,1)*(3*gm(2,3)*dgm(1,2)+12*gm(2,2)&
 392 &   *dgm(1,3))+gm(1,3)*(12*gm(2,2)*dgm(1,1)+6*gm(1,2)*dgm(1,2)+12*gm(1,1)&
 393 &   *dgm(2,2))+3*gm(1,2)*(gm(2,3)*dgm(1,1)+gm(1,1)*dgm(2,3))
 394    cm(1,7,4)=-3*gm(2,2)**2*dgm(1,3)+9*gm(1,2)*gm(2,3)*dgm(2,2)+gm(2,2)&
 395 &   *(9*gm(2,3)*dgm(1,2)-6*gm(1,3)*dgm(2,2)+9*gm(1,2)*dgm(2,3))
 396    cm(1,8,4)=3*gm(2,3)**2*dgm(1,3)+gm(3,3)*(12*gm(1,3)*dgm(2,2)+3*gm(1,2)&
 397 &   *dgm(2,3))+gm(2,3)*(3*gm(3,3)*dgm(1,2)+6*gm(1,3)*dgm(2,3)+3*gm(1,2)&
 398 &   *dgm(3,3))+12*gm(2,2)*(gm(3,3)*dgm(1,3)+gm(1,3)*dgm(3,3))
 399    cm(1,9,4)=3*gm(2,3)**2*dgm(1,2)+12*gm(1,2)*gm(3,3)*dgm(2,2)+gm(2,3)&
 400 &   *(3*gm(2,2)*dgm(1,3)+3*gm(1,3)*dgm(2,2)+6*gm(1,2)*dgm(2,3))+gm(2,2)&
 401 &   *(12*gm(3,3)*dgm(1,2)+3*gm(1,3)*dgm(2,3)+12*gm(1,2)*dgm(3,3))
 402    cm(1,10,4)=-3*gm(3,3)**2*dgm(1,2)+9*gm(1,3)*gm(2,3)*dgm(3,3)+gm(3,3)&
 403 &   *(9*gm(2,3)*dgm(1,3)+9*gm(1,3)*dgm(2,3)-6*gm(1,2)*dgm(3,3))
 404    cm(1,1,5)=gm(1,1)*(6*gm(1,3)*dgm(1,1)+3*gm(1,1)*dgm(1,3))
 405    cm(1,2,5)=1.5d0*gm(1,2)**2*dgm(1,3)+gm(1,1)*(12*gm(2,3)*dgm(1,2)&
 406 &   -4.5d0*gm(2,2)*dgm(1,3))+gm(1,3)*(-4.5d0*gm(2,2)*dgm(1,1)+3*gm(1,2)&
 407 &   *dgm(1,2)-4.5d0*gm(1,1)*dgm(2,2))+12*gm(1,2)*(gm(2,3)*dgm(1,1)&
 408 &   +gm(1,1)*dgm(2,3))
 409    cm(1,3,5)=4.5d0*gm(1,3)**2*dgm(1,3)+7.5d0*gm(1,1)*gm(3,3)*dgm(1,3)&
 410 &   +7.5d0*gm(1,3)*(gm(3,3)*dgm(1,1)+gm(1,1)*dgm(3,3))
 411    cm(1,4,5)=3*gm(1,3)**2*dgm(1,2)+gm(1,1)*(12*gm(3,3)*dgm(1,2)+3*gm(2,3)&
 412 &   *dgm(1,3))+gm(1,3)*(3*gm(2,3)*dgm(1,1)+6*gm(1,2)*dgm(1,3)+3*gm(1,1)&
 413 &   *dgm(2,3))+12*gm(1,2)*(gm(3,3)*dgm(1,1)+gm(1,1)*dgm(3,3))
 414    cm(1,5,5)=3*gm(1,3)**2*dgm(1,1)+6*gm(1,1)*gm(1,3)*dgm(1,3)+gm(1,1)&
 415 &   *(12*gm(3,3)*dgm(1,1)+6*gm(1,1)*dgm(3,3))
 416    cm(1,6,5)=3*gm(1,2)*(gm(1,3)*dgm(1,1)+gm(1,1)*dgm(1,3))+gm(1,1)&
 417 &   *(12*gm(2,3)*dgm(1,1)+3*gm(1,3)*dgm(1,2)+6*gm(1,1)*dgm(2,3))
 418    cm(1,7,5)=-1.5d0*gm(1,1)*gm(2,3)*dgm(2,2)+gm(1,2)*(15*gm(2,3)&
 419 &   *dgm(1,2)-3*gm(1,3)*dgm(2,2))+7.5d0*gm(1,2)**2*dgm(2,3)+gm(2,2)&
 420 &   *(-1.5d0*gm(2,3)*dgm(1,1)-3*gm(1,3)*dgm(1,2)-3*gm(1,2)*dgm(1,3)&
 421 &   -1.5d0*gm(1,1)*dgm(2,3))
 422    cm(1,8,5)=1.5d0*gm(1,3)**2*dgm(2,3)+gm(3,3)*(12*gm(1,2)*dgm(1,3)&
 423 &   -4.5d0*gm(1,1)*dgm(2,3))+gm(2,3)*(-4.5d0*gm(3,3)*dgm(1,1)+3*gm(1,3)&
 424 &   *dgm(1,3)-4.5d0*gm(1,1)*dgm(3,3))+12*gm(1,3)*(gm(3,3)*dgm(1,2)&
 425 &   +gm(1,2)*dgm(3,3))
 426    cm(1,9,5)=-3*gm(2,3)**2*dgm(1,1)+15*gm(1,2)*gm(3,3)*dgm(1,2)-3*gm(1,3)&
 427 &   **2*dgm(2,2)-1.5d0*gm(1,1)*gm(3,3)*dgm(2,2)+9*gm(1,2)*gm(1,3)&
 428 &   *dgm(2,3)+gm(2,3)*(9*gm(1,3)*dgm(1,2)+9*gm(1,2)*dgm(1,3)-6*gm(1,1)&
 429 &   *dgm(2,3))+7.5d0*gm(1,2)**2*dgm(3,3)+gm(2,2)*(-1.5d0*gm(3,3)&
 430 &   *dgm(1,1)-6*gm(1,3)*dgm(1,3)-1.5d0*gm(1,1)*dgm(3,3))
 431    cm(1,10,5)=-1.5d0*gm(3,3)**2*dgm(1,1)+4.5d0*gm(1,3)**2*dgm(3,3)&
 432 &   +gm(3,3)*(9*gm(1,3)*dgm(1,3)-3*gm(1,1)*dgm(3,3))
 433    cm(1,1,6)=gm(1,1)*(6*gm(1,2)*dgm(1,1)+3*gm(1,1)*dgm(1,2))
 434    cm(1,2,6)=4.5d0*gm(1,2)**2*dgm(1,2)+7.5d0*gm(1,1)*gm(2,2)*dgm(1,2)&
 435 &   +7.5d0*gm(1,2)*(gm(2,2)*dgm(1,1)+gm(1,1)*dgm(2,2))
 436    cm(1,3,6)=1.5d0*gm(1,3)**2*dgm(1,2)+gm(1,1)*(-4.5d0*gm(3,3)*dgm(1,2)&
 437 &   +12*gm(2,3)*dgm(1,3))+gm(1,3)*(12*gm(2,3)*dgm(1,1)+3*gm(1,2)&
 438 &   *dgm(1,3)+12*gm(1,1)*dgm(2,3))-4.5d0*gm(1,2)*(gm(3,3)*dgm(1,1)&
 439 &   +gm(1,1)*dgm(3,3))
 440    cm(1,4,6)=3*gm(1,2)**2*dgm(1,3)+gm(1,1)*(3*gm(2,3)*dgm(1,2)+12*gm(2,2)&
 441 &   *dgm(1,3))+gm(1,3)*(12*gm(2,2)*dgm(1,1)+6*gm(1,2)*dgm(1,2)+12*gm(1,1)&
 442 &   *dgm(2,2))+3*gm(1,2)*(gm(2,3)*dgm(1,1)+gm(1,1)*dgm(2,3))
 443    cm(1,5,6)=3*gm(1,2)*(gm(1,3)*dgm(1,1)+gm(1,1)*dgm(1,3))+gm(1,1)&
 444 &   *(12*gm(2,3)*dgm(1,1)+3*gm(1,3)*dgm(1,2)+6*gm(1,1)*dgm(2,3))
 445    cm(1,6,6)=3*gm(1,2)**2*dgm(1,1)+6*gm(1,1)*gm(1,2)*dgm(1,2)+gm(1,1)&
 446 &   *(12*gm(2,2)*dgm(1,1)+6*gm(1,1)*dgm(2,2))
 447    cm(1,7,6)=-1.5d0*gm(2,2)**2*dgm(1,1)+4.5d0*gm(1,2)**2*dgm(2,2)&
 448 &   +gm(2,2)*(9*gm(1,2)*dgm(1,2)-3*gm(1,1)*dgm(2,2))
 449    cm(1,8,6)=-3*gm(2,3)**2*dgm(1,1)-6*gm(1,2)*gm(3,3)*dgm(1,2)+7.5d0*gm(1,3)&
 450 &   **2*dgm(2,2)-1.5d0*gm(1,1)*gm(3,3)*dgm(2,2)+9*gm(1,2)*gm(1,3)&
 451 &   *dgm(2,3)+gm(2,3)*(9*gm(1,3)*dgm(1,2)+9*gm(1,2)*dgm(1,3)-6*gm(1,1)&
 452 &   *dgm(2,3))-3*gm(1,2)**2*dgm(3,3)+gm(2,2)*(-1.5d0*gm(3,3)*dgm(1,1)&
 453 &   +15*gm(1,3)*dgm(1,3)-1.5d0*gm(1,1)*dgm(3,3))
 454    cm(1,9,6)=-4.5d0*gm(1,1)*gm(2,3)*dgm(2,2)+gm(1,2)*(3*gm(2,3)*dgm(1,2)&
 455 &   +12*gm(1,3)*dgm(2,2))+1.5d0*gm(1,2)**2*dgm(2,3)+gm(2,2)*(-4.5d0*gm(2,3)&
 456 &   *dgm(1,1)+12*gm(1,3)*dgm(1,2)+12*gm(1,2)*dgm(1,3)-4.5d0*gm(1,1)&
 457 &   *dgm(2,3))
 458    cm(1,10,6)=7.5d0*gm(1,3)**2*dgm(2,3)+gm(3,3)*(-3*gm(1,2)*dgm(1,3)&
 459 &   -1.5d0*gm(1,1)*dgm(2,3))+gm(2,3)*(-1.5d0*gm(3,3)*dgm(1,1)+15*gm(1,3)&
 460 &   *dgm(1,3)-1.5d0*gm(1,1)*dgm(3,3))-3*gm(1,3)*(gm(3,3)*dgm(1,2)&
 461 &   +gm(1,2)*dgm(3,3))
 462    cm(1,1,7)=7.5d0*gm(1,2)**2*dgm(1,2)-1.5d0*gm(1,1)*gm(2,2)*dgm(1,2)&
 463 &   -1.5d0*gm(1,2)*(gm(2,2)*dgm(1,1)+gm(1,1)*dgm(2,2))
 464    cm(1,2,7)=gm(2,2)*(3*gm(2,2)*dgm(1,2)+6*gm(1,2)*dgm(2,2))
 465    cm(1,3,7)=7.5d0*gm(2,3)**2*dgm(1,2)-1.5d0*gm(1,2)*gm(3,3)*dgm(2,2)&
 466 &   +gm(2,3)*(-3*gm(2,2)*dgm(1,3)-3*gm(1,3)*dgm(2,2)+15*gm(1,2)*dgm(2,3))&
 467 &   +gm(2,2)*(-1.5d0*gm(3,3)*dgm(1,2)-3*gm(1,3)*dgm(2,3)-1.5d0*gm(1,2)&
 468 &   *dgm(3,3))
 469    cm(1,4,7)=-3*gm(2,2)**2*dgm(1,3)+9*gm(1,2)*gm(2,3)*dgm(2,2)+gm(2,2)&
 470 &   *(9*gm(2,3)*dgm(1,2)-6*gm(1,3)*dgm(2,2)+9*gm(1,2)*dgm(2,3))
 471    cm(1,5,7)=-1.5d0*gm(1,1)*gm(2,3)*dgm(2,2)+gm(1,2)*(15*gm(2,3)&
 472 &   *dgm(1,2)-3*gm(1,3)*dgm(2,2))+7.5d0*gm(1,2)**2*dgm(2,3)+gm(2,2)&
 473 &   *(-1.5d0*gm(2,3)*dgm(1,1)-3*gm(1,3)*dgm(1,2)-3*gm(1,2)*dgm(1,3)&
 474 &   -1.5d0*gm(1,1)*dgm(2,3))
 475    cm(1,6,7)=-1.5d0*gm(2,2)**2*dgm(1,1)+4.5d0*gm(1,2)**2*dgm(2,2)&
 476 &   +gm(2,2)*(9*gm(1,2)*dgm(1,2)-3*gm(1,1)*dgm(2,2))
 477    cm(1,7,7)=3*gm(2,2)**2*dgm(2,2)
 478    cm(1,8,7)=4.5d0*gm(2,3)**2*dgm(2,2)+9*gm(2,2)*gm(2,3)*dgm(2,3)&
 479 &   +gm(2,2)*(-3*gm(3,3)*dgm(2,2)-1.5d0*gm(2,2)*dgm(3,3))
 480    cm(1,9,7)=gm(2,2)*(6*gm(2,3)*dgm(2,2)+3*gm(2,2)*dgm(2,3))
 481    cm(1,10,7)=7.5d0*gm(2,3)**2*dgm(2,3)-1.5d0*gm(2,2)*gm(3,3)*dgm(2,3)&
 482 &   -1.5d0*gm(2,3)*(gm(3,3)*dgm(2,2)+gm(2,2)*dgm(3,3))
 483    cm(1,1,8)=7.5d0*gm(1,3)**2*dgm(1,2)+gm(1,1)*(-1.5d0*gm(3,3)*dgm(1,2)&
 484 &   -3*gm(2,3)*dgm(1,3))+gm(1,3)*(-3*gm(2,3)*dgm(1,1)+15*gm(1,2)&
 485 &   *dgm(1,3)-3*gm(1,1)*dgm(2,3))-1.5d0*gm(1,2)*(gm(3,3)*dgm(1,1)&
 486 &   +gm(1,1)*dgm(3,3))
 487    cm(1,2,8)=1.5d0*gm(2,3)**2*dgm(1,2)-4.5d0*gm(1,2)*gm(3,3)*dgm(2,2)&
 488 &   +gm(2,3)*(12*gm(2,2)*dgm(1,3)+12*gm(1,3)*dgm(2,2)+3*gm(1,2)*dgm(2,3))&
 489 &   +gm(2,2)*(-4.5d0*gm(3,3)*dgm(1,2)+12*gm(1,3)*dgm(2,3)-4.5d0*gm(1,2)&
 490 &   *dgm(3,3))
 491    cm(1,3,8)=6*gm(3,3)**2*dgm(1,2)+3*gm(1,3)*gm(2,3)*dgm(3,3)+gm(3,3)&
 492 &   *(3*gm(2,3)*dgm(1,3)+3*gm(1,3)*dgm(2,3)+12*gm(1,2)*dgm(3,3))
 493    cm(1,4,8)=3*gm(2,3)**2*dgm(1,3)+gm(3,3)*(12*gm(1,3)*dgm(2,2)+3*gm(1,2)&
 494 &   *dgm(2,3))+gm(2,3)*(3*gm(3,3)*dgm(1,2)+6*gm(1,3)*dgm(2,3)+3*gm(1,2)&
 495 &   *dgm(3,3))+12*gm(2,2)*(gm(3,3)*dgm(1,3)+gm(1,3)*dgm(3,3))
 496    cm(1,5,8)=1.5d0*gm(1,3)**2*dgm(2,3)+gm(3,3)*(12*gm(1,2)*dgm(1,3)&
 497 &   -4.5d0*gm(1,1)*dgm(2,3))+gm(2,3)*(-4.5d0*gm(3,3)*dgm(1,1)+3*gm(1,3)&
 498 &   *dgm(1,3)-4.5d0*gm(1,1)*dgm(3,3))+12*gm(1,3)*(gm(3,3)*dgm(1,2)&
 499 &   +gm(1,2)*dgm(3,3))
 500    cm(1,6,8)=-3*gm(2,3)**2*dgm(1,1)-6*gm(1,2)*gm(3,3)*dgm(1,2)+7.5d0*gm(1,3)&
 501 &   **2*dgm(2,2)-1.5d0*gm(1,1)*gm(3,3)*dgm(2,2)+9*gm(1,2)*gm(1,3)&
 502 &   *dgm(2,3)+gm(2,3)*(9*gm(1,3)*dgm(1,2)+9*gm(1,2)*dgm(1,3)-6*gm(1,1)&
 503 &   *dgm(2,3))-3*gm(1,2)**2*dgm(3,3)+gm(2,2)*(-1.5d0*gm(3,3)*dgm(1,1)&
 504 &   +15*gm(1,3)*dgm(1,3)-1.5d0*gm(1,1)*dgm(3,3))
 505    cm(1,7,8)=4.5d0*gm(2,3)**2*dgm(2,2)+9*gm(2,2)*gm(2,3)*dgm(2,3)&
 506 &   +gm(2,2)*(-3*gm(3,3)*dgm(2,2)-1.5d0*gm(2,2)*dgm(3,3))
 507    cm(1,8,8)=6*gm(3,3)**2*dgm(2,2)+3*gm(2,3)**2*dgm(3,3)+gm(3,3)&
 508 &   *(6*gm(2,3)*dgm(2,3)+12*gm(2,2)*dgm(3,3))
 509    cm(1,9,8)=4.5d0*gm(2,3)**2*dgm(2,3)+7.5d0*gm(2,2)*gm(3,3)*dgm(2,3)&
 510 &   +7.5d0*gm(2,3)*(gm(3,3)*dgm(2,2)+gm(2,2)*dgm(3,3))
 511    cm(1,10,8)=gm(3,3)*(3*gm(3,3)*dgm(2,3)+6*gm(2,3)*dgm(3,3))
 512    cm(1,1,9)=7.5d0*gm(1,2)**2*dgm(1,3)+gm(1,1)*(-3*gm(2,3)*dgm(1,2)&
 513 &   -1.5d0*gm(2,2)*dgm(1,3))+gm(1,3)*(-1.5d0*gm(2,2)*dgm(1,1)+15*gm(1,2)&
 514 &   *dgm(1,2)-1.5d0*gm(1,1)*dgm(2,2))-3*gm(1,2)*(gm(2,3)*dgm(1,1)&
 515 &   +gm(1,1)*dgm(2,3))
 516    cm(1,2,9)=6*gm(2,2)**2*dgm(1,3)+3*gm(1,2)*gm(2,3)*dgm(2,2)+gm(2,2)&
 517 &   *(3*gm(2,3)*dgm(1,2)+12*gm(1,3)*dgm(2,2)+3*gm(1,2)*dgm(2,3))
 518    cm(1,3,9)=1.5d0*gm(2,3)**2*dgm(1,3)+gm(3,3)*(-4.5d0*gm(1,3)*dgm(2,2)&
 519 &   +12*gm(1,2)*dgm(2,3))+gm(2,3)*(12*gm(3,3)*dgm(1,2)+3*gm(1,3)&
 520 &   *dgm(2,3)+12*gm(1,2)*dgm(3,3))-4.5d0*gm(2,2)*(gm(3,3)*dgm(1,3)&
 521 &   +gm(1,3)*dgm(3,3))
 522    cm(1,4,9)=3*gm(2,3)**2*dgm(1,2)+12*gm(1,2)*gm(3,3)*dgm(2,2)+gm(2,3)&
 523 &   *(3*gm(2,2)*dgm(1,3)+3*gm(1,3)*dgm(2,2)+6*gm(1,2)*dgm(2,3))+gm(2,2)&
 524 &   *(12*gm(3,3)*dgm(1,2)+3*gm(1,3)*dgm(2,3)+12*gm(1,2)*dgm(3,3))
 525    cm(1,5,9)=-3*gm(2,3)**2*dgm(1,1)+15*gm(1,2)*gm(3,3)*dgm(1,2)-3*gm(1,3)&
 526 &   **2*dgm(2,2)-1.5d0*gm(1,1)*gm(3,3)*dgm(2,2)+9*gm(1,2)*gm(1,3)&
 527 &   *dgm(2,3)+gm(2,3)*(9*gm(1,3)*dgm(1,2)+9*gm(1,2)*dgm(1,3)-6*gm(1,1)&
 528 &   *dgm(2,3))+7.5d0*gm(1,2)**2*dgm(3,3)+gm(2,2)*(-1.5d0*gm(3,3)&
 529 &   *dgm(1,1)-6*gm(1,3)*dgm(1,3)-1.5d0*gm(1,1)*dgm(3,3))
 530    cm(1,6,9)=-4.5d0*gm(1,1)*gm(2,3)*dgm(2,2)+gm(1,2)*(3*gm(2,3)*dgm(1,2)&
 531 &   +12*gm(1,3)*dgm(2,2))+1.5d0*gm(1,2)**2*dgm(2,3)+gm(2,2)*(-4.5d0*gm(2,3)&
 532 &   *dgm(1,1)+12*gm(1,3)*dgm(1,2)+12*gm(1,2)*dgm(1,3)-4.5d0*gm(1,1)&
 533 &   *dgm(2,3))
 534    cm(1,7,9)=gm(2,2)*(6*gm(2,3)*dgm(2,2)+3*gm(2,2)*dgm(2,3))
 535    cm(1,8,9)=4.5d0*gm(2,3)**2*dgm(2,3)+7.5d0*gm(2,2)*gm(3,3)*dgm(2,3)&
 536 &   +7.5d0*gm(2,3)*(gm(3,3)*dgm(2,2)+gm(2,2)*dgm(3,3))
 537    cm(1,9,9)=3*gm(2,3)**2*dgm(2,2)+6*gm(2,2)*gm(2,3)*dgm(2,3)+gm(2,2)&
 538 &   *(12*gm(3,3)*dgm(2,2)+6*gm(2,2)*dgm(3,3))
 539    cm(1,10,9)=-1.5d0*gm(3,3)**2*dgm(2,2)+4.5d0*gm(2,3)**2*dgm(3,3)&
 540 &   +gm(3,3)*(9*gm(2,3)*dgm(2,3)-3*gm(2,2)*dgm(3,3))
 541    cm(1,1,10)=7.5d0*gm(1,3)**2*dgm(1,3)-1.5d0*gm(1,1)*gm(3,3)*dgm(1,3)&
 542 &   -1.5d0*gm(1,3)*(gm(3,3)*dgm(1,1)+gm(1,1)*dgm(3,3))
 543    cm(1,2,10)=7.5d0*gm(2,3)**2*dgm(1,3)+gm(3,3)*(-1.5d0*gm(1,3)*dgm(2,2)&
 544 &   -3*gm(1,2)*dgm(2,3))+gm(2,3)*(-3*gm(3,3)*dgm(1,2)+15*gm(1,3)&
 545 &   *dgm(2,3)-3*gm(1,2)*dgm(3,3))-1.5d0*gm(2,2)*(gm(3,3)*dgm(1,3)&
 546 &   +gm(1,3)*dgm(3,3))
 547    cm(1,3,10)=gm(3,3)*(3*gm(3,3)*dgm(1,3)+6*gm(1,3)*dgm(3,3))
 548    cm(1,4,10)=-3*gm(3,3)**2*dgm(1,2)+9*gm(1,3)*gm(2,3)*dgm(3,3)+gm(3,3)&
 549 &   *(9*gm(2,3)*dgm(1,3)+9*gm(1,3)*dgm(2,3)-6*gm(1,2)*dgm(3,3))
 550    cm(1,5,10)=-1.5d0*gm(3,3)**2*dgm(1,1)+4.5d0*gm(1,3)**2*dgm(3,3)&
 551 &   +gm(3,3)*(9*gm(1,3)*dgm(1,3)-3*gm(1,1)*dgm(3,3))
 552    cm(1,6,10)=7.5d0*gm(1,3)**2*dgm(2,3)+gm(3,3)*(-3*gm(1,2)*dgm(1,3)&
 553 &   -1.5d0*gm(1,1)*dgm(2,3))+gm(2,3)*(-1.5d0*gm(3,3)*dgm(1,1)+15*gm(1,3)&
 554 &   *dgm(1,3)-1.5d0*gm(1,1)*dgm(3,3))-3*gm(1,3)*(gm(3,3)*dgm(1,2)&
 555 &   +gm(1,2)*dgm(3,3))
 556    cm(1,7,10)=7.5d0*gm(2,3)**2*dgm(2,3)-1.5d0*gm(2,2)*gm(3,3)*dgm(2,3)&
 557 &   -1.5d0*gm(2,3)*(gm(3,3)*dgm(2,2)+gm(2,2)*dgm(3,3))
 558    cm(1,8,10)=gm(3,3)*(3*gm(3,3)*dgm(2,3)+6*gm(2,3)*dgm(3,3))
 559    cm(1,9,10)=-1.5d0*gm(3,3)**2*dgm(2,2)+4.5d0*gm(2,3)**2*dgm(3,3)&
 560 &   +gm(3,3)*(9*gm(2,3)*dgm(2,3)-3*gm(2,2)*dgm(3,3))
 561    cm(1,10,10)=3*gm(3,3)**2*dgm(3,3)
 562    cm(2,1,2)=gm(1,1)*(6*gm(1,2)*dgm(1,1)+3*gm(1,1)*dgm(1,2))
 563    cm(2,2,2)=4.5d0*gm(1,2)**2*dgm(1,2)+7.5d0*gm(1,1)*gm(2,2)*dgm(1,2)&
 564 &   +7.5d0*gm(1,2)*(gm(2,2)*dgm(1,1)+gm(1,1)*dgm(2,2))
 565    cm(2,3,2)=1.5d0*gm(1,3)**2*dgm(1,2)+gm(1,1)*(-4.5d0*gm(3,3)*dgm(1,2)&
 566 &   +12*gm(2,3)*dgm(1,3))+gm(1,3)*(12*gm(2,3)*dgm(1,1)+3*gm(1,2)&
 567 &   *dgm(1,3)+12*gm(1,1)*dgm(2,3))-4.5d0*gm(1,2)*(gm(3,3)*dgm(1,1)&
 568 &   +gm(1,1)*dgm(3,3))
 569    cm(2,4,2)=3*gm(1,2)**2*dgm(1,3)+gm(1,1)*(3*gm(2,3)*dgm(1,2)+12*gm(2,2)&
 570 &   *dgm(1,3))+gm(1,3)*(12*gm(2,2)*dgm(1,1)+6*gm(1,2)*dgm(1,2)+12*gm(1,1)&
 571 &   *dgm(2,2))+3*gm(1,2)*(gm(2,3)*dgm(1,1)+gm(1,1)*dgm(2,3))
 572    cm(2,5,2)=3*gm(1,2)*(gm(1,3)*dgm(1,1)+gm(1,1)*dgm(1,3))+gm(1,1)&
 573 &   *(12*gm(2,3)*dgm(1,1)+3*gm(1,3)*dgm(1,2)+6*gm(1,1)*dgm(2,3))
 574    cm(2,6,2)=3*gm(1,2)**2*dgm(1,1)+6*gm(1,1)*gm(1,2)*dgm(1,2)+gm(1,1)&
 575 &   *(12*gm(2,2)*dgm(1,1)+6*gm(1,1)*dgm(2,2))
 576    cm(2,7,2)=-1.5d0*gm(2,2)**2*dgm(1,1)+4.5d0*gm(1,2)**2*dgm(2,2)&
 577 &   +gm(2,2)*(9*gm(1,2)*dgm(1,2)-3*gm(1,1)*dgm(2,2))
 578    cm(2,8,2)=-3*gm(2,3)**2*dgm(1,1)-6*gm(1,2)*gm(3,3)*dgm(1,2)+7.5d0*gm(1,3)&
 579 &   **2*dgm(2,2)-1.5d0*gm(1,1)*gm(3,3)*dgm(2,2)+9*gm(1,2)*gm(1,3)&
 580 &   *dgm(2,3)+gm(2,3)*(9*gm(1,3)*dgm(1,2)+9*gm(1,2)*dgm(1,3)-6*gm(1,1)&
 581 &   *dgm(2,3))-3*gm(1,2)**2*dgm(3,3)+gm(2,2)*(-1.5d0*gm(3,3)*dgm(1,1)&
 582 &   +15*gm(1,3)*dgm(1,3)-1.5d0*gm(1,1)*dgm(3,3))
 583    cm(2,9,2)=-4.5d0*gm(1,1)*gm(2,3)*dgm(2,2)+gm(1,2)*(3*gm(2,3)*dgm(1,2)&
 584 &   +12*gm(1,3)*dgm(2,2))+1.5d0*gm(1,2)**2*dgm(2,3)+gm(2,2)*(-4.5d0*gm(2,3)&
 585 &   *dgm(1,1)+12*gm(1,3)*dgm(1,2)+12*gm(1,2)*dgm(1,3)-4.5d0*gm(1,1)&
 586 &   *dgm(2,3))
 587    cm(2,10,2)=7.5d0*gm(1,3)**2*dgm(2,3)+gm(3,3)*(-3*gm(1,2)*dgm(1,3)&
 588 &   -1.5d0*gm(1,1)*dgm(2,3))+gm(2,3)*(-1.5d0*gm(3,3)*dgm(1,1)+15*gm(1,3)&
 589 &   *dgm(1,3)-1.5d0*gm(1,1)*dgm(3,3))-3*gm(1,3)*(gm(3,3)*dgm(1,2)&
 590 &   +gm(1,2)*dgm(3,3))
 591    cm(2,1,4)=gm(1,1)*(6*gm(1,3)*dgm(1,1)+3*gm(1,1)*dgm(1,3))
 592    cm(2,2,4)=1.5d0*gm(1,2)**2*dgm(1,3)+gm(1,1)*(12*gm(2,3)*dgm(1,2)&
 593 &   -4.5d0*gm(2,2)*dgm(1,3))+gm(1,3)*(-4.5d0*gm(2,2)*dgm(1,1)+3*gm(1,2)&
 594 &   *dgm(1,2)-4.5d0*gm(1,1)*dgm(2,2))+12*gm(1,2)*(gm(2,3)*dgm(1,1)&
 595 &   +gm(1,1)*dgm(2,3))
 596    cm(2,3,4)=4.5d0*gm(1,3)**2*dgm(1,3)+7.5d0*gm(1,1)*gm(3,3)*dgm(1,3)&
 597 &   +7.5d0*gm(1,3)*(gm(3,3)*dgm(1,1)+gm(1,1)*dgm(3,3))
 598    cm(2,4,4)=3*gm(1,3)**2*dgm(1,2)+gm(1,1)*(12*gm(3,3)*dgm(1,2)+3*gm(2,3)&
 599 &   *dgm(1,3))+gm(1,3)*(3*gm(2,3)*dgm(1,1)+6*gm(1,2)*dgm(1,3)+3*gm(1,1)&
 600 &   *dgm(2,3))+12*gm(1,2)*(gm(3,3)*dgm(1,1)+gm(1,1)*dgm(3,3))
 601    cm(2,5,4)=3*gm(1,3)**2*dgm(1,1)+6*gm(1,1)*gm(1,3)*dgm(1,3)+gm(1,1)&
 602 &   *(12*gm(3,3)*dgm(1,1)+6*gm(1,1)*dgm(3,3))
 603    cm(2,6,4)=3*gm(1,2)*(gm(1,3)*dgm(1,1)+gm(1,1)*dgm(1,3))+gm(1,1)&
 604 &   *(12*gm(2,3)*dgm(1,1)+3*gm(1,3)*dgm(1,2)+6*gm(1,1)*dgm(2,3))
 605    cm(2,7,4)=-1.5d0*gm(1,1)*gm(2,3)*dgm(2,2)+gm(1,2)*(15*gm(2,3)&
 606 &   *dgm(1,2)-3*gm(1,3)*dgm(2,2))+7.5d0*gm(1,2)**2*dgm(2,3)+gm(2,2)&
 607 &   *(-1.5d0*gm(2,3)*dgm(1,1)-3*gm(1,3)*dgm(1,2)-3*gm(1,2)*dgm(1,3)&
 608 &   -1.5d0*gm(1,1)*dgm(2,3))
 609    cm(2,8,4)=1.5d0*gm(1,3)**2*dgm(2,3)+gm(3,3)*(12*gm(1,2)*dgm(1,3)&
 610 &   -4.5d0*gm(1,1)*dgm(2,3))+gm(2,3)*(-4.5d0*gm(3,3)*dgm(1,1)+3*gm(1,3)&
 611 &   *dgm(1,3)-4.5d0*gm(1,1)*dgm(3,3))+12*gm(1,3)*(gm(3,3)*dgm(1,2)&
 612 &   +gm(1,2)*dgm(3,3))
 613    cm(2,9,4)=-3*gm(2,3)**2*dgm(1,1)+15*gm(1,2)*gm(3,3)*dgm(1,2)-3*gm(1,3)&
 614 &   **2*dgm(2,2)-1.5d0*gm(1,1)*gm(3,3)*dgm(2,2)+9*gm(1,2)*gm(1,3)&
 615 &   *dgm(2,3)+gm(2,3)*(9*gm(1,3)*dgm(1,2)+9*gm(1,2)*dgm(1,3)-6*gm(1,1)&
 616 &   *dgm(2,3))+7.5d0*gm(1,2)**2*dgm(3,3)+gm(2,2)*(-1.5d0*gm(3,3)&
 617 &   *dgm(1,1)-6*gm(1,3)*dgm(1,3)-1.5d0*gm(1,1)*dgm(3,3))
 618    cm(2,10,4)=-1.5d0*gm(3,3)**2*dgm(1,1)+4.5d0*gm(1,3)**2*dgm(3,3)&
 619 &   +gm(3,3)*(9*gm(1,3)*dgm(1,3)-3*gm(1,1)*dgm(3,3))
 620    cm(2,1,6)=3*gm(1,1)**2*dgm(1,1)
 621    cm(2,2,6)=4.5d0*gm(1,2)**2*dgm(1,1)+9*gm(1,1)*gm(1,2)*dgm(1,2)&
 622 &   +gm(1,1)*(-3*gm(2,2)*dgm(1,1)-1.5d0*gm(1,1)*dgm(2,2))
 623    cm(2,3,6)=4.5d0*gm(1,3)**2*dgm(1,1)+9*gm(1,1)*gm(1,3)*dgm(1,3)&
 624 &   +gm(1,1)*(-3*gm(3,3)*dgm(1,1)-1.5d0*gm(1,1)*dgm(3,3))
 625    cm(2,4,6)=9*gm(1,2)*(gm(1,3)*dgm(1,1)+gm(1,1)*dgm(1,3))+gm(1,1)&
 626 &   *(-6*gm(2,3)*dgm(1,1)+9*gm(1,3)*dgm(1,2)-3*gm(1,1)*dgm(2,3))
 627    cm(2,5,6)=gm(1,1)*(6*gm(1,3)*dgm(1,1)+3*gm(1,1)*dgm(1,3))
 628    cm(2,6,6)=gm(1,1)*(6*gm(1,2)*dgm(1,1)+3*gm(1,1)*dgm(1,2))
 629    cm(2,7,6)=7.5d0*gm(1,2)**2*dgm(1,2)-1.5d0*gm(1,1)*gm(2,2)*dgm(1,2)&
 630 &   -1.5d0*gm(1,2)*(gm(2,2)*dgm(1,1)+gm(1,1)*dgm(2,2))
 631    cm(2,8,6)=7.5d0*gm(1,3)**2*dgm(1,2)+gm(1,1)*(-1.5d0*gm(3,3)*dgm(1,2)&
 632 &   -3*gm(2,3)*dgm(1,3))+gm(1,3)*(-3*gm(2,3)*dgm(1,1)+15*gm(1,2)&
 633 &   *dgm(1,3)-3*gm(1,1)*dgm(2,3))-1.5d0*gm(1,2)*(gm(3,3)*dgm(1,1)&
 634 &   +gm(1,1)*dgm(3,3))
 635    cm(2,9,6)=7.5d0*gm(1,2)**2*dgm(1,3)+gm(1,1)*(-3*gm(2,3)*dgm(1,2)&
 636 &   -1.5d0*gm(2,2)*dgm(1,3))+gm(1,3)*(-1.5d0*gm(2,2)*dgm(1,1)+15*gm(1,2)&
 637 &   *dgm(1,2)-1.5d0*gm(1,1)*dgm(2,2))-3*gm(1,2)*(gm(2,3)*dgm(1,1)&
 638 &   +gm(1,1)*dgm(2,3))
 639    cm(2,10,6)=7.5d0*gm(1,3)**2*dgm(1,3)-1.5d0*gm(1,1)*gm(3,3)*dgm(1,3)&
 640 &   -1.5d0*gm(1,3)*(gm(3,3)*dgm(1,1)+gm(1,1)*dgm(3,3))
 641    cm(2,1,7)=4.5d0*gm(1,2)**2*dgm(1,1)+9*gm(1,1)*gm(1,2)*dgm(1,2)&
 642 &   +gm(1,1)*(-3*gm(2,2)*dgm(1,1)-1.5d0*gm(1,1)*dgm(2,2))
 643    cm(2,2,7)=6*gm(2,2)**2*dgm(1,1)+3*gm(1,2)**2*dgm(2,2)+gm(2,2)&
 644 &   *(6*gm(1,2)*dgm(1,2)+12*gm(1,1)*dgm(2,2))
 645    cm(2,3,7)=7.5d0*gm(2,3)**2*dgm(1,1)-6*gm(1,2)*gm(3,3)*dgm(1,2)&
 646 &   -3*gm(1,3)**2*dgm(2,2)-1.5d0*gm(1,1)*gm(3,3)*dgm(2,2)+9*gm(1,2)&
 647 &   *gm(1,3)*dgm(2,3)+gm(2,3)*(9*gm(1,3)*dgm(1,2)+9*gm(1,2)*dgm(1,3)&
 648 &   +15*gm(1,1)*dgm(2,3))-3*gm(1,2)**2*dgm(3,3)+gm(2,2)*(-1.5d0*gm(3,3)&
 649 &   *dgm(1,1)-6*gm(1,3)*dgm(1,3)-1.5d0*gm(1,1)*dgm(3,3))
 650    cm(2,4,7)=12*gm(1,1)*gm(2,3)*dgm(2,2)+gm(1,2)*(6*gm(2,3)*dgm(1,2)&
 651 &   +3*gm(1,3)*dgm(2,2))+3*gm(1,2)**2*dgm(2,3)+gm(2,2)*(12*gm(2,3)&
 652 &   *dgm(1,1)+3*gm(1,3)*dgm(1,2)+3*gm(1,2)*dgm(1,3)+12*gm(1,1)*dgm(2,3))
 653    cm(2,5,7)=1.5d0*gm(1,2)**2*dgm(1,3)+gm(1,1)*(12*gm(2,3)*dgm(1,2)&
 654 &   -4.5d0*gm(2,2)*dgm(1,3))+gm(1,3)*(-4.5d0*gm(2,2)*dgm(1,1)+3*gm(1,2)&
 655 &   *dgm(1,2)-4.5d0*gm(1,1)*dgm(2,2))+12*gm(1,2)*(gm(2,3)*dgm(1,1)&
 656 &   +gm(1,1)*dgm(2,3))
 657    cm(2,6,7)=4.5d0*gm(1,2)**2*dgm(1,2)+7.5d0*gm(1,1)*gm(2,2)*dgm(1,2)&
 658 &   +7.5d0*gm(1,2)*(gm(2,2)*dgm(1,1)+gm(1,1)*dgm(2,2))
 659    cm(2,7,7)=gm(2,2)*(3*gm(2,2)*dgm(1,2)+6*gm(1,2)*dgm(2,2))
 660    cm(2,8,7)=1.5d0*gm(2,3)**2*dgm(1,2)-4.5d0*gm(1,2)*gm(3,3)*dgm(2,2)&
 661 &   +gm(2,3)*(12*gm(2,2)*dgm(1,3)+12*gm(1,3)*dgm(2,2)+3*gm(1,2)*dgm(2,3))&
 662 &   +gm(2,2)*(-4.5d0*gm(3,3)*dgm(1,2)+12*gm(1,3)*dgm(2,3)-4.5d0*gm(1,2)&
 663 &   *dgm(3,3))
 664    cm(2,9,7)=6*gm(2,2)**2*dgm(1,3)+3*gm(1,2)*gm(2,3)*dgm(2,2)+gm(2,2)&
 665 &   *(3*gm(2,3)*dgm(1,2)+12*gm(1,3)*dgm(2,2)+3*gm(1,2)*dgm(2,3))
 666    cm(2,10,7)=7.5d0*gm(2,3)**2*dgm(1,3)+gm(3,3)*(-1.5d0*gm(1,3)*dgm(2,2)&
 667 &   -3*gm(1,2)*dgm(2,3))+gm(2,3)*(-3*gm(3,3)*dgm(1,2)+15*gm(1,3)&
 668 &   *dgm(2,3)-3*gm(1,2)*dgm(3,3))-1.5d0*gm(2,2)*(gm(3,3)*dgm(1,3)&
 669 &   +gm(1,3)*dgm(3,3))
 670    cm(2,1,8)=4.5d0*gm(1,3)**2*dgm(1,1)+9*gm(1,1)*gm(1,3)*dgm(1,3)&
 671 &   +gm(1,1)*(-3*gm(3,3)*dgm(1,1)-1.5d0*gm(1,1)*dgm(3,3))
 672    cm(2,2,8)=7.5d0*gm(2,3)**2*dgm(1,1)-6*gm(1,2)*gm(3,3)*dgm(1,2)&
 673 &   -3*gm(1,3)**2*dgm(2,2)-1.5d0*gm(1,1)*gm(3,3)*dgm(2,2)+9*gm(1,2)&
 674 &   *gm(1,3)*dgm(2,3)+gm(2,3)*(9*gm(1,3)*dgm(1,2)+9*gm(1,2)*dgm(1,3)&
 675 &   +15*gm(1,1)*dgm(2,3))-3*gm(1,2)**2*dgm(3,3)+gm(2,2)*(-1.5d0*gm(3,3)&
 676 &   *dgm(1,1)-6*gm(1,3)*dgm(1,3)-1.5d0*gm(1,1)*dgm(3,3))
 677    cm(2,3,8)=6*gm(3,3)**2*dgm(1,1)+3*gm(1,3)**2*dgm(3,3)+gm(3,3)&
 678 &   *(6*gm(1,3)*dgm(1,3)+12*gm(1,1)*dgm(3,3))
 679    cm(2,4,8)=3*gm(1,3)**2*dgm(2,3)+gm(3,3)*(3*gm(1,2)*dgm(1,3)+12*gm(1,1)&
 680 &   *dgm(2,3))+gm(2,3)*(12*gm(3,3)*dgm(1,1)+6*gm(1,3)*dgm(1,3)+12*gm(1,1)&
 681 &   *dgm(3,3))+3*gm(1,3)*(gm(3,3)*dgm(1,2)+gm(1,2)*dgm(3,3))
 682    cm(2,5,8)=4.5d0*gm(1,3)**2*dgm(1,3)+7.5d0*gm(1,1)*gm(3,3)*dgm(1,3)&
 683 &   +7.5d0*gm(1,3)*(gm(3,3)*dgm(1,1)+gm(1,1)*dgm(3,3))
 684    cm(2,6,8)=1.5d0*gm(1,3)**2*dgm(1,2)+gm(1,1)*(-4.5d0*gm(3,3)*dgm(1,2)&
 685 &   +12*gm(2,3)*dgm(1,3))+gm(1,3)*(12*gm(2,3)*dgm(1,1)+3*gm(1,2)&
 686 &   *dgm(1,3)+12*gm(1,1)*dgm(2,3))-4.5d0*gm(1,2)*(gm(3,3)*dgm(1,1)&
 687 &   +gm(1,1)*dgm(3,3))
 688    cm(2,7,8)=7.5d0*gm(2,3)**2*dgm(1,2)-1.5d0*gm(1,2)*gm(3,3)*dgm(2,2)&
 689 &   +gm(2,3)*(-3*gm(2,2)*dgm(1,3)-3*gm(1,3)*dgm(2,2)+15*gm(1,2)*dgm(2,3))&
 690 &   +gm(2,2)*(-1.5d0*gm(3,3)*dgm(1,2)-3*gm(1,3)*dgm(2,3)-1.5d0*gm(1,2)&
 691 &   *dgm(3,3))
 692    cm(2,8,8)=6*gm(3,3)**2*dgm(1,2)+3*gm(1,3)*gm(2,3)*dgm(3,3)+gm(3,3)&
 693 &   *(3*gm(2,3)*dgm(1,3)+3*gm(1,3)*dgm(2,3)+12*gm(1,2)*dgm(3,3))
 694    cm(2,9,8)=1.5d0*gm(2,3)**2*dgm(1,3)+gm(3,3)*(-4.5d0*gm(1,3)*dgm(2,2)&
 695 &   +12*gm(1,2)*dgm(2,3))+gm(2,3)*(12*gm(3,3)*dgm(1,2)+3*gm(1,3)&
 696 &   *dgm(2,3)+12*gm(1,2)*dgm(3,3))-4.5d0*gm(2,2)*(gm(3,3)*dgm(1,3)&
 697 &   +gm(1,3)*dgm(3,3))
 698    cm(2,10,8)=gm(3,3)*(3*gm(3,3)*dgm(1,3)+6*gm(1,3)*dgm(3,3))
 699    cm(2,1,9)=9*gm(1,2)*(gm(1,3)*dgm(1,1)+gm(1,1)*dgm(1,3))+gm(1,1)&
 700 &   *(-6*gm(2,3)*dgm(1,1)+9*gm(1,3)*dgm(1,2)-3*gm(1,1)*dgm(2,3))
 701    cm(2,2,9)=12*gm(1,1)*gm(2,3)*dgm(2,2)+gm(1,2)*(6*gm(2,3)*dgm(1,2)&
 702 &   +3*gm(1,3)*dgm(2,2))+3*gm(1,2)**2*dgm(2,3)+gm(2,2)*(12*gm(2,3)&
 703 &   *dgm(1,1)+3*gm(1,3)*dgm(1,2)+3*gm(1,2)*dgm(1,3)+12*gm(1,1)*dgm(2,3))
 704    cm(2,3,9)=3*gm(1,3)**2*dgm(2,3)+gm(3,3)*(3*gm(1,2)*dgm(1,3)+12*gm(1,1)&
 705 &   *dgm(2,3))+gm(2,3)*(12*gm(3,3)*dgm(1,1)+6*gm(1,3)*dgm(1,3)+12*gm(1,1)&
 706 &   *dgm(3,3))+3*gm(1,3)*(gm(3,3)*dgm(1,2)+gm(1,2)*dgm(3,3))
 707    cm(2,4,9)=9*gm(2,3)**2*dgm(1,1)+18*gm(1,2)*gm(3,3)*dgm(1,2)+9*gm(1,3)&
 708 &   **2*dgm(2,2)+15*gm(1,1)*gm(3,3)*dgm(2,2)-6*gm(1,2)*gm(1,3)*dgm(2,3)&
 709 &   +gm(2,3)*(-6*gm(1,3)*dgm(1,2)-6*gm(1,2)*dgm(1,3)+18*gm(1,1)*dgm(2,3))&
 710 &   +9*gm(1,2)**2*dgm(3,3)+gm(2,2)*(15*gm(3,3)*dgm(1,1)+18*gm(1,3)&
 711 &   *dgm(1,3)+15*gm(1,1)*dgm(3,3))
 712    cm(2,5,9)=3*gm(1,3)**2*dgm(1,2)+gm(1,1)*(12*gm(3,3)*dgm(1,2)+3*gm(2,3)&
 713 &   *dgm(1,3))+gm(1,3)*(3*gm(2,3)*dgm(1,1)+6*gm(1,2)*dgm(1,3)+3*gm(1,1)&
 714 &   *dgm(2,3))+12*gm(1,2)*(gm(3,3)*dgm(1,1)+gm(1,1)*dgm(3,3))
 715    cm(2,6,9)=3*gm(1,2)**2*dgm(1,3)+gm(1,1)*(3*gm(2,3)*dgm(1,2)+12*gm(2,2)&
 716 &   *dgm(1,3))+gm(1,3)*(12*gm(2,2)*dgm(1,1)+6*gm(1,2)*dgm(1,2)+12*gm(1,1)&
 717 &   *dgm(2,2))+3*gm(1,2)*(gm(2,3)*dgm(1,1)+gm(1,1)*dgm(2,3))
 718    cm(2,7,9)=-3*gm(2,2)**2*dgm(1,3)+9*gm(1,2)*gm(2,3)*dgm(2,2)+gm(2,2)&
 719 &   *(9*gm(2,3)*dgm(1,2)-6*gm(1,3)*dgm(2,2)+9*gm(1,2)*dgm(2,3))
 720    cm(2,8,9)=3*gm(2,3)**2*dgm(1,3)+gm(3,3)*(12*gm(1,3)*dgm(2,2)+3*gm(1,2)&
 721 &   *dgm(2,3))+gm(2,3)*(3*gm(3,3)*dgm(1,2)+6*gm(1,3)*dgm(2,3)+3*gm(1,2)&
 722 &   *dgm(3,3))+12*gm(2,2)*(gm(3,3)*dgm(1,3)+gm(1,3)*dgm(3,3))
 723    cm(2,9,9)=3*gm(2,3)**2*dgm(1,2)+12*gm(1,2)*gm(3,3)*dgm(2,2)+gm(2,3)&
 724 &   *(3*gm(2,2)*dgm(1,3)+3*gm(1,3)*dgm(2,2)+6*gm(1,2)*dgm(2,3))+gm(2,2)&
 725 &   *(12*gm(3,3)*dgm(1,2)+3*gm(1,3)*dgm(2,3)+12*gm(1,2)*dgm(3,3))
 726    cm(2,10,9)=-3*gm(3,3)**2*dgm(1,2)+9*gm(1,3)*gm(2,3)*dgm(3,3)+gm(3,3)&
 727 &   *(9*gm(2,3)*dgm(1,3)+9*gm(1,3)*dgm(2,3)-6*gm(1,2)*dgm(3,3))
 728    cm(2,1,11)=7.5d0*gm(1,2)**2*dgm(1,2)-1.5d0*gm(1,1)*gm(2,2)*dgm(1,2)&
 729 &   -1.5d0*gm(1,2)*(gm(2,2)*dgm(1,1)+gm(1,1)*dgm(2,2))
 730    cm(2,2,11)=gm(2,2)*(3*gm(2,2)*dgm(1,2)+6*gm(1,2)*dgm(2,2))
 731    cm(2,3,11)=7.5d0*gm(2,3)**2*dgm(1,2)-1.5d0*gm(1,2)*gm(3,3)*dgm(2,2)&
 732 &   +gm(2,3)*(-3*gm(2,2)*dgm(1,3)-3*gm(1,3)*dgm(2,2)+15*gm(1,2)*dgm(2,3))&
 733 &   +gm(2,2)*(-1.5d0*gm(3,3)*dgm(1,2)-3*gm(1,3)*dgm(2,3)-1.5d0*gm(1,2)&
 734 &   *dgm(3,3))
 735    cm(2,4,11)=-3*gm(2,2)**2*dgm(1,3)+9*gm(1,2)*gm(2,3)*dgm(2,2)+gm(2,2)&
 736 &   *(9*gm(2,3)*dgm(1,2)-6*gm(1,3)*dgm(2,2)+9*gm(1,2)*dgm(2,3))
 737    cm(2,5,11)=-1.5d0*gm(1,1)*gm(2,3)*dgm(2,2)+gm(1,2)*(15*gm(2,3)&
 738 &   *dgm(1,2)-3*gm(1,3)*dgm(2,2))+7.5d0*gm(1,2)**2*dgm(2,3)+gm(2,2)&
 739 &   *(-1.5d0*gm(2,3)*dgm(1,1)-3*gm(1,3)*dgm(1,2)-3*gm(1,2)*dgm(1,3)&
 740 &   -1.5d0*gm(1,1)*dgm(2,3))
 741    cm(2,6,11)=-1.5d0*gm(2,2)**2*dgm(1,1)+4.5d0*gm(1,2)**2*dgm(2,2)&
 742 &   +gm(2,2)*(9*gm(1,2)*dgm(1,2)-3*gm(1,1)*dgm(2,2))
 743    cm(2,7,11)=3*gm(2,2)**2*dgm(2,2)
 744    cm(2,8,11)=4.5d0*gm(2,3)**2*dgm(2,2)+9*gm(2,2)*gm(2,3)*dgm(2,3)&
 745 &   +gm(2,2)*(-3*gm(3,3)*dgm(2,2)-1.5d0*gm(2,2)*dgm(3,3))
 746    cm(2,9,11)=gm(2,2)*(6*gm(2,3)*dgm(2,2)+3*gm(2,2)*dgm(2,3))
 747    cm(2,10,11)=7.5d0*gm(2,3)**2*dgm(2,3)-1.5d0*gm(2,2)*gm(3,3)*dgm(2,3)&
 748 &   -1.5d0*gm(2,3)*(gm(3,3)*dgm(2,2)+gm(2,2)*dgm(3,3))
 749    cm(2,1,12)=7.5d0*gm(1,3)**2*dgm(1,2)+gm(1,1)*(-1.5d0*gm(3,3)*dgm(1,2)&
 750 &   -3*gm(2,3)*dgm(1,3))+gm(1,3)*(-3*gm(2,3)*dgm(1,1)+15*gm(1,2)&
 751 &   *dgm(1,3)-3*gm(1,1)*dgm(2,3))-1.5d0*gm(1,2)*(gm(3,3)*dgm(1,1)&
 752 &   +gm(1,1)*dgm(3,3))
 753    cm(2,2,12)=1.5d0*gm(2,3)**2*dgm(1,2)-4.5d0*gm(1,2)*gm(3,3)*dgm(2,2)&
 754 &   +gm(2,3)*(12*gm(2,2)*dgm(1,3)+12*gm(1,3)*dgm(2,2)+3*gm(1,2)*dgm(2,3))&
 755 &   +gm(2,2)*(-4.5d0*gm(3,3)*dgm(1,2)+12*gm(1,3)*dgm(2,3)-4.5d0*gm(1,2)&
 756 &   *dgm(3,3))
 757    cm(2,3,12)=6*gm(3,3)**2*dgm(1,2)+3*gm(1,3)*gm(2,3)*dgm(3,3)+gm(3,3)&
 758 &   *(3*gm(2,3)*dgm(1,3)+3*gm(1,3)*dgm(2,3)+12*gm(1,2)*dgm(3,3))
 759    cm(2,4,12)=3*gm(2,3)**2*dgm(1,3)+gm(3,3)*(12*gm(1,3)*dgm(2,2)&
 760 &   +3*gm(1,2)*dgm(2,3))+gm(2,3)*(3*gm(3,3)*dgm(1,2)+6*gm(1,3)*dgm(2,3)&
 761 &   +3*gm(1,2)*dgm(3,3))+12*gm(2,2)*(gm(3,3)*dgm(1,3)+gm(1,3)*dgm(3,3))
 762    cm(2,5,12)=1.5d0*gm(1,3)**2*dgm(2,3)+gm(3,3)*(12*gm(1,2)*dgm(1,3)&
 763 &   -4.5d0*gm(1,1)*dgm(2,3))+gm(2,3)*(-4.5d0*gm(3,3)*dgm(1,1)+3*gm(1,3)&
 764 &   *dgm(1,3)-4.5d0*gm(1,1)*dgm(3,3))+12*gm(1,3)*(gm(3,3)*dgm(1,2)&
 765 &   +gm(1,2)*dgm(3,3))
 766    cm(2,6,12)=-3*gm(2,3)**2*dgm(1,1)-6*gm(1,2)*gm(3,3)*dgm(1,2)+7.5d0*gm(1,3)&
 767 &   **2*dgm(2,2)-1.5d0*gm(1,1)*gm(3,3)*dgm(2,2)+9*gm(1,2)*gm(1,3)&
 768 &   *dgm(2,3)+gm(2,3)*(9*gm(1,3)*dgm(1,2)+9*gm(1,2)*dgm(1,3)-6*gm(1,1)&
 769 &   *dgm(2,3))-3*gm(1,2)**2*dgm(3,3)+gm(2,2)*(-1.5d0*gm(3,3)*dgm(1,1)&
 770 &   +15*gm(1,3)*dgm(1,3)-1.5d0*gm(1,1)*dgm(3,3))
 771    cm(2,7,12)=4.5d0*gm(2,3)**2*dgm(2,2)+9*gm(2,2)*gm(2,3)*dgm(2,3)&
 772 &   +gm(2,2)*(-3*gm(3,3)*dgm(2,2)-1.5d0*gm(2,2)*dgm(3,3))
 773    cm(2,8,12)=6*gm(3,3)**2*dgm(2,2)+3*gm(2,3)**2*dgm(3,3)+gm(3,3)&
 774 &   *(6*gm(2,3)*dgm(2,3)+12*gm(2,2)*dgm(3,3))
 775    cm(2,9,12)=4.5d0*gm(2,3)**2*dgm(2,3)+7.5d0*gm(2,2)*gm(3,3)*dgm(2,3)&
 776 &   +7.5d0*gm(2,3)*(gm(3,3)*dgm(2,2)+gm(2,2)*dgm(3,3))
 777    cm(2,10,12)=gm(3,3)*(3*gm(3,3)*dgm(2,3)+6*gm(2,3)*dgm(3,3))
 778    cm(2,1,13)=7.5d0*gm(1,2)**2*dgm(1,3)+gm(1,1)*(-3*gm(2,3)*dgm(1,2)&
 779 &   -1.5d0*gm(2,2)*dgm(1,3))+gm(1,3)*(-1.5d0*gm(2,2)*dgm(1,1)+15*gm(1,2)&
 780 &   *dgm(1,2)-1.5d0*gm(1,1)*dgm(2,2))-3*gm(1,2)*(gm(2,3)*dgm(1,1)&
 781 &   +gm(1,1)*dgm(2,3))
 782    cm(2,2,13)=6*gm(2,2)**2*dgm(1,3)+3*gm(1,2)*gm(2,3)*dgm(2,2)+gm(2,2)&
 783 &   *(3*gm(2,3)*dgm(1,2)+12*gm(1,3)*dgm(2,2)+3*gm(1,2)*dgm(2,3))
 784    cm(2,3,13)=1.5d0*gm(2,3)**2*dgm(1,3)+gm(3,3)*(-4.5d0*gm(1,3)*dgm(2,2)&
 785 &   +12*gm(1,2)*dgm(2,3))+gm(2,3)*(12*gm(3,3)*dgm(1,2)+3*gm(1,3)&
 786 &   *dgm(2,3)+12*gm(1,2)*dgm(3,3))-4.5d0*gm(2,2)*(gm(3,3)*dgm(1,3)&
 787 &   +gm(1,3)*dgm(3,3))
 788    cm(2,4,13)=3*gm(2,3)**2*dgm(1,2)+12*gm(1,2)*gm(3,3)*dgm(2,2)+gm(2,3)&
 789 &   *(3*gm(2,2)*dgm(1,3)+3*gm(1,3)*dgm(2,2)+6*gm(1,2)*dgm(2,3))+gm(2,2)&
 790 &   *(12*gm(3,3)*dgm(1,2)+3*gm(1,3)*dgm(2,3)+12*gm(1,2)*dgm(3,3))
 791    cm(2,5,13)=-3*gm(2,3)**2*dgm(1,1)+15*gm(1,2)*gm(3,3)*dgm(1,2)&
 792 &   -3*gm(1,3)**2*dgm(2,2)-1.5d0*gm(1,1)*gm(3,3)*dgm(2,2)+9*gm(1,2)&
 793 &   *gm(1,3)*dgm(2,3)+gm(2,3)*(9*gm(1,3)*dgm(1,2)+9*gm(1,2)*dgm(1,3)&
 794 &   -6*gm(1,1)*dgm(2,3))+7.5d0*gm(1,2)**2*dgm(3,3)+gm(2,2)*(-1.5d0*gm(3,3)&
 795 &   *dgm(1,1)-6*gm(1,3)*dgm(1,3)-1.5d0*gm(1,1)*dgm(3,3))
 796    cm(2,6,13)=-4.5d0*gm(1,1)*gm(2,3)*dgm(2,2)+gm(1,2)*(3*gm(2,3)&
 797 &   *dgm(1,2)+12*gm(1,3)*dgm(2,2))+1.5d0*gm(1,2)**2*dgm(2,3)+gm(2,2)&
 798 &   *(-4.5d0*gm(2,3)*dgm(1,1)+12*gm(1,3)*dgm(1,2)+12*gm(1,2)*dgm(1,3)&
 799 &   -4.5d0*gm(1,1)*dgm(2,3))
 800    cm(2,7,13)=gm(2,2)*(6*gm(2,3)*dgm(2,2)+3*gm(2,2)*dgm(2,3))
 801    cm(2,8,13)=4.5d0*gm(2,3)**2*dgm(2,3)+7.5d0*gm(2,2)*gm(3,3)*dgm(2,3)&
 802 &   +7.5d0*gm(2,3)*(gm(3,3)*dgm(2,2)+gm(2,2)*dgm(3,3))
 803    cm(2,9,13)=3*gm(2,3)**2*dgm(2,2)+6*gm(2,2)*gm(2,3)*dgm(2,3)+gm(2,2)&
 804 &   *(12*gm(3,3)*dgm(2,2)+6*gm(2,2)*dgm(3,3))
 805    cm(2,10,13)=-1.5d0*gm(3,3)**2*dgm(2,2)+4.5d0*gm(2,3)**2*dgm(3,3)&
 806 &   +gm(3,3)*(9*gm(2,3)*dgm(2,3)-3*gm(2,2)*dgm(3,3))
 807    cm(2,1,14)=7.5d0*gm(1,3)**2*dgm(1,3)-1.5d0*gm(1,1)*gm(3,3)*dgm(1,3)&
 808 &   -1.5d0*gm(1,3)*(gm(3,3)*dgm(1,1)+gm(1,1)*dgm(3,3))
 809    cm(2,2,14)=7.5d0*gm(2,3)**2*dgm(1,3)+gm(3,3)*(-1.5d0*gm(1,3)*dgm(2,2)&
 810 &   -3*gm(1,2)*dgm(2,3))+gm(2,3)*(-3*gm(3,3)*dgm(1,2)+15*gm(1,3)&
 811 &   *dgm(2,3)-3*gm(1,2)*dgm(3,3))-1.5d0*gm(2,2)*(gm(3,3)*dgm(1,3)&
 812 &   +gm(1,3)*dgm(3,3))
 813    cm(2,3,14)=gm(3,3)*(3*gm(3,3)*dgm(1,3)+6*gm(1,3)*dgm(3,3))
 814    cm(2,4,14)=-3*gm(3,3)**2*dgm(1,2)+9*gm(1,3)*gm(2,3)*dgm(3,3)+gm(3,3)&
 815 &   *(9*gm(2,3)*dgm(1,3)+9*gm(1,3)*dgm(2,3)-6*gm(1,2)*dgm(3,3))
 816    cm(2,5,14)=-1.5d0*gm(3,3)**2*dgm(1,1)+4.5d0*gm(1,3)**2*dgm(3,3)&
 817 &   +gm(3,3)*(9*gm(1,3)*dgm(1,3)-3*gm(1,1)*dgm(3,3))
 818    cm(2,6,14)=7.5d0*gm(1,3)**2*dgm(2,3)+gm(3,3)*(-3*gm(1,2)*dgm(1,3)&
 819 &   -1.5d0*gm(1,1)*dgm(2,3))+gm(2,3)*(-1.5d0*gm(3,3)*dgm(1,1)+15*gm(1,3)&
 820 &   *dgm(1,3)-1.5d0*gm(1,1)*dgm(3,3))-3*gm(1,3)*(gm(3,3)*dgm(1,2)&
 821 &   +gm(1,2)*dgm(3,3))
 822    cm(2,7,14)=7.5d0*gm(2,3)**2*dgm(2,3)-1.5d0*gm(2,2)*gm(3,3)*dgm(2,3)&
 823 &   -1.5d0*gm(2,3)*(gm(3,3)*dgm(2,2)+gm(2,2)*dgm(3,3))
 824    cm(2,8,14)=gm(3,3)*(3*gm(3,3)*dgm(2,3)+6*gm(2,3)*dgm(3,3))
 825    cm(2,9,14)=-1.5d0*gm(3,3)**2*dgm(2,2)+4.5d0*gm(2,3)**2*dgm(3,3)&
 826 &   +gm(3,3)*(9*gm(2,3)*dgm(2,3)-3*gm(2,2)*dgm(3,3))
 827    cm(2,10,14)=3*gm(3,3)**2*dgm(3,3)
 828    cm(3,1,3)=gm(1,1)*(6*gm(1,3)*dgm(1,1)+3*gm(1,1)*dgm(1,3))
 829    cm(3,2,3)=1.5d0*gm(1,2)**2*dgm(1,3)+gm(1,1)*(12*gm(2,3)*dgm(1,2)&
 830 &   -4.5d0*gm(2,2)*dgm(1,3))+gm(1,3)*(-4.5d0*gm(2,2)*dgm(1,1)+3*gm(1,2)&
 831 &   *dgm(1,2)-4.5d0*gm(1,1)*dgm(2,2))+12*gm(1,2)*(gm(2,3)*dgm(1,1)&
 832 &   +gm(1,1)*dgm(2,3))
 833    cm(3,3,3)=4.5d0*gm(1,3)**2*dgm(1,3)+7.5d0*gm(1,1)*gm(3,3)*dgm(1,3)&
 834 &   +7.5d0*gm(1,3)*(gm(3,3)*dgm(1,1)+gm(1,1)*dgm(3,3))
 835    cm(3,4,3)=3*gm(1,3)**2*dgm(1,2)+gm(1,1)*(12*gm(3,3)*dgm(1,2)+3*gm(2,3)&
 836 &   *dgm(1,3))+gm(1,3)*(3*gm(2,3)*dgm(1,1)+6*gm(1,2)*dgm(1,3)+3*gm(1,1)&
 837 &   *dgm(2,3))+12*gm(1,2)*(gm(3,3)*dgm(1,1)+gm(1,1)*dgm(3,3))
 838    cm(3,5,3)=3*gm(1,3)**2*dgm(1,1)+6*gm(1,1)*gm(1,3)*dgm(1,3)+gm(1,1)&
 839 &   *(12*gm(3,3)*dgm(1,1)+6*gm(1,1)*dgm(3,3))
 840    cm(3,6,3)=3*gm(1,2)*(gm(1,3)*dgm(1,1)+gm(1,1)*dgm(1,3))+gm(1,1)&
 841 &   *(12*gm(2,3)*dgm(1,1)+3*gm(1,3)*dgm(1,2)+6*gm(1,1)*dgm(2,3))
 842    cm(3,7,3)=-1.5d0*gm(1,1)*gm(2,3)*dgm(2,2)+gm(1,2)*(15*gm(2,3)&
 843 &   *dgm(1,2)-3*gm(1,3)*dgm(2,2))+7.5d0*gm(1,2)**2*dgm(2,3)+gm(2,2)&
 844 &   *(-1.5d0*gm(2,3)*dgm(1,1)-3*gm(1,3)*dgm(1,2)-3*gm(1,2)*dgm(1,3)&
 845 &   -1.5d0*gm(1,1)*dgm(2,3))
 846    cm(3,8,3)=1.5d0*gm(1,3)**2*dgm(2,3)+gm(3,3)*(12*gm(1,2)*dgm(1,3)&
 847 &   -4.5d0*gm(1,1)*dgm(2,3))+gm(2,3)*(-4.5d0*gm(3,3)*dgm(1,1)+3*gm(1,3)&
 848 &   *dgm(1,3)-4.5d0*gm(1,1)*dgm(3,3))+12*gm(1,3)*(gm(3,3)*dgm(1,2)&
 849 &   +gm(1,2)*dgm(3,3))
 850    cm(3,9,3)=-3*gm(2,3)**2*dgm(1,1)+15*gm(1,2)*gm(3,3)*dgm(1,2)-3*gm(1,3)&
 851 &   **2*dgm(2,2)-1.5d0*gm(1,1)*gm(3,3)*dgm(2,2)+9*gm(1,2)*gm(1,3)&
 852 &   *dgm(2,3)+gm(2,3)*(9*gm(1,3)*dgm(1,2)+9*gm(1,2)*dgm(1,3)-6*gm(1,1)&
 853 &   *dgm(2,3))+7.5d0*gm(1,2)**2*dgm(3,3)+gm(2,2)*(-1.5d0*gm(3,3)&
 854 &   *dgm(1,1)-6*gm(1,3)*dgm(1,3)-1.5d0*gm(1,1)*dgm(3,3))
 855    cm(3,10,3)=-1.5d0*gm(3,3)**2*dgm(1,1)+4.5d0*gm(1,3)**2*dgm(3,3)&
 856 &   +gm(3,3)*(9*gm(1,3)*dgm(1,3)-3*gm(1,1)*dgm(3,3))
 857    cm(3,1,4)=gm(1,1)*(6*gm(1,2)*dgm(1,1)+3*gm(1,1)*dgm(1,2))
 858    cm(3,2,4)=4.5d0*gm(1,2)**2*dgm(1,2)+7.5d0*gm(1,1)*gm(2,2)*dgm(1,2)&
 859 &   +7.5d0*gm(1,2)*(gm(2,2)*dgm(1,1)+gm(1,1)*dgm(2,2))
 860    cm(3,3,4)=1.5d0*gm(1,3)**2*dgm(1,2)+gm(1,1)*(-4.5d0*gm(3,3)*dgm(1,2)&
 861 &   +12*gm(2,3)*dgm(1,3))+gm(1,3)*(12*gm(2,3)*dgm(1,1)+3*gm(1,2)&
 862 &   *dgm(1,3)+12*gm(1,1)*dgm(2,3))-4.5d0*gm(1,2)*(gm(3,3)*dgm(1,1)&
 863 &   +gm(1,1)*dgm(3,3))
 864    cm(3,4,4)=3*gm(1,2)**2*dgm(1,3)+gm(1,1)*(3*gm(2,3)*dgm(1,2)+12*gm(2,2)&
 865 &   *dgm(1,3))+gm(1,3)*(12*gm(2,2)*dgm(1,1)+6*gm(1,2)*dgm(1,2)+12*gm(1,1)&
 866 &   *dgm(2,2))+3*gm(1,2)*(gm(2,3)*dgm(1,1)+gm(1,1)*dgm(2,3))
 867    cm(3,5,4)=3*gm(1,2)*(gm(1,3)*dgm(1,1)+gm(1,1)*dgm(1,3))+gm(1,1)&
 868 &   *(12*gm(2,3)*dgm(1,1)+3*gm(1,3)*dgm(1,2)+6*gm(1,1)*dgm(2,3))
 869    cm(3,6,4)=3*gm(1,2)**2*dgm(1,1)+6*gm(1,1)*gm(1,2)*dgm(1,2)+gm(1,1)&
 870 &   *(12*gm(2,2)*dgm(1,1)+6*gm(1,1)*dgm(2,2))
 871    cm(3,7,4)=-1.5d0*gm(2,2)**2*dgm(1,1)+4.5d0*gm(1,2)**2*dgm(2,2)&
 872 &   +gm(2,2)*(9*gm(1,2)*dgm(1,2)-3*gm(1,1)*dgm(2,2))
 873    cm(3,8,4)=-3*gm(2,3)**2*dgm(1,1)-6*gm(1,2)*gm(3,3)*dgm(1,2)+7.5d0*gm(1,3)&
 874 &   **2*dgm(2,2)-1.5d0*gm(1,1)*gm(3,3)*dgm(2,2)+9*gm(1,2)*gm(1,3)&
 875 &   *dgm(2,3)+gm(2,3)*(9*gm(1,3)*dgm(1,2)+9*gm(1,2)*dgm(1,3)-6*gm(1,1)&
 876 &   *dgm(2,3))-3*gm(1,2)**2*dgm(3,3)+gm(2,2)*(-1.5d0*gm(3,3)*dgm(1,1)&
 877 &   +15*gm(1,3)*dgm(1,3)-1.5d0*gm(1,1)*dgm(3,3))
 878    cm(3,9,4)=-4.5d0*gm(1,1)*gm(2,3)*dgm(2,2)+gm(1,2)*(3*gm(2,3)*dgm(1,2)&
 879 &   +12*gm(1,3)*dgm(2,2))+1.5d0*gm(1,2)**2*dgm(2,3)+gm(2,2)*(-4.5d0*gm(2,3)&
 880 &   *dgm(1,1)+12*gm(1,3)*dgm(1,2)+12*gm(1,2)*dgm(1,3)-4.5d0*gm(1,1)&
 881 &   *dgm(2,3))
 882    cm(3,10,4)=7.5d0*gm(1,3)**2*dgm(2,3)+gm(3,3)*(-3*gm(1,2)*dgm(1,3)&
 883 &   -1.5d0*gm(1,1)*dgm(2,3))+gm(2,3)*(-1.5d0*gm(3,3)*dgm(1,1)+15*gm(1,3)&
 884 &   *dgm(1,3)-1.5d0*gm(1,1)*dgm(3,3))-3*gm(1,3)*(gm(3,3)*dgm(1,2)&
 885 &   +gm(1,2)*dgm(3,3))
 886    cm(3,1,5)=3*gm(1,1)**2*dgm(1,1)
 887    cm(3,2,5)=4.5d0*gm(1,2)**2*dgm(1,1)+9*gm(1,1)*gm(1,2)*dgm(1,2)&
 888 &   +gm(1,1)*(-3*gm(2,2)*dgm(1,1)-1.5d0*gm(1,1)*dgm(2,2))
 889    cm(3,3,5)=4.5d0*gm(1,3)**2*dgm(1,1)+9*gm(1,1)*gm(1,3)*dgm(1,3)&
 890 &   +gm(1,1)*(-3*gm(3,3)*dgm(1,1)-1.5d0*gm(1,1)*dgm(3,3))
 891    cm(3,4,5)=9*gm(1,2)*(gm(1,3)*dgm(1,1)+gm(1,1)*dgm(1,3))+gm(1,1)&
 892 &   *(-6*gm(2,3)*dgm(1,1)+9*gm(1,3)*dgm(1,2)-3*gm(1,1)*dgm(2,3))
 893    cm(3,5,5)=gm(1,1)*(6*gm(1,3)*dgm(1,1)+3*gm(1,1)*dgm(1,3))
 894    cm(3,6,5)=gm(1,1)*(6*gm(1,2)*dgm(1,1)+3*gm(1,1)*dgm(1,2))
 895    cm(3,7,5)=7.5d0*gm(1,2)**2*dgm(1,2)-1.5d0*gm(1,1)*gm(2,2)*dgm(1,2)&
 896 &   -1.5d0*gm(1,2)*(gm(2,2)*dgm(1,1)+gm(1,1)*dgm(2,2))
 897    cm(3,8,5)=7.5d0*gm(1,3)**2*dgm(1,2)+gm(1,1)*(-1.5d0*gm(3,3)*dgm(1,2)&
 898 &   -3*gm(2,3)*dgm(1,3))+gm(1,3)*(-3*gm(2,3)*dgm(1,1)+15*gm(1,2)&
 899 &   *dgm(1,3)-3*gm(1,1)*dgm(2,3))-1.5d0*gm(1,2)*(gm(3,3)*dgm(1,1)&
 900 &   +gm(1,1)*dgm(3,3))
 901    cm(3,9,5)=7.5d0*gm(1,2)**2*dgm(1,3)+gm(1,1)*(-3*gm(2,3)*dgm(1,2)&
 902 &   -1.5d0*gm(2,2)*dgm(1,3))+gm(1,3)*(-1.5d0*gm(2,2)*dgm(1,1)+15*gm(1,2)&
 903 &   *dgm(1,2)-1.5d0*gm(1,1)*dgm(2,2))-3*gm(1,2)*(gm(2,3)*dgm(1,1)&
 904 &   +gm(1,1)*dgm(2,3))
 905    cm(3,10,5)=7.5d0*gm(1,3)**2*dgm(1,3)-1.5d0*gm(1,1)*gm(3,3)*dgm(1,3)&
 906 &   -1.5d0*gm(1,3)*(gm(3,3)*dgm(1,1)+gm(1,1)*dgm(3,3))
 907    cm(3,1,8)=9*gm(1,2)*(gm(1,3)*dgm(1,1)+gm(1,1)*dgm(1,3))+gm(1,1)&
 908 &   *(-6*gm(2,3)*dgm(1,1)+9*gm(1,3)*dgm(1,2)-3*gm(1,1)*dgm(2,3))
 909    cm(3,2,8)=12*gm(1,1)*gm(2,3)*dgm(2,2)+gm(1,2)*(6*gm(2,3)*dgm(1,2)&
 910 &   +3*gm(1,3)*dgm(2,2))+3*gm(1,2)**2*dgm(2,3)+gm(2,2)*(12*gm(2,3)&
 911 &   *dgm(1,1)+3*gm(1,3)*dgm(1,2)+3*gm(1,2)*dgm(1,3)+12*gm(1,1)*dgm(2,3))
 912    cm(3,3,8)=3*gm(1,3)**2*dgm(2,3)+gm(3,3)*(3*gm(1,2)*dgm(1,3)+12*gm(1,1)&
 913 &   *dgm(2,3))+gm(2,3)*(12*gm(3,3)*dgm(1,1)+6*gm(1,3)*dgm(1,3)+12*gm(1,1)&
 914 &   *dgm(3,3))+3*gm(1,3)*(gm(3,3)*dgm(1,2)+gm(1,2)*dgm(3,3))
 915    cm(3,4,8)=9*gm(2,3)**2*dgm(1,1)+18*gm(1,2)*gm(3,3)*dgm(1,2)+9*gm(1,3)&
 916 &   **2*dgm(2,2)+15*gm(1,1)*gm(3,3)*dgm(2,2)-6*gm(1,2)*gm(1,3)*dgm(2,3)&
 917 &   +gm(2,3)*(-6*gm(1,3)*dgm(1,2)-6*gm(1,2)*dgm(1,3)+18*gm(1,1)*dgm(2,3))&
 918 &   +9*gm(1,2)**2*dgm(3,3)+gm(2,2)*(15*gm(3,3)*dgm(1,1)+18*gm(1,3)&
 919 &   *dgm(1,3)+15*gm(1,1)*dgm(3,3))
 920    cm(3,5,8)=3*gm(1,3)**2*dgm(1,2)+gm(1,1)*(12*gm(3,3)*dgm(1,2)+3*gm(2,3)&
 921 &   *dgm(1,3))+gm(1,3)*(3*gm(2,3)*dgm(1,1)+6*gm(1,2)*dgm(1,3)+3*gm(1,1)&
 922 &   *dgm(2,3))+12*gm(1,2)*(gm(3,3)*dgm(1,1)+gm(1,1)*dgm(3,3))
 923    cm(3,6,8)=3*gm(1,2)**2*dgm(1,3)+gm(1,1)*(3*gm(2,3)*dgm(1,2)+12*gm(2,2)&
 924 &   *dgm(1,3))+gm(1,3)*(12*gm(2,2)*dgm(1,1)+6*gm(1,2)*dgm(1,2)+12*gm(1,1)&
 925 &   *dgm(2,2))+3*gm(1,2)*(gm(2,3)*dgm(1,1)+gm(1,1)*dgm(2,3))
 926    cm(3,7,8)=-3*gm(2,2)**2*dgm(1,3)+9*gm(1,2)*gm(2,3)*dgm(2,2)+gm(2,2)&
 927 &   *(9*gm(2,3)*dgm(1,2)-6*gm(1,3)*dgm(2,2)+9*gm(1,2)*dgm(2,3))
 928    cm(3,8,8)=3*gm(2,3)**2*dgm(1,3)+gm(3,3)*(12*gm(1,3)*dgm(2,2)+3*gm(1,2)&
 929 &   *dgm(2,3))+gm(2,3)*(3*gm(3,3)*dgm(1,2)+6*gm(1,3)*dgm(2,3)+3*gm(1,2)&
 930 &   *dgm(3,3))+12*gm(2,2)*(gm(3,3)*dgm(1,3)+gm(1,3)*dgm(3,3))
 931    cm(3,9,8)=3*gm(2,3)**2*dgm(1,2)+12*gm(1,2)*gm(3,3)*dgm(2,2)+gm(2,3)&
 932 &   *(3*gm(2,2)*dgm(1,3)+3*gm(1,3)*dgm(2,2)+6*gm(1,2)*dgm(2,3))+gm(2,2)&
 933 &   *(12*gm(3,3)*dgm(1,2)+3*gm(1,3)*dgm(2,3)+12*gm(1,2)*dgm(3,3))
 934    cm(3,10,8)=-3*gm(3,3)**2*dgm(1,2)+9*gm(1,3)*gm(2,3)*dgm(3,3)+gm(3,3)&
 935 &   *(9*gm(2,3)*dgm(1,3)+9*gm(1,3)*dgm(2,3)-6*gm(1,2)*dgm(3,3))
 936    cm(3,1,9)=4.5d0*gm(1,2)**2*dgm(1,1)+9*gm(1,1)*gm(1,2)*dgm(1,2)&
 937 &   +gm(1,1)*(-3*gm(2,2)*dgm(1,1)-1.5d0*gm(1,1)*dgm(2,2))
 938    cm(3,2,9)=6*gm(2,2)**2*dgm(1,1)+3*gm(1,2)**2*dgm(2,2)+gm(2,2)&
 939 &   *(6*gm(1,2)*dgm(1,2)+12*gm(1,1)*dgm(2,2))
 940    cm(3,3,9)=7.5d0*gm(2,3)**2*dgm(1,1)-6*gm(1,2)*gm(3,3)*dgm(1,2)&
 941 &   -3*gm(1,3)**2*dgm(2,2)-1.5d0*gm(1,1)*gm(3,3)*dgm(2,2)+9*gm(1,2)&
 942 &   *gm(1,3)*dgm(2,3)+gm(2,3)*(9*gm(1,3)*dgm(1,2)+9*gm(1,2)*dgm(1,3)&
 943 &   +15*gm(1,1)*dgm(2,3))-3*gm(1,2)**2*dgm(3,3)+gm(2,2)*(-1.5d0*gm(3,3)&
 944 &   *dgm(1,1)-6*gm(1,3)*dgm(1,3)-1.5d0*gm(1,1)*dgm(3,3))
 945    cm(3,4,9)=12*gm(1,1)*gm(2,3)*dgm(2,2)+gm(1,2)*(6*gm(2,3)*dgm(1,2)&
 946 &   +3*gm(1,3)*dgm(2,2))+3*gm(1,2)**2*dgm(2,3)+gm(2,2)*(12*gm(2,3)&
 947 &   *dgm(1,1)+3*gm(1,3)*dgm(1,2)+3*gm(1,2)*dgm(1,3)+12*gm(1,1)*dgm(2,3))
 948    cm(3,5,9)=1.5d0*gm(1,2)**2*dgm(1,3)+gm(1,1)*(12*gm(2,3)*dgm(1,2)&
 949 &   -4.5d0*gm(2,2)*dgm(1,3))+gm(1,3)*(-4.5d0*gm(2,2)*dgm(1,1)+3*gm(1,2)&
 950 &   *dgm(1,2)-4.5d0*gm(1,1)*dgm(2,2))+12*gm(1,2)*(gm(2,3)*dgm(1,1)&
 951 &   +gm(1,1)*dgm(2,3))
 952    cm(3,6,9)=4.5d0*gm(1,2)**2*dgm(1,2)+7.5d0*gm(1,1)*gm(2,2)*dgm(1,2)&
 953 &   +7.5d0*gm(1,2)*(gm(2,2)*dgm(1,1)+gm(1,1)*dgm(2,2))
 954    cm(3,7,9)=gm(2,2)*(3*gm(2,2)*dgm(1,2)+6*gm(1,2)*dgm(2,2))
 955    cm(3,8,9)=1.5d0*gm(2,3)**2*dgm(1,2)-4.5d0*gm(1,2)*gm(3,3)*dgm(2,2)&
 956 &   +gm(2,3)*(12*gm(2,2)*dgm(1,3)+12*gm(1,3)*dgm(2,2)+3*gm(1,2)*dgm(2,3))&
 957 &   +gm(2,2)*(-4.5d0*gm(3,3)*dgm(1,2)+12*gm(1,3)*dgm(2,3)-4.5d0*gm(1,2)&
 958 &   *dgm(3,3))
 959    cm(3,9,9)=6*gm(2,2)**2*dgm(1,3)+3*gm(1,2)*gm(2,3)*dgm(2,2)+gm(2,2)&
 960 &   *(3*gm(2,3)*dgm(1,2)+12*gm(1,3)*dgm(2,2)+3*gm(1,2)*dgm(2,3))
 961    cm(3,10,9)=7.5d0*gm(2,3)**2*dgm(1,3)+gm(3,3)*(-1.5d0*gm(1,3)*dgm(2,2)&
 962 &   -3*gm(1,2)*dgm(2,3))+gm(2,3)*(-3*gm(3,3)*dgm(1,2)+15*gm(1,3)&
 963 &   *dgm(2,3)-3*gm(1,2)*dgm(3,3))-1.5d0*gm(2,2)*(gm(3,3)*dgm(1,3)&
 964 &   +gm(1,3)*dgm(3,3))
 965    cm(3,1,10)=4.5d0*gm(1,3)**2*dgm(1,1)+9*gm(1,1)*gm(1,3)*dgm(1,3)&
 966 &   +gm(1,1)*(-3*gm(3,3)*dgm(1,1)-1.5d0*gm(1,1)*dgm(3,3))
 967    cm(3,2,10)=7.5d0*gm(2,3)**2*dgm(1,1)-6*gm(1,2)*gm(3,3)*dgm(1,2)&
 968 &   -3*gm(1,3)**2*dgm(2,2)-1.5d0*gm(1,1)*gm(3,3)*dgm(2,2)+9*gm(1,2)&
 969 &   *gm(1,3)*dgm(2,3)+gm(2,3)*(9*gm(1,3)*dgm(1,2)+9*gm(1,2)*dgm(1,3)&
 970 &   +15*gm(1,1)*dgm(2,3))-3*gm(1,2)**2*dgm(3,3)+gm(2,2)*(-1.5d0*gm(3,3)&
 971 &   *dgm(1,1)-6*gm(1,3)*dgm(1,3)-1.5d0*gm(1,1)*dgm(3,3))
 972    cm(3,3,10)=6*gm(3,3)**2*dgm(1,1)+3*gm(1,3)**2*dgm(3,3)+gm(3,3)&
 973 &   *(6*gm(1,3)*dgm(1,3)+12*gm(1,1)*dgm(3,3))
 974    cm(3,4,10)=3*gm(1,3)**2*dgm(2,3)+gm(3,3)*(3*gm(1,2)*dgm(1,3)+12*gm(1,1)&
 975 &   *dgm(2,3))+gm(2,3)*(12*gm(3,3)*dgm(1,1)+6*gm(1,3)*dgm(1,3)+12*gm(1,1)&
 976 &   *dgm(3,3))+3*gm(1,3)*(gm(3,3)*dgm(1,2)+gm(1,2)*dgm(3,3))
 977    cm(3,5,10)=4.5d0*gm(1,3)**2*dgm(1,3)+7.5d0*gm(1,1)*gm(3,3)*dgm(1,3)&
 978 &   +7.5d0*gm(1,3)*(gm(3,3)*dgm(1,1)+gm(1,1)*dgm(3,3))
 979    cm(3,6,10)=1.5d0*gm(1,3)**2*dgm(1,2)+gm(1,1)*(-4.5d0*gm(3,3)*dgm(1,2)&
 980 &   +12*gm(2,3)*dgm(1,3))+gm(1,3)*(12*gm(2,3)*dgm(1,1)+3*gm(1,2)&
 981 &   *dgm(1,3)+12*gm(1,1)*dgm(2,3))-4.5d0*gm(1,2)*(gm(3,3)*dgm(1,1)&
 982 &   +gm(1,1)*dgm(3,3))
 983    cm(3,7,10)=7.5d0*gm(2,3)**2*dgm(1,2)-1.5d0*gm(1,2)*gm(3,3)*dgm(2,2)&
 984 &   +gm(2,3)*(-3*gm(2,2)*dgm(1,3)-3*gm(1,3)*dgm(2,2)+15*gm(1,2)*dgm(2,3))&
 985 &   +gm(2,2)*(-1.5d0*gm(3,3)*dgm(1,2)-3*gm(1,3)*dgm(2,3)-1.5d0*gm(1,2)&
 986 &   *dgm(3,3))
 987    cm(3,8,10)=6*gm(3,3)**2*dgm(1,2)+3*gm(1,3)*gm(2,3)*dgm(3,3)+gm(3,3)&
 988 &   *(3*gm(2,3)*dgm(1,3)+3*gm(1,3)*dgm(2,3)+12*gm(1,2)*dgm(3,3))
 989    cm(3,9,10)=1.5d0*gm(2,3)**2*dgm(1,3)+gm(3,3)*(-4.5d0*gm(1,3)*dgm(2,2)&
 990 &   +12*gm(1,2)*dgm(2,3))+gm(2,3)*(12*gm(3,3)*dgm(1,2)+3*gm(1,3)&
 991 &   *dgm(2,3)+12*gm(1,2)*dgm(3,3))-4.5d0*gm(2,2)*(gm(3,3)*dgm(1,3)&
 992 &   +gm(1,3)*dgm(3,3))
 993    cm(3,10,10)=gm(3,3)*(3*gm(3,3)*dgm(1,3)+6*gm(1,3)*dgm(3,3))
 994    cm(3,1,12)=7.5d0*gm(1,2)**2*dgm(1,3)+gm(1,1)*(-3*gm(2,3)*dgm(1,2)&
 995 &   -1.5d0*gm(2,2)*dgm(1,3))+gm(1,3)*(-1.5d0*gm(2,2)*dgm(1,1)+15*gm(1,2)&
 996 &   *dgm(1,2)-1.5d0*gm(1,1)*dgm(2,2))-3*gm(1,2)*(gm(2,3)*dgm(1,1)&
 997 &   +gm(1,1)*dgm(2,3))
 998    cm(3,2,12)=6*gm(2,2)**2*dgm(1,3)+3*gm(1,2)*gm(2,3)*dgm(2,2)+gm(2,2)&
 999 &   *(3*gm(2,3)*dgm(1,2)+12*gm(1,3)*dgm(2,2)+3*gm(1,2)*dgm(2,3))
1000    cm(3,3,12)=1.5d0*gm(2,3)**2*dgm(1,3)+gm(3,3)*(-4.5d0*gm(1,3)*dgm(2,2)&
1001 &   +12*gm(1,2)*dgm(2,3))+gm(2,3)*(12*gm(3,3)*dgm(1,2)+3*gm(1,3)&
1002 &   *dgm(2,3)+12*gm(1,2)*dgm(3,3))-4.5d0*gm(2,2)*(gm(3,3)*dgm(1,3)&
1003 &   +gm(1,3)*dgm(3,3))
1004    cm(3,4,12)=3*gm(2,3)**2*dgm(1,2)+12*gm(1,2)*gm(3,3)*dgm(2,2)+gm(2,3)&
1005 &   *(3*gm(2,2)*dgm(1,3)+3*gm(1,3)*dgm(2,2)+6*gm(1,2)*dgm(2,3))+gm(2,2)&
1006 &   *(12*gm(3,3)*dgm(1,2)+3*gm(1,3)*dgm(2,3)+12*gm(1,2)*dgm(3,3))
1007    cm(3,5,12)=-3*gm(2,3)**2*dgm(1,1)+15*gm(1,2)*gm(3,3)*dgm(1,2)&
1008 &   -3*gm(1,3)**2*dgm(2,2)-1.5d0*gm(1,1)*gm(3,3)*dgm(2,2)+9*gm(1,2)&
1009 &   *gm(1,3)*dgm(2,3)+gm(2,3)*(9*gm(1,3)*dgm(1,2)+9*gm(1,2)*dgm(1,3)&
1010 &   -6*gm(1,1)*dgm(2,3))+7.5d0*gm(1,2)**2*dgm(3,3)+gm(2,2)*(-1.5d0*gm(3,3)&
1011 &   *dgm(1,1)-6*gm(1,3)*dgm(1,3)-1.5d0*gm(1,1)*dgm(3,3))
1012    cm(3,6,12)=-4.5d0*gm(1,1)*gm(2,3)*dgm(2,2)+gm(1,2)*(3*gm(2,3)&
1013 &   *dgm(1,2)+12*gm(1,3)*dgm(2,2))+1.5d0*gm(1,2)**2*dgm(2,3)+gm(2,2)&
1014 &   *(-4.5d0*gm(2,3)*dgm(1,1)+12*gm(1,3)*dgm(1,2)+12*gm(1,2)*dgm(1,3)&
1015 &   -4.5d0*gm(1,1)*dgm(2,3))
1016    cm(3,7,12)=gm(2,2)*(6*gm(2,3)*dgm(2,2)+3*gm(2,2)*dgm(2,3))
1017    cm(3,8,12)=4.5d0*gm(2,3)**2*dgm(2,3)+7.5d0*gm(2,2)*gm(3,3)*dgm(2,3)&
1018 &   +7.5d0*gm(2,3)*(gm(3,3)*dgm(2,2)+gm(2,2)*dgm(3,3))
1019    cm(3,9,12)=3*gm(2,3)**2*dgm(2,2)+6*gm(2,2)*gm(2,3)*dgm(2,3)+gm(2,2)&
1020 &   *(12*gm(3,3)*dgm(2,2)+6*gm(2,2)*dgm(3,3))
1021    cm(3,10,12)=-1.5d0*gm(3,3)**2*dgm(2,2)+4.5d0*gm(2,3)**2*dgm(3,3)&
1022 &   +gm(3,3)*(9*gm(2,3)*dgm(2,3)-3*gm(2,2)*dgm(3,3))
1023    cm(3,1,13)=7.5d0*gm(1,2)**2*dgm(1,2)-1.5d0*gm(1,1)*gm(2,2)*dgm(1,2)&
1024 &   -1.5d0*gm(1,2)*(gm(2,2)*dgm(1,1)+gm(1,1)*dgm(2,2))
1025    cm(3,2,13)=gm(2,2)*(3*gm(2,2)*dgm(1,2)+6*gm(1,2)*dgm(2,2))
1026    cm(3,3,13)=7.5d0*gm(2,3)**2*dgm(1,2)-1.5d0*gm(1,2)*gm(3,3)*dgm(2,2)&
1027 &   +gm(2,3)*(-3*gm(2,2)*dgm(1,3)-3*gm(1,3)*dgm(2,2)+15*gm(1,2)*dgm(2,3))&
1028 &   +gm(2,2)*(-1.5d0*gm(3,3)*dgm(1,2)-3*gm(1,3)*dgm(2,3)-1.5d0*gm(1,2)&
1029 &   *dgm(3,3))
1030    cm(3,4,13)=-3*gm(2,2)**2*dgm(1,3)+9*gm(1,2)*gm(2,3)*dgm(2,2)+gm(2,2)&
1031 &   *(9*gm(2,3)*dgm(1,2)-6*gm(1,3)*dgm(2,2)+9*gm(1,2)*dgm(2,3))
1032    cm(3,5,13)=-1.5d0*gm(1,1)*gm(2,3)*dgm(2,2)+gm(1,2)*(15*gm(2,3)&
1033 &   *dgm(1,2)-3*gm(1,3)*dgm(2,2))+7.5d0*gm(1,2)**2*dgm(2,3)+gm(2,2)&
1034 &   *(-1.5d0*gm(2,3)*dgm(1,1)-3*gm(1,3)*dgm(1,2)-3*gm(1,2)*dgm(1,3)&
1035 &   -1.5d0*gm(1,1)*dgm(2,3))
1036    cm(3,6,13)=-1.5d0*gm(2,2)**2*dgm(1,1)+4.5d0*gm(1,2)**2*dgm(2,2)&
1037 &   +gm(2,2)*(9*gm(1,2)*dgm(1,2)-3*gm(1,1)*dgm(2,2))
1038    cm(3,7,13)=3*gm(2,2)**2*dgm(2,2)
1039    cm(3,8,13)=4.5d0*gm(2,3)**2*dgm(2,2)+9*gm(2,2)*gm(2,3)*dgm(2,3)&
1040 &   +gm(2,2)*(-3*gm(3,3)*dgm(2,2)-1.5d0*gm(2,2)*dgm(3,3))
1041    cm(3,9,13)=gm(2,2)*(6*gm(2,3)*dgm(2,2)+3*gm(2,2)*dgm(2,3))
1042    cm(3,10,13)=7.5d0*gm(2,3)**2*dgm(2,3)-1.5d0*gm(2,2)*gm(3,3)*dgm(2,3)&
1043 &   -1.5d0*gm(2,3)*(gm(3,3)*dgm(2,2)+gm(2,2)*dgm(3,3))
1044    cm(3,1,14)=7.5d0*gm(1,3)**2*dgm(1,2)+gm(1,1)*(-1.5d0*gm(3,3)*dgm(1,2)&
1045 &   -3*gm(2,3)*dgm(1,3))+gm(1,3)*(-3*gm(2,3)*dgm(1,1)+15*gm(1,2)&
1046 &   *dgm(1,3)-3*gm(1,1)*dgm(2,3))-1.5d0*gm(1,2)*(gm(3,3)*dgm(1,1)&
1047 &   +gm(1,1)*dgm(3,3))
1048    cm(3,2,14)=1.5d0*gm(2,3)**2*dgm(1,2)-4.5d0*gm(1,2)*gm(3,3)*dgm(2,2)&
1049 &   +gm(2,3)*(12*gm(2,2)*dgm(1,3)+12*gm(1,3)*dgm(2,2)+3*gm(1,2)*dgm(2,3))&
1050 &   +gm(2,2)*(-4.5d0*gm(3,3)*dgm(1,2)+12*gm(1,3)*dgm(2,3)-4.5d0*gm(1,2)&
1051 &   *dgm(3,3))
1052    cm(3,3,14)=6*gm(3,3)**2*dgm(1,2)+3*gm(1,3)*gm(2,3)*dgm(3,3)+gm(3,3)&
1053 &   *(3*gm(2,3)*dgm(1,3)+3*gm(1,3)*dgm(2,3)+12*gm(1,2)*dgm(3,3))
1054    cm(3,4,14)=3*gm(2,3)**2*dgm(1,3)+gm(3,3)*(12*gm(1,3)*dgm(2,2)&
1055 &   +3*gm(1,2)*dgm(2,3))+gm(2,3)*(3*gm(3,3)*dgm(1,2)+6*gm(1,3)*dgm(2,3)&
1056 &   +3*gm(1,2)*dgm(3,3))+12*gm(2,2)*(gm(3,3)*dgm(1,3)+gm(1,3)*dgm(3,3))
1057    cm(3,5,14)=1.5d0*gm(1,3)**2*dgm(2,3)+gm(3,3)*(12*gm(1,2)*dgm(1,3)&
1058 &   -4.5d0*gm(1,1)*dgm(2,3))+gm(2,3)*(-4.5d0*gm(3,3)*dgm(1,1)+3*gm(1,3)&
1059 &   *dgm(1,3)-4.5d0*gm(1,1)*dgm(3,3))+12*gm(1,3)*(gm(3,3)*dgm(1,2)&
1060 &   +gm(1,2)*dgm(3,3))
1061    cm(3,6,14)=-3*gm(2,3)**2*dgm(1,1)-6*gm(1,2)*gm(3,3)*dgm(1,2)+7.5d0*gm(1,3)&
1062 &   **2*dgm(2,2)-1.5d0*gm(1,1)*gm(3,3)*dgm(2,2)+9*gm(1,2)*gm(1,3)&
1063 &   *dgm(2,3)+gm(2,3)*(9*gm(1,3)*dgm(1,2)+9*gm(1,2)*dgm(1,3)-6*gm(1,1)&
1064 &   *dgm(2,3))-3*gm(1,2)**2*dgm(3,3)+gm(2,2)*(-1.5d0*gm(3,3)*dgm(1,1)&
1065 &   +15*gm(1,3)*dgm(1,3)-1.5d0*gm(1,1)*dgm(3,3))
1066    cm(3,7,14)=4.5d0*gm(2,3)**2*dgm(2,2)+9*gm(2,2)*gm(2,3)*dgm(2,3)&
1067 &   +gm(2,2)*(-3*gm(3,3)*dgm(2,2)-1.5d0*gm(2,2)*dgm(3,3))
1068    cm(3,8,14)=6*gm(3,3)**2*dgm(2,2)+3*gm(2,3)**2*dgm(3,3)+gm(3,3)&
1069 &   *(6*gm(2,3)*dgm(2,3)+12*gm(2,2)*dgm(3,3))
1070    cm(3,9,14)=4.5d0*gm(2,3)**2*dgm(2,3)+7.5d0*gm(2,2)*gm(3,3)*dgm(2,3)&
1071 &   +7.5d0*gm(2,3)*(gm(3,3)*dgm(2,2)+gm(2,2)*dgm(3,3))
1072    cm(3,10,14)=gm(3,3)*(3*gm(3,3)*dgm(2,3)+6*gm(2,3)*dgm(3,3))
1073    cm(3,1,15)=7.5d0*gm(1,3)**2*dgm(1,3)-1.5d0*gm(1,1)*gm(3,3)*dgm(1,3)&
1074 &   -1.5d0*gm(1,3)*(gm(3,3)*dgm(1,1)+gm(1,1)*dgm(3,3))
1075    cm(3,2,15)=7.5d0*gm(2,3)**2*dgm(1,3)+gm(3,3)*(-1.5d0*gm(1,3)*dgm(2,2)&
1076 &   -3*gm(1,2)*dgm(2,3))+gm(2,3)*(-3*gm(3,3)*dgm(1,2)+15*gm(1,3)&
1077 &   *dgm(2,3)-3*gm(1,2)*dgm(3,3))-1.5d0*gm(2,2)*(gm(3,3)*dgm(1,3)&
1078 &   +gm(1,3)*dgm(3,3))
1079    cm(3,3,15)=gm(3,3)*(3*gm(3,3)*dgm(1,3)+6*gm(1,3)*dgm(3,3))
1080    cm(3,4,15)=-3*gm(3,3)**2*dgm(1,2)+9*gm(1,3)*gm(2,3)*dgm(3,3)+gm(3,3)&
1081 &   *(9*gm(2,3)*dgm(1,3)+9*gm(1,3)*dgm(2,3)-6*gm(1,2)*dgm(3,3))
1082    cm(3,5,15)=-1.5d0*gm(3,3)**2*dgm(1,1)+4.5d0*gm(1,3)**2*dgm(3,3)&
1083 &   +gm(3,3)*(9*gm(1,3)*dgm(1,3)-3*gm(1,1)*dgm(3,3))
1084    cm(3,6,15)=7.5d0*gm(1,3)**2*dgm(2,3)+gm(3,3)*(-3*gm(1,2)*dgm(1,3)&
1085 &   -1.5d0*gm(1,1)*dgm(2,3))+gm(2,3)*(-1.5d0*gm(3,3)*dgm(1,1)+15*gm(1,3)&
1086 &   *dgm(1,3)-1.5d0*gm(1,1)*dgm(3,3))-3*gm(1,3)*(gm(3,3)*dgm(1,2)&
1087 &   +gm(1,2)*dgm(3,3))
1088    cm(3,7,15)=7.5d0*gm(2,3)**2*dgm(2,3)-1.5d0*gm(2,2)*gm(3,3)*dgm(2,3)&
1089 &   -1.5d0*gm(2,3)*(gm(3,3)*dgm(2,2)+gm(2,2)*dgm(3,3))
1090    cm(3,8,15)=gm(3,3)*(3*gm(3,3)*dgm(2,3)+6*gm(2,3)*dgm(3,3))
1091    cm(3,9,15)=-1.5d0*gm(3,3)**2*dgm(2,2)+4.5d0*gm(2,3)**2*dgm(3,3)&
1092 &   +gm(3,3)*(9*gm(2,3)*dgm(2,3)-3*gm(2,2)*dgm(3,3))
1093    cm(3,10,15)=3*gm(3,3)**2*dgm(3,3)
1094  end if
1095 !
1096 !Section for volume derivative term for diagonal strains.  This is
1097 !equivalent to calculating the reduced gradients but it is easier
1098 !to embed it here as a separate calculation than cope with the logic
1099 !and ffkg index rearrangements necessary to do a separate force
1100 !calculation.
1101 !
1102  if(istr<=3) then
1103    if(rank==0)then
1104      cm(1,1,1)=cm(1,1,1)-(1)
1105      cm(2,1,2)=cm(2,1,2)-(1)
1106      cm(3,1,3)=cm(3,1,3)-(1)
1107    elseif(rank==1)then
1108      cm(1,1,1)=cm(1,1,1)-(gm(1,1))
1109      cm(1,2,1)=cm(1,2,1)-(gm(1,2))
1110      cm(1,3,1)=cm(1,3,1)-(gm(1,3))
1111      cm(1,1,5)=cm(1,1,5)-(gm(1,3))
1112      cm(1,2,5)=cm(1,2,5)-(gm(2,3))
1113      cm(1,3,5)=cm(1,3,5)-(gm(3,3))
1114      cm(1,1,6)=cm(1,1,6)-(gm(1,2))
1115      cm(1,2,6)=cm(1,2,6)-(gm(2,2))
1116      cm(1,3,6)=cm(1,3,6)-(gm(2,3))
1117      cm(2,1,2)=cm(2,1,2)-(gm(1,2))
1118      cm(2,2,2)=cm(2,2,2)-(gm(2,2))
1119      cm(2,3,2)=cm(2,3,2)-(gm(2,3))
1120      cm(2,1,4)=cm(2,1,4)-(gm(1,3))
1121      cm(2,2,4)=cm(2,2,4)-(gm(2,3))
1122      cm(2,3,4)=cm(2,3,4)-(gm(3,3))
1123      cm(2,1,6)=cm(2,1,6)-(gm(1,1))
1124      cm(2,2,6)=cm(2,2,6)-(gm(1,2))
1125      cm(2,3,6)=cm(2,3,6)-(gm(1,3))
1126      cm(3,1,3)=cm(3,1,3)-(gm(1,3))
1127      cm(3,2,3)=cm(3,2,3)-(gm(2,3))
1128      cm(3,3,3)=cm(3,3,3)-(gm(3,3))
1129      cm(3,1,4)=cm(3,1,4)-(gm(1,2))
1130      cm(3,2,4)=cm(3,2,4)-(gm(2,2))
1131      cm(3,3,4)=cm(3,3,4)-(gm(2,3))
1132      cm(3,1,5)=cm(3,1,5)-(gm(1,1))
1133      cm(3,2,5)=cm(3,2,5)-(gm(1,2))
1134      cm(3,3,5)=cm(3,3,5)-(gm(1,3))
1135    elseif(rank==2)then
1136      cm(1,1,1)=cm(1,1,1)-(gm(1,1)**2)
1137      cm(1,2,1)=cm(1,2,1)-(1.5d0*gm(1,2)**2-0.5d0*gm(1,1)*gm(2,2))
1138      cm(1,3,1)=cm(1,3,1)-(1.5d0*gm(1,3)**2-0.5d0*gm(1,1)*gm(3,3))
1139      cm(1,4,1)=cm(1,4,1)-(3*gm(1,2)*gm(1,3)-gm(1,1)*gm(2,3))
1140      cm(1,5,1)=cm(1,5,1)-(2*gm(1,1)*gm(1,3))
1141      cm(1,6,1)=cm(1,6,1)-(2*gm(1,1)*gm(1,2))
1142      cm(1,1,2)=cm(1,1,2)-(1.5d0*gm(1,2)**2-0.5d0*gm(1,1)*gm(2,2))
1143      cm(1,2,2)=cm(1,2,2)-(gm(2,2)**2)
1144      cm(1,3,2)=cm(1,3,2)-(1.5d0*gm(2,3)**2-0.5d0*gm(2,2)*gm(3,3))
1145      cm(1,4,2)=cm(1,4,2)-(2*gm(2,2)*gm(2,3))
1146      cm(1,5,2)=cm(1,5,2)-(-gm(1,3)*gm(2,2)+3*gm(1,2)*gm(2,3))
1147      cm(1,6,2)=cm(1,6,2)-(2*gm(1,2)*gm(2,2))
1148      cm(1,1,3)=cm(1,1,3)-(1.5d0*gm(1,3)**2-0.5d0*gm(1,1)*gm(3,3))
1149      cm(1,2,3)=cm(1,2,3)-(1.5d0*gm(2,3)**2-0.5d0*gm(2,2)*gm(3,3))
1150      cm(1,3,3)=cm(1,3,3)-(gm(3,3)**2)
1151      cm(1,4,3)=cm(1,4,3)-(2*gm(2,3)*gm(3,3))
1152      cm(1,5,3)=cm(1,5,3)-(2*gm(1,3)*gm(3,3))
1153      cm(1,6,3)=cm(1,6,3)-(3*gm(1,3)*gm(2,3)-gm(1,2)*gm(3,3))
1154      cm(1,1,4)=cm(1,1,4)-(3*gm(1,2)*gm(1,3)-gm(1,1)*gm(2,3))
1155      cm(1,2,4)=cm(1,2,4)-(2*gm(2,2)*gm(2,3))
1156      cm(1,3,4)=cm(1,3,4)-(2*gm(2,3)*gm(3,3))
1157      cm(1,4,4)=cm(1,4,4)-(gm(2,3)**2+3*gm(2,2)*gm(3,3))
1158      cm(1,5,4)=cm(1,5,4)-(gm(1,3)*gm(2,3)+3*gm(1,2)*gm(3,3))
1159      cm(1,6,4)=cm(1,6,4)-(3*gm(1,3)*gm(2,2)+gm(1,2)*gm(2,3))
1160      cm(1,1,5)=cm(1,1,5)-(2*gm(1,1)*gm(1,3))
1161      cm(1,2,5)=cm(1,2,5)-(-gm(1,3)*gm(2,2)+3*gm(1,2)*gm(2,3))
1162      cm(1,3,5)=cm(1,3,5)-(2*gm(1,3)*gm(3,3))
1163      cm(1,4,5)=cm(1,4,5)-(gm(1,3)*gm(2,3)+3*gm(1,2)*gm(3,3))
1164      cm(1,5,5)=cm(1,5,5)-(gm(1,3)**2+3*gm(1,1)*gm(3,3))
1165      cm(1,6,5)=cm(1,6,5)-(gm(1,2)*gm(1,3)+3*gm(1,1)*gm(2,3))
1166      cm(1,1,6)=cm(1,1,6)-(2*gm(1,1)*gm(1,2))
1167      cm(1,2,6)=cm(1,2,6)-(2*gm(1,2)*gm(2,2))
1168      cm(1,3,6)=cm(1,3,6)-(3*gm(1,3)*gm(2,3)-gm(1,2)*gm(3,3))
1169      cm(1,4,6)=cm(1,4,6)-(3*gm(1,3)*gm(2,2)+gm(1,2)*gm(2,3))
1170      cm(1,5,6)=cm(1,5,6)-(gm(1,2)*gm(1,3)+3*gm(1,1)*gm(2,3))
1171      cm(1,6,6)=cm(1,6,6)-(gm(1,2)**2+3*gm(1,1)*gm(2,2))
1172      cm(2,1,2)=cm(2,1,2)-(2*gm(1,1)*gm(1,2))
1173      cm(2,2,2)=cm(2,2,2)-(2*gm(1,2)*gm(2,2))
1174      cm(2,3,2)=cm(2,3,2)-(3*gm(1,3)*gm(2,3)-gm(1,2)*gm(3,3))
1175      cm(2,4,2)=cm(2,4,2)-(3*gm(1,3)*gm(2,2)+gm(1,2)*gm(2,3))
1176      cm(2,5,2)=cm(2,5,2)-(gm(1,2)*gm(1,3)+3*gm(1,1)*gm(2,3))
1177      cm(2,6,2)=cm(2,6,2)-(gm(1,2)**2+3*gm(1,1)*gm(2,2))
1178      cm(2,1,4)=cm(2,1,4)-(2*gm(1,1)*gm(1,3))
1179      cm(2,2,4)=cm(2,2,4)-(-gm(1,3)*gm(2,2)+3*gm(1,2)*gm(2,3))
1180      cm(2,3,4)=cm(2,3,4)-(2*gm(1,3)*gm(3,3))
1181      cm(2,4,4)=cm(2,4,4)-(gm(1,3)*gm(2,3)+3*gm(1,2)*gm(3,3))
1182      cm(2,5,4)=cm(2,5,4)-(gm(1,3)**2+3*gm(1,1)*gm(3,3))
1183      cm(2,6,4)=cm(2,6,4)-(gm(1,2)*gm(1,3)+3*gm(1,1)*gm(2,3))
1184      cm(2,1,6)=cm(2,1,6)-(gm(1,1)**2)
1185      cm(2,2,6)=cm(2,2,6)-(1.5d0*gm(1,2)**2-0.5d0*gm(1,1)*gm(2,2))
1186      cm(2,3,6)=cm(2,3,6)-(1.5d0*gm(1,3)**2-0.5d0*gm(1,1)*gm(3,3))
1187      cm(2,4,6)=cm(2,4,6)-(3*gm(1,2)*gm(1,3)-gm(1,1)*gm(2,3))
1188      cm(2,5,6)=cm(2,5,6)-(2*gm(1,1)*gm(1,3))
1189      cm(2,6,6)=cm(2,6,6)-(2*gm(1,1)*gm(1,2))
1190      cm(2,1,7)=cm(2,1,7)-(1.5d0*gm(1,2)**2-0.5d0*gm(1,1)*gm(2,2))
1191      cm(2,2,7)=cm(2,2,7)-(gm(2,2)**2)
1192      cm(2,3,7)=cm(2,3,7)-(1.5d0*gm(2,3)**2-0.5d0*gm(2,2)*gm(3,3))
1193      cm(2,4,7)=cm(2,4,7)-(2*gm(2,2)*gm(2,3))
1194      cm(2,5,7)=cm(2,5,7)-(-gm(1,3)*gm(2,2)+3*gm(1,2)*gm(2,3))
1195      cm(2,6,7)=cm(2,6,7)-(2*gm(1,2)*gm(2,2))
1196      cm(2,1,8)=cm(2,1,8)-(1.5d0*gm(1,3)**2-0.5d0*gm(1,1)*gm(3,3))
1197      cm(2,2,8)=cm(2,2,8)-(1.5d0*gm(2,3)**2-0.5d0*gm(2,2)*gm(3,3))
1198      cm(2,3,8)=cm(2,3,8)-(gm(3,3)**2)
1199      cm(2,4,8)=cm(2,4,8)-(2*gm(2,3)*gm(3,3))
1200      cm(2,5,8)=cm(2,5,8)-(2*gm(1,3)*gm(3,3))
1201      cm(2,6,8)=cm(2,6,8)-(3*gm(1,3)*gm(2,3)-gm(1,2)*gm(3,3))
1202      cm(2,1,9)=cm(2,1,9)-(3*gm(1,2)*gm(1,3)-gm(1,1)*gm(2,3))
1203      cm(2,2,9)=cm(2,2,9)-(2*gm(2,2)*gm(2,3))
1204      cm(2,3,9)=cm(2,3,9)-(2*gm(2,3)*gm(3,3))
1205      cm(2,4,9)=cm(2,4,9)-(gm(2,3)**2+3*gm(2,2)*gm(3,3))
1206      cm(2,5,9)=cm(2,5,9)-(gm(1,3)*gm(2,3)+3*gm(1,2)*gm(3,3))
1207      cm(2,6,9)=cm(2,6,9)-(3*gm(1,3)*gm(2,2)+gm(1,2)*gm(2,3))
1208      cm(3,1,3)=cm(3,1,3)-(2*gm(1,1)*gm(1,3))
1209      cm(3,2,3)=cm(3,2,3)-(-gm(1,3)*gm(2,2)+3*gm(1,2)*gm(2,3))
1210      cm(3,3,3)=cm(3,3,3)-(2*gm(1,3)*gm(3,3))
1211      cm(3,4,3)=cm(3,4,3)-(gm(1,3)*gm(2,3)+3*gm(1,2)*gm(3,3))
1212      cm(3,5,3)=cm(3,5,3)-(gm(1,3)**2+3*gm(1,1)*gm(3,3))
1213      cm(3,6,3)=cm(3,6,3)-(gm(1,2)*gm(1,3)+3*gm(1,1)*gm(2,3))
1214      cm(3,1,4)=cm(3,1,4)-(2*gm(1,1)*gm(1,2))
1215      cm(3,2,4)=cm(3,2,4)-(2*gm(1,2)*gm(2,2))
1216      cm(3,3,4)=cm(3,3,4)-(3*gm(1,3)*gm(2,3)-gm(1,2)*gm(3,3))
1217      cm(3,4,4)=cm(3,4,4)-(3*gm(1,3)*gm(2,2)+gm(1,2)*gm(2,3))
1218      cm(3,5,4)=cm(3,5,4)-(gm(1,2)*gm(1,3)+3*gm(1,1)*gm(2,3))
1219      cm(3,6,4)=cm(3,6,4)-(gm(1,2)**2+3*gm(1,1)*gm(2,2))
1220      cm(3,1,5)=cm(3,1,5)-(gm(1,1)**2)
1221      cm(3,2,5)=cm(3,2,5)-(1.5d0*gm(1,2)**2-0.5d0*gm(1,1)*gm(2,2))
1222      cm(3,3,5)=cm(3,3,5)-(1.5d0*gm(1,3)**2-0.5d0*gm(1,1)*gm(3,3))
1223      cm(3,4,5)=cm(3,4,5)-(3*gm(1,2)*gm(1,3)-gm(1,1)*gm(2,3))
1224      cm(3,5,5)=cm(3,5,5)-(2*gm(1,1)*gm(1,3))
1225      cm(3,6,5)=cm(3,6,5)-(2*gm(1,1)*gm(1,2))
1226      cm(3,1,8)=cm(3,1,8)-(3*gm(1,2)*gm(1,3)-gm(1,1)*gm(2,3))
1227      cm(3,2,8)=cm(3,2,8)-(2*gm(2,2)*gm(2,3))
1228      cm(3,3,8)=cm(3,3,8)-(2*gm(2,3)*gm(3,3))
1229      cm(3,4,8)=cm(3,4,8)-(gm(2,3)**2+3*gm(2,2)*gm(3,3))
1230      cm(3,5,8)=cm(3,5,8)-(gm(1,3)*gm(2,3)+3*gm(1,2)*gm(3,3))
1231      cm(3,6,8)=cm(3,6,8)-(3*gm(1,3)*gm(2,2)+gm(1,2)*gm(2,3))
1232      cm(3,1,9)=cm(3,1,9)-(1.5d0*gm(1,2)**2-0.5d0*gm(1,1)*gm(2,2))
1233      cm(3,2,9)=cm(3,2,9)-(gm(2,2)**2)
1234      cm(3,3,9)=cm(3,3,9)-(1.5d0*gm(2,3)**2-0.5d0*gm(2,2)*gm(3,3))
1235      cm(3,4,9)=cm(3,4,9)-(2*gm(2,2)*gm(2,3))
1236      cm(3,5,9)=cm(3,5,9)-(-gm(1,3)*gm(2,2)+3*gm(1,2)*gm(2,3))
1237      cm(3,6,9)=cm(3,6,9)-(2*gm(1,2)*gm(2,2))
1238      cm(3,1,10)=cm(3,1,10)-(1.5d0*gm(1,3)**2-0.5d0*gm(1,1)*gm(3,3))
1239      cm(3,2,10)=cm(3,2,10)-(1.5d0*gm(2,3)**2-0.5d0*gm(2,2)*gm(3,3))
1240      cm(3,3,10)=cm(3,3,10)-(gm(3,3)**2)
1241      cm(3,4,10)=cm(3,4,10)-(2*gm(2,3)*gm(3,3))
1242      cm(3,5,10)=cm(3,5,10)-(2*gm(1,3)*gm(3,3))
1243      cm(3,6,10)=cm(3,6,10)-(3*gm(1,3)*gm(2,3)-gm(1,2)*gm(3,3))
1244    elseif(rank==3)then
1245      cm(1,1,1)=cm(1,1,1)-(gm(1,1)**3)
1246      cm(1,2,1)=cm(1,2,1)-(gm(1,1)*(4.5d0*gm(1,2)**2-1.5d0*gm(1,1)*gm(2,2)))
1247      cm(1,3,1)=cm(1,3,1)-(gm(1,1)*(4.5d0*gm(1,3)**2-1.5d0*gm(1,1)*gm(3,3)))
1248      cm(1,4,1)=cm(1,4,1)-(gm(1,1)*(9*gm(1,2)*gm(1,3)-3*gm(1,1)*gm(2,3)))
1249      cm(1,5,1)=cm(1,5,1)-(3*gm(1,1)**2*gm(1,3))
1250      cm(1,6,1)=cm(1,6,1)-(3*gm(1,1)**2*gm(1,2))
1251      cm(1,7,1)=cm(1,7,1)-(2.5d0*gm(1,2)**3-1.5d0*gm(1,1)*gm(1,2)*gm(2,2))
1252      cm(1,8,1)=cm(1,8,1)-(-3*gm(1,1)*gm(1,3)*gm(2,3)+gm(1,2)*(7.5d0*gm(1,3)&
1253 &     **2-1.5d0*gm(1,1)*gm(3,3)))
1254      cm(1,9,1)=cm(1,9,1)-(7.5d0*gm(1,2)**2*gm(1,3)-1.5d0*gm(1,1)*gm(1,3)&
1255 &     *gm(2,2)-3*gm(1,1)*gm(1,2)*gm(2,3))
1256      cm(1,10,1)=cm(1,10,1)-(2.5d0*gm(1,3)**3-1.5d0*gm(1,1)*gm(1,3)&
1257 &     *gm(3,3))
1258      cm(1,1,2)=cm(1,1,2)-(gm(1,1)*(4.5d0*gm(1,2)**2-1.5d0*gm(1,1)*gm(2,2)))
1259      cm(1,2,2)=cm(1,2,2)-(gm(2,2)*(3*gm(1,2)**2+6*gm(1,1)*gm(2,2)))
1260      cm(1,3,2)=cm(1,3,2)-(-3*gm(1,3)**2*gm(2,2)+9*gm(1,2)*gm(1,3)*gm(2,3)&
1261 &     -3*gm(1,2)**2*gm(3,3)+gm(1,1)*(7.5d0*gm(2,3)**2-1.5d0*gm(2,2)&
1262 &     *gm(3,3)))
1263      cm(1,4,2)=cm(1,4,2)-(3*gm(1,2)*gm(1,3)*gm(2,2)+3*gm(1,2)**2*gm(2,3)&
1264 &     +12*gm(1,1)*gm(2,2)*gm(2,3))
1265      cm(1,5,2)=cm(1,5,2)-(1.5d0*gm(1,2)**2*gm(1,3)-4.5d0*gm(1,1)*gm(1,3)&
1266 &     *gm(2,2)+12*gm(1,1)*gm(1,2)*gm(2,3))
1267      cm(1,6,2)=cm(1,6,2)-(1.5d0*gm(1,2)**3+7.5d0*gm(1,1)*gm(1,2)*gm(2,2))
1268      cm(1,7,2)=cm(1,7,2)-(3*gm(1,2)*gm(2,2)**2)
1269      cm(1,8,2)=cm(1,8,2)-(12*gm(1,3)*gm(2,2)*gm(2,3)+gm(1,2)*(1.5d0*gm(2,3)&
1270 &     **2-4.5d0*gm(2,2)*gm(3,3)))
1271      cm(1,9,2)=cm(1,9,2)-(gm(2,2)*(6*gm(1,3)*gm(2,2)+3*gm(1,2)*gm(2,3)))
1272      cm(1,10,2)=cm(1,10,2)-(-3*gm(1,2)*gm(2,3)*gm(3,3)+gm(1,3)*(7.5d0*gm(2,3)&
1273 &     **2-1.5d0*gm(2,2)*gm(3,3)))
1274      cm(1,1,3)=cm(1,1,3)-(gm(1,1)*(4.5d0*gm(1,3)**2-1.5d0*gm(1,1)*gm(3,3)))
1275      cm(1,2,3)=cm(1,2,3)-(-3*gm(1,3)**2*gm(2,2)+9*gm(1,2)*gm(1,3)*gm(2,3)&
1276 &     -3*gm(1,2)**2*gm(3,3)+gm(1,1)*(7.5d0*gm(2,3)**2-1.5d0*gm(2,2)&
1277 &     *gm(3,3)))
1278      cm(1,3,3)=cm(1,3,3)-(gm(3,3)*(3*gm(1,3)**2+6*gm(1,1)*gm(3,3)))
1279      cm(1,4,3)=cm(1,4,3)-(3*gm(1,3)**2*gm(2,3)+3*gm(1,2)*gm(1,3)*gm(3,3)&
1280 &     +12*gm(1,1)*gm(2,3)*gm(3,3))
1281      cm(1,5,3)=cm(1,5,3)-(1.5d0*gm(1,3)**3+7.5d0*gm(1,1)*gm(1,3)*gm(3,3))
1282      cm(1,6,3)=cm(1,6,3)-(12*gm(1,1)*gm(1,3)*gm(2,3)+gm(1,2)*(1.5d0*gm(1,3)&
1283 &     **2-4.5d0*gm(1,1)*gm(3,3)))
1284      cm(1,7,3)=cm(1,7,3)-(-3*gm(1,3)*gm(2,2)*gm(2,3)+gm(1,2)*(7.5d0*gm(2,3)&
1285 &     **2-1.5d0*gm(2,2)*gm(3,3)))
1286      cm(1,8,3)=cm(1,8,3)-(gm(3,3)*(3*gm(1,3)*gm(2,3)+6*gm(1,2)*gm(3,3)))
1287      cm(1,9,3)=cm(1,9,3)-(12*gm(1,2)*gm(2,3)*gm(3,3)+gm(1,3)*(1.5d0*gm(2,3)&
1288 &     **2-4.5d0*gm(2,2)*gm(3,3)))
1289      cm(1,10,3)=cm(1,10,3)-(3*gm(1,3)*gm(3,3)**2)
1290      cm(1,1,4)=cm(1,1,4)-(gm(1,1)*(9*gm(1,2)*gm(1,3)-3*gm(1,1)*gm(2,3)))
1291      cm(1,2,4)=cm(1,2,4)-(3*gm(1,2)*gm(1,3)*gm(2,2)+3*gm(1,2)**2*gm(2,3)&
1292 &     +12*gm(1,1)*gm(2,2)*gm(2,3))
1293      cm(1,3,4)=cm(1,3,4)-(3*gm(1,3)**2*gm(2,3)+3*gm(1,2)*gm(1,3)*gm(3,3)&
1294 &     +12*gm(1,1)*gm(2,3)*gm(3,3))
1295      cm(1,4,4)=cm(1,4,4)-(9*gm(1,3)**2*gm(2,2)-6*gm(1,2)*gm(1,3)*gm(2,3)&
1296 &     +9*gm(1,2)**2*gm(3,3)+gm(1,1)*(9*gm(2,3)**2+15*gm(2,2)*gm(3,3)))
1297      cm(1,5,4)=cm(1,5,4)-(3*gm(1,1)*gm(1,3)*gm(2,3)+gm(1,2)*(3*gm(1,3)&
1298 &     **2+12*gm(1,1)*gm(3,3)))
1299      cm(1,6,4)=cm(1,6,4)-(3*gm(1,2)**2*gm(1,3)+12*gm(1,1)*gm(1,3)*gm(2,2)&
1300 &     +3*gm(1,1)*gm(1,2)*gm(2,3))
1301      cm(1,7,4)=cm(1,7,4)-(gm(2,2)*(-3*gm(1,3)*gm(2,2)+9*gm(1,2)*gm(2,3)))
1302      cm(1,8,4)=cm(1,8,4)-(3*gm(1,2)*gm(2,3)*gm(3,3)+gm(1,3)*(3*gm(2,3)&
1303 &     **2+12*gm(2,2)*gm(3,3)))
1304      cm(1,9,4)=cm(1,9,4)-(3*gm(1,3)*gm(2,2)*gm(2,3)+gm(1,2)*(3*gm(2,3)&
1305 &     **2+12*gm(2,2)*gm(3,3)))
1306      cm(1,10,4)=cm(1,10,4)-(gm(3,3)*(9*gm(1,3)*gm(2,3)-3*gm(1,2)*gm(3,3)))
1307      cm(1,1,5)=cm(1,1,5)-(3*gm(1,1)**2*gm(1,3))
1308      cm(1,2,5)=cm(1,2,5)-(1.5d0*gm(1,2)**2*gm(1,3)-4.5d0*gm(1,1)*gm(1,3)&
1309 &     *gm(2,2)+12*gm(1,1)*gm(1,2)*gm(2,3))
1310      cm(1,3,5)=cm(1,3,5)-(1.5d0*gm(1,3)**3+7.5d0*gm(1,1)*gm(1,3)*gm(3,3))
1311      cm(1,4,5)=cm(1,4,5)-(3*gm(1,1)*gm(1,3)*gm(2,3)+gm(1,2)*(3*gm(1,3)&
1312 &     **2+12*gm(1,1)*gm(3,3)))
1313      cm(1,5,5)=cm(1,5,5)-(gm(1,1)*(3*gm(1,3)**2+6*gm(1,1)*gm(3,3)))
1314      cm(1,6,5)=cm(1,6,5)-(gm(1,1)*(3*gm(1,2)*gm(1,3)+6*gm(1,1)*gm(2,3)))
1315      cm(1,7,5)=cm(1,7,5)-(-3*gm(1,2)*gm(1,3)*gm(2,2)+7.5d0*gm(1,2)&
1316 &     **2*gm(2,3)-1.5d0*gm(1,1)*gm(2,2)*gm(2,3))
1317      cm(1,8,5)=cm(1,8,5)-(1.5d0*gm(1,3)**2*gm(2,3)+12*gm(1,2)*gm(1,3)&
1318 &     *gm(3,3)-4.5d0*gm(1,1)*gm(2,3)*gm(3,3))
1319      cm(1,9,5)=cm(1,9,5)-(-3*gm(1,3)**2*gm(2,2)+9*gm(1,2)*gm(1,3)*gm(2,3)&
1320 &     +7.5d0*gm(1,2)**2*gm(3,3)+gm(1,1)*(-3*gm(2,3)**2-1.5d0*gm(2,2)&
1321 &     *gm(3,3)))
1322      cm(1,10,5)=cm(1,10,5)-(gm(3,3)*(4.5d0*gm(1,3)**2-1.5d0*gm(1,1)&
1323 &     *gm(3,3)))
1324      cm(1,1,6)=cm(1,1,6)-(3*gm(1,1)**2*gm(1,2))
1325      cm(1,2,6)=cm(1,2,6)-(1.5d0*gm(1,2)**3+7.5d0*gm(1,1)*gm(1,2)*gm(2,2))
1326      cm(1,3,6)=cm(1,3,6)-(12*gm(1,1)*gm(1,3)*gm(2,3)+gm(1,2)*(1.5d0*gm(1,3)&
1327 &     **2-4.5d0*gm(1,1)*gm(3,3)))
1328      cm(1,4,6)=cm(1,4,6)-(3*gm(1,2)**2*gm(1,3)+12*gm(1,1)*gm(1,3)*gm(2,2)&
1329 &     +3*gm(1,1)*gm(1,2)*gm(2,3))
1330      cm(1,5,6)=cm(1,5,6)-(gm(1,1)*(3*gm(1,2)*gm(1,3)+6*gm(1,1)*gm(2,3)))
1331      cm(1,6,6)=cm(1,6,6)-(gm(1,1)*(3*gm(1,2)**2+6*gm(1,1)*gm(2,2)))
1332      cm(1,7,6)=cm(1,7,6)-(gm(2,2)*(4.5d0*gm(1,2)**2-1.5d0*gm(1,1)*gm(2,2)))
1333      cm(1,8,6)=cm(1,8,6)-(7.5d0*gm(1,3)**2*gm(2,2)+9*gm(1,2)*gm(1,3)&
1334 &     *gm(2,3)-3*gm(1,2)**2*gm(3,3)+gm(1,1)*(-3*gm(2,3)**2-1.5d0*gm(2,2)&
1335 &     *gm(3,3)))
1336      cm(1,9,6)=cm(1,9,6)-(12*gm(1,2)*gm(1,3)*gm(2,2)+1.5d0*gm(1,2)&
1337 &     **2*gm(2,3)-4.5d0*gm(1,1)*gm(2,2)*gm(2,3))
1338      cm(1,10,6)=cm(1,10,6)-(7.5d0*gm(1,3)**2*gm(2,3)-3*gm(1,2)*gm(1,3)&
1339 &     *gm(3,3)-1.5d0*gm(1,1)*gm(2,3)*gm(3,3))
1340      cm(1,1,7)=cm(1,1,7)-(2.5d0*gm(1,2)**3-1.5d0*gm(1,1)*gm(1,2)*gm(2,2))
1341      cm(1,2,7)=cm(1,2,7)-(3*gm(1,2)*gm(2,2)**2)
1342      cm(1,3,7)=cm(1,3,7)-(-3*gm(1,3)*gm(2,2)*gm(2,3)+gm(1,2)*(7.5d0*gm(2,3)&
1343 &     **2-1.5d0*gm(2,2)*gm(3,3)))
1344      cm(1,4,7)=cm(1,4,7)-(gm(2,2)*(-3*gm(1,3)*gm(2,2)+9*gm(1,2)*gm(2,3)))
1345      cm(1,5,7)=cm(1,5,7)-(-3*gm(1,2)*gm(1,3)*gm(2,2)+7.5d0*gm(1,2)&
1346 &     **2*gm(2,3)-1.5d0*gm(1,1)*gm(2,2)*gm(2,3))
1347      cm(1,6,7)=cm(1,6,7)-(gm(2,2)*(4.5d0*gm(1,2)**2-1.5d0*gm(1,1)*gm(2,2)))
1348      cm(1,7,7)=cm(1,7,7)-(gm(2,2)**3)
1349      cm(1,8,7)=cm(1,8,7)-(gm(2,2)*(4.5d0*gm(2,3)**2-1.5d0*gm(2,2)*gm(3,3)))
1350      cm(1,9,7)=cm(1,9,7)-(3*gm(2,2)**2*gm(2,3))
1351      cm(1,10,7)=cm(1,10,7)-(2.5d0*gm(2,3)**3-1.5d0*gm(2,2)*gm(2,3)&
1352 &     *gm(3,3))
1353      cm(1,1,8)=cm(1,1,8)-(-3*gm(1,1)*gm(1,3)*gm(2,3)+gm(1,2)*(7.5d0*gm(1,3)&
1354 &     **2-1.5d0*gm(1,1)*gm(3,3)))
1355      cm(1,2,8)=cm(1,2,8)-(12*gm(1,3)*gm(2,2)*gm(2,3)+gm(1,2)*(1.5d0*gm(2,3)&
1356 &     **2-4.5d0*gm(2,2)*gm(3,3)))
1357      cm(1,3,8)=cm(1,3,8)-(gm(3,3)*(3*gm(1,3)*gm(2,3)+6*gm(1,2)*gm(3,3)))
1358      cm(1,4,8)=cm(1,4,8)-(3*gm(1,2)*gm(2,3)*gm(3,3)+gm(1,3)*(3*gm(2,3)&
1359 &     **2+12*gm(2,2)*gm(3,3)))
1360      cm(1,5,8)=cm(1,5,8)-(1.5d0*gm(1,3)**2*gm(2,3)+12*gm(1,2)*gm(1,3)&
1361 &     *gm(3,3)-4.5d0*gm(1,1)*gm(2,3)*gm(3,3))
1362      cm(1,6,8)=cm(1,6,8)-(7.5d0*gm(1,3)**2*gm(2,2)+9*gm(1,2)*gm(1,3)&
1363 &     *gm(2,3)-3*gm(1,2)**2*gm(3,3)+gm(1,1)*(-3*gm(2,3)**2-1.5d0*gm(2,2)&
1364 &     *gm(3,3)))
1365      cm(1,7,8)=cm(1,7,8)-(gm(2,2)*(4.5d0*gm(2,3)**2-1.5d0*gm(2,2)*gm(3,3)))
1366      cm(1,8,8)=cm(1,8,8)-(gm(3,3)*(3*gm(2,3)**2+6*gm(2,2)*gm(3,3)))
1367      cm(1,9,8)=cm(1,9,8)-(1.5d0*gm(2,3)**3+7.5d0*gm(2,2)*gm(2,3)*gm(3,3))
1368      cm(1,10,8)=cm(1,10,8)-(3*gm(2,3)*gm(3,3)**2)
1369      cm(1,1,9)=cm(1,1,9)-(7.5d0*gm(1,2)**2*gm(1,3)-1.5d0*gm(1,1)*gm(1,3)&
1370 &     *gm(2,2)-3*gm(1,1)*gm(1,2)*gm(2,3))
1371      cm(1,2,9)=cm(1,2,9)-(gm(2,2)*(6*gm(1,3)*gm(2,2)+3*gm(1,2)*gm(2,3)))
1372      cm(1,3,9)=cm(1,3,9)-(12*gm(1,2)*gm(2,3)*gm(3,3)+gm(1,3)*(1.5d0*gm(2,3)&
1373 &     **2-4.5d0*gm(2,2)*gm(3,3)))
1374      cm(1,4,9)=cm(1,4,9)-(3*gm(1,3)*gm(2,2)*gm(2,3)+gm(1,2)*(3*gm(2,3)&
1375 &     **2+12*gm(2,2)*gm(3,3)))
1376      cm(1,5,9)=cm(1,5,9)-(-3*gm(1,3)**2*gm(2,2)+9*gm(1,2)*gm(1,3)*gm(2,3)&
1377 &     +7.5d0*gm(1,2)**2*gm(3,3)+gm(1,1)*(-3*gm(2,3)**2-1.5d0*gm(2,2)&
1378 &     *gm(3,3)))
1379      cm(1,6,9)=cm(1,6,9)-(12*gm(1,2)*gm(1,3)*gm(2,2)+1.5d0*gm(1,2)&
1380 &     **2*gm(2,3)-4.5d0*gm(1,1)*gm(2,2)*gm(2,3))
1381      cm(1,7,9)=cm(1,7,9)-(3*gm(2,2)**2*gm(2,3))
1382      cm(1,8,9)=cm(1,8,9)-(1.5d0*gm(2,3)**3+7.5d0*gm(2,2)*gm(2,3)*gm(3,3))
1383      cm(1,9,9)=cm(1,9,9)-(gm(2,2)*(3*gm(2,3)**2+6*gm(2,2)*gm(3,3)))
1384      cm(1,10,9)=cm(1,10,9)-(gm(3,3)*(4.5d0*gm(2,3)**2-1.5d0*gm(2,2)&
1385 &     *gm(3,3)))
1386      cm(1,1,10)=cm(1,1,10)-(2.5d0*gm(1,3)**3-1.5d0*gm(1,1)*gm(1,3)&
1387 &     *gm(3,3))
1388      cm(1,2,10)=cm(1,2,10)-(-3*gm(1,2)*gm(2,3)*gm(3,3)+gm(1,3)*(7.5d0*gm(2,3)&
1389 &     **2-1.5d0*gm(2,2)*gm(3,3)))
1390      cm(1,3,10)=cm(1,3,10)-(3*gm(1,3)*gm(3,3)**2)
1391      cm(1,4,10)=cm(1,4,10)-(gm(3,3)*(9*gm(1,3)*gm(2,3)-3*gm(1,2)*gm(3,3)))
1392      cm(1,5,10)=cm(1,5,10)-(gm(3,3)*(4.5d0*gm(1,3)**2-1.5d0*gm(1,1)&
1393 &     *gm(3,3)))
1394      cm(1,6,10)=cm(1,6,10)-(7.5d0*gm(1,3)**2*gm(2,3)-3*gm(1,2)*gm(1,3)&
1395 &     *gm(3,3)-1.5d0*gm(1,1)*gm(2,3)*gm(3,3))
1396      cm(1,7,10)=cm(1,7,10)-(2.5d0*gm(2,3)**3-1.5d0*gm(2,2)*gm(2,3)&
1397 &     *gm(3,3))
1398      cm(1,8,10)=cm(1,8,10)-(3*gm(2,3)*gm(3,3)**2)
1399      cm(1,9,10)=cm(1,9,10)-(gm(3,3)*(4.5d0*gm(2,3)**2-1.5d0*gm(2,2)&
1400 &     *gm(3,3)))
1401      cm(1,10,10)=cm(1,10,10)-(gm(3,3)**3)
1402      cm(2,1,2)=cm(2,1,2)-(3*gm(1,1)**2*gm(1,2))
1403      cm(2,2,2)=cm(2,2,2)-(1.5d0*gm(1,2)**3+7.5d0*gm(1,1)*gm(1,2)*gm(2,2))
1404      cm(2,3,2)=cm(2,3,2)-(12*gm(1,1)*gm(1,3)*gm(2,3)+gm(1,2)*(1.5d0*gm(1,3)&
1405 &     **2-4.5d0*gm(1,1)*gm(3,3)))
1406      cm(2,4,2)=cm(2,4,2)-(3*gm(1,2)**2*gm(1,3)+12*gm(1,1)*gm(1,3)*gm(2,2)&
1407 &     +3*gm(1,1)*gm(1,2)*gm(2,3))
1408      cm(2,5,2)=cm(2,5,2)-(gm(1,1)*(3*gm(1,2)*gm(1,3)+6*gm(1,1)*gm(2,3)))
1409      cm(2,6,2)=cm(2,6,2)-(gm(1,1)*(3*gm(1,2)**2+6*gm(1,1)*gm(2,2)))
1410      cm(2,7,2)=cm(2,7,2)-(gm(2,2)*(4.5d0*gm(1,2)**2-1.5d0*gm(1,1)*gm(2,2)))
1411      cm(2,8,2)=cm(2,8,2)-(7.5d0*gm(1,3)**2*gm(2,2)+9*gm(1,2)*gm(1,3)&
1412 &     *gm(2,3)-3*gm(1,2)**2*gm(3,3)+gm(1,1)*(-3*gm(2,3)**2-1.5d0*gm(2,2)&
1413 &     *gm(3,3)))
1414      cm(2,9,2)=cm(2,9,2)-(12*gm(1,2)*gm(1,3)*gm(2,2)+1.5d0*gm(1,2)&
1415 &     **2*gm(2,3)-4.5d0*gm(1,1)*gm(2,2)*gm(2,3))
1416      cm(2,10,2)=cm(2,10,2)-(7.5d0*gm(1,3)**2*gm(2,3)-3*gm(1,2)*gm(1,3)&
1417 &     *gm(3,3)-1.5d0*gm(1,1)*gm(2,3)*gm(3,3))
1418      cm(2,1,4)=cm(2,1,4)-(3*gm(1,1)**2*gm(1,3))
1419      cm(2,2,4)=cm(2,2,4)-(1.5d0*gm(1,2)**2*gm(1,3)-4.5d0*gm(1,1)*gm(1,3)&
1420 &     *gm(2,2)+12*gm(1,1)*gm(1,2)*gm(2,3))
1421      cm(2,3,4)=cm(2,3,4)-(1.5d0*gm(1,3)**3+7.5d0*gm(1,1)*gm(1,3)*gm(3,3))
1422      cm(2,4,4)=cm(2,4,4)-(3*gm(1,1)*gm(1,3)*gm(2,3)+gm(1,2)*(3*gm(1,3)&
1423 &     **2+12*gm(1,1)*gm(3,3)))
1424      cm(2,5,4)=cm(2,5,4)-(gm(1,1)*(3*gm(1,3)**2+6*gm(1,1)*gm(3,3)))
1425      cm(2,6,4)=cm(2,6,4)-(gm(1,1)*(3*gm(1,2)*gm(1,3)+6*gm(1,1)*gm(2,3)))
1426      cm(2,7,4)=cm(2,7,4)-(-3*gm(1,2)*gm(1,3)*gm(2,2)+7.5d0*gm(1,2)&
1427 &     **2*gm(2,3)-1.5d0*gm(1,1)*gm(2,2)*gm(2,3))
1428      cm(2,8,4)=cm(2,8,4)-(1.5d0*gm(1,3)**2*gm(2,3)+12*gm(1,2)*gm(1,3)&
1429 &     *gm(3,3)-4.5d0*gm(1,1)*gm(2,3)*gm(3,3))
1430      cm(2,9,4)=cm(2,9,4)-(-3*gm(1,3)**2*gm(2,2)+9*gm(1,2)*gm(1,3)*gm(2,3)&
1431 &     +7.5d0*gm(1,2)**2*gm(3,3)+gm(1,1)*(-3*gm(2,3)**2-1.5d0*gm(2,2)&
1432 &     *gm(3,3)))
1433      cm(2,10,4)=cm(2,10,4)-(gm(3,3)*(4.5d0*gm(1,3)**2-1.5d0*gm(1,1)&
1434 &     *gm(3,3)))
1435      cm(2,1,6)=cm(2,1,6)-(gm(1,1)**3)
1436      cm(2,2,6)=cm(2,2,6)-(gm(1,1)*(4.5d0*gm(1,2)**2-1.5d0*gm(1,1)*gm(2,2)))
1437      cm(2,3,6)=cm(2,3,6)-(gm(1,1)*(4.5d0*gm(1,3)**2-1.5d0*gm(1,1)*gm(3,3)))
1438      cm(2,4,6)=cm(2,4,6)-(gm(1,1)*(9*gm(1,2)*gm(1,3)-3*gm(1,1)*gm(2,3)))
1439      cm(2,5,6)=cm(2,5,6)-(3*gm(1,1)**2*gm(1,3))
1440      cm(2,6,6)=cm(2,6,6)-(3*gm(1,1)**2*gm(1,2))
1441      cm(2,7,6)=cm(2,7,6)-(2.5d0*gm(1,2)**3-1.5d0*gm(1,1)*gm(1,2)*gm(2,2))
1442      cm(2,8,6)=cm(2,8,6)-(-3*gm(1,1)*gm(1,3)*gm(2,3)+gm(1,2)*(7.5d0*gm(1,3)&
1443 &     **2-1.5d0*gm(1,1)*gm(3,3)))
1444      cm(2,9,6)=cm(2,9,6)-(7.5d0*gm(1,2)**2*gm(1,3)-1.5d0*gm(1,1)*gm(1,3)&
1445 &     *gm(2,2)-3*gm(1,1)*gm(1,2)*gm(2,3))
1446      cm(2,10,6)=cm(2,10,6)-(2.5d0*gm(1,3)**3-1.5d0*gm(1,1)*gm(1,3)&
1447 &     *gm(3,3))
1448      cm(2,1,7)=cm(2,1,7)-(gm(1,1)*(4.5d0*gm(1,2)**2-1.5d0*gm(1,1)*gm(2,2)))
1449      cm(2,2,7)=cm(2,2,7)-(gm(2,2)*(3*gm(1,2)**2+6*gm(1,1)*gm(2,2)))
1450      cm(2,3,7)=cm(2,3,7)-(-3*gm(1,3)**2*gm(2,2)+9*gm(1,2)*gm(1,3)*gm(2,3)&
1451 &     -3*gm(1,2)**2*gm(3,3)+gm(1,1)*(7.5d0*gm(2,3)**2-1.5d0*gm(2,2)&
1452 &     *gm(3,3)))
1453      cm(2,4,7)=cm(2,4,7)-(3*gm(1,2)*gm(1,3)*gm(2,2)+3*gm(1,2)**2*gm(2,3)&
1454 &     +12*gm(1,1)*gm(2,2)*gm(2,3))
1455      cm(2,5,7)=cm(2,5,7)-(1.5d0*gm(1,2)**2*gm(1,3)-4.5d0*gm(1,1)*gm(1,3)&
1456 &     *gm(2,2)+12*gm(1,1)*gm(1,2)*gm(2,3))
1457      cm(2,6,7)=cm(2,6,7)-(1.5d0*gm(1,2)**3+7.5d0*gm(1,1)*gm(1,2)*gm(2,2))
1458      cm(2,7,7)=cm(2,7,7)-(3*gm(1,2)*gm(2,2)**2)
1459      cm(2,8,7)=cm(2,8,7)-(12*gm(1,3)*gm(2,2)*gm(2,3)+gm(1,2)*(1.5d0*gm(2,3)&
1460 &     **2-4.5d0*gm(2,2)*gm(3,3)))
1461      cm(2,9,7)=cm(2,9,7)-(gm(2,2)*(6*gm(1,3)*gm(2,2)+3*gm(1,2)*gm(2,3)))
1462      cm(2,10,7)=cm(2,10,7)-(-3*gm(1,2)*gm(2,3)*gm(3,3)+gm(1,3)*(7.5d0*gm(2,3)&
1463 &     **2-1.5d0*gm(2,2)*gm(3,3)))
1464      cm(2,1,8)=cm(2,1,8)-(gm(1,1)*(4.5d0*gm(1,3)**2-1.5d0*gm(1,1)*gm(3,3)))
1465      cm(2,2,8)=cm(2,2,8)-(-3*gm(1,3)**2*gm(2,2)+9*gm(1,2)*gm(1,3)*gm(2,3)&
1466 &     -3*gm(1,2)**2*gm(3,3)+gm(1,1)*(7.5d0*gm(2,3)**2-1.5d0*gm(2,2)&
1467 &     *gm(3,3)))
1468      cm(2,3,8)=cm(2,3,8)-(gm(3,3)*(3*gm(1,3)**2+6*gm(1,1)*gm(3,3)))
1469      cm(2,4,8)=cm(2,4,8)-(3*gm(1,3)**2*gm(2,3)+3*gm(1,2)*gm(1,3)*gm(3,3)&
1470 &     +12*gm(1,1)*gm(2,3)*gm(3,3))
1471      cm(2,5,8)=cm(2,5,8)-(1.5d0*gm(1,3)**3+7.5d0*gm(1,1)*gm(1,3)*gm(3,3))
1472      cm(2,6,8)=cm(2,6,8)-(12*gm(1,1)*gm(1,3)*gm(2,3)+gm(1,2)*(1.5d0*gm(1,3)&
1473 &     **2-4.5d0*gm(1,1)*gm(3,3)))
1474      cm(2,7,8)=cm(2,7,8)-(-3*gm(1,3)*gm(2,2)*gm(2,3)+gm(1,2)*(7.5d0*gm(2,3)&
1475 &     **2-1.5d0*gm(2,2)*gm(3,3)))
1476      cm(2,8,8)=cm(2,8,8)-(gm(3,3)*(3*gm(1,3)*gm(2,3)+6*gm(1,2)*gm(3,3)))
1477      cm(2,9,8)=cm(2,9,8)-(12*gm(1,2)*gm(2,3)*gm(3,3)+gm(1,3)*(1.5d0*gm(2,3)&
1478 &     **2-4.5d0*gm(2,2)*gm(3,3)))
1479      cm(2,10,8)=cm(2,10,8)-(3*gm(1,3)*gm(3,3)**2)
1480      cm(2,1,9)=cm(2,1,9)-(gm(1,1)*(9*gm(1,2)*gm(1,3)-3*gm(1,1)*gm(2,3)))
1481      cm(2,2,9)=cm(2,2,9)-(3*gm(1,2)*gm(1,3)*gm(2,2)+3*gm(1,2)**2*gm(2,3)&
1482 &     +12*gm(1,1)*gm(2,2)*gm(2,3))
1483      cm(2,3,9)=cm(2,3,9)-(3*gm(1,3)**2*gm(2,3)+3*gm(1,2)*gm(1,3)*gm(3,3)&
1484 &     +12*gm(1,1)*gm(2,3)*gm(3,3))
1485      cm(2,4,9)=cm(2,4,9)-(9*gm(1,3)**2*gm(2,2)-6*gm(1,2)*gm(1,3)*gm(2,3)&
1486 &     +9*gm(1,2)**2*gm(3,3)+gm(1,1)*(9*gm(2,3)**2+15*gm(2,2)*gm(3,3)))
1487      cm(2,5,9)=cm(2,5,9)-(3*gm(1,1)*gm(1,3)*gm(2,3)+gm(1,2)*(3*gm(1,3)&
1488 &     **2+12*gm(1,1)*gm(3,3)))
1489      cm(2,6,9)=cm(2,6,9)-(3*gm(1,2)**2*gm(1,3)+12*gm(1,1)*gm(1,3)*gm(2,2)&
1490 &     +3*gm(1,1)*gm(1,2)*gm(2,3))
1491      cm(2,7,9)=cm(2,7,9)-(gm(2,2)*(-3*gm(1,3)*gm(2,2)+9*gm(1,2)*gm(2,3)))
1492      cm(2,8,9)=cm(2,8,9)-(3*gm(1,2)*gm(2,3)*gm(3,3)+gm(1,3)*(3*gm(2,3)&
1493 &     **2+12*gm(2,2)*gm(3,3)))
1494      cm(2,9,9)=cm(2,9,9)-(3*gm(1,3)*gm(2,2)*gm(2,3)+gm(1,2)*(3*gm(2,3)&
1495 &     **2+12*gm(2,2)*gm(3,3)))
1496      cm(2,10,9)=cm(2,10,9)-(gm(3,3)*(9*gm(1,3)*gm(2,3)-3*gm(1,2)*gm(3,3)))
1497      cm(2,1,11)=cm(2,1,11)-(2.5d0*gm(1,2)**3-1.5d0*gm(1,1)*gm(1,2)&
1498 &     *gm(2,2))
1499      cm(2,2,11)=cm(2,2,11)-(3*gm(1,2)*gm(2,2)**2)
1500      cm(2,3,11)=cm(2,3,11)-(-3*gm(1,3)*gm(2,2)*gm(2,3)+gm(1,2)*(7.5d0*gm(2,3)&
1501 &     **2-1.5d0*gm(2,2)*gm(3,3)))
1502      cm(2,4,11)=cm(2,4,11)-(gm(2,2)*(-3*gm(1,3)*gm(2,2)+9*gm(1,2)*gm(2,3)))
1503      cm(2,5,11)=cm(2,5,11)-(-3*gm(1,2)*gm(1,3)*gm(2,2)+7.5d0*gm(1,2)&
1504 &     **2*gm(2,3)-1.5d0*gm(1,1)*gm(2,2)*gm(2,3))
1505      cm(2,6,11)=cm(2,6,11)-(gm(2,2)*(4.5d0*gm(1,2)**2-1.5d0*gm(1,1)&
1506 &     *gm(2,2)))
1507      cm(2,7,11)=cm(2,7,11)-(gm(2,2)**3)
1508      cm(2,8,11)=cm(2,8,11)-(gm(2,2)*(4.5d0*gm(2,3)**2-1.5d0*gm(2,2)&
1509 &     *gm(3,3)))
1510      cm(2,9,11)=cm(2,9,11)-(3*gm(2,2)**2*gm(2,3))
1511      cm(2,10,11)=cm(2,10,11)-(2.5d0*gm(2,3)**3-1.5d0*gm(2,2)*gm(2,3)&
1512 &     *gm(3,3))
1513      cm(2,1,12)=cm(2,1,12)-(-3*gm(1,1)*gm(1,3)*gm(2,3)+gm(1,2)*(7.5d0*gm(1,3)&
1514 &     **2-1.5d0*gm(1,1)*gm(3,3)))
1515      cm(2,2,12)=cm(2,2,12)-(12*gm(1,3)*gm(2,2)*gm(2,3)+gm(1,2)*(1.5d0*gm(2,3)&
1516 &     **2-4.5d0*gm(2,2)*gm(3,3)))
1517      cm(2,3,12)=cm(2,3,12)-(gm(3,3)*(3*gm(1,3)*gm(2,3)+6*gm(1,2)*gm(3,3)))
1518      cm(2,4,12)=cm(2,4,12)-(3*gm(1,2)*gm(2,3)*gm(3,3)+gm(1,3)*(3*gm(2,3)&
1519 &     **2+12*gm(2,2)*gm(3,3)))
1520      cm(2,5,12)=cm(2,5,12)-(1.5d0*gm(1,3)**2*gm(2,3)+12*gm(1,2)*gm(1,3)&
1521 &     *gm(3,3)-4.5d0*gm(1,1)*gm(2,3)*gm(3,3))
1522      cm(2,6,12)=cm(2,6,12)-(7.5d0*gm(1,3)**2*gm(2,2)+9*gm(1,2)*gm(1,3)&
1523 &     *gm(2,3)-3*gm(1,2)**2*gm(3,3)+gm(1,1)*(-3*gm(2,3)**2-1.5d0*gm(2,2)&
1524 &     *gm(3,3)))
1525      cm(2,7,12)=cm(2,7,12)-(gm(2,2)*(4.5d0*gm(2,3)**2-1.5d0*gm(2,2)&
1526 &     *gm(3,3)))
1527      cm(2,8,12)=cm(2,8,12)-(gm(3,3)*(3*gm(2,3)**2+6*gm(2,2)*gm(3,3)))
1528      cm(2,9,12)=cm(2,9,12)-(1.5d0*gm(2,3)**3+7.5d0*gm(2,2)*gm(2,3)&
1529 &     *gm(3,3))
1530      cm(2,10,12)=cm(2,10,12)-(3*gm(2,3)*gm(3,3)**2)
1531      cm(2,1,13)=cm(2,1,13)-(7.5d0*gm(1,2)**2*gm(1,3)-1.5d0*gm(1,1)&
1532 &     *gm(1,3)*gm(2,2)-3*gm(1,1)*gm(1,2)*gm(2,3))
1533      cm(2,2,13)=cm(2,2,13)-(gm(2,2)*(6*gm(1,3)*gm(2,2)+3*gm(1,2)*gm(2,3)))
1534      cm(2,3,13)=cm(2,3,13)-(12*gm(1,2)*gm(2,3)*gm(3,3)+gm(1,3)*(1.5d0*gm(2,3)&
1535 &     **2-4.5d0*gm(2,2)*gm(3,3)))
1536      cm(2,4,13)=cm(2,4,13)-(3*gm(1,3)*gm(2,2)*gm(2,3)+gm(1,2)*(3*gm(2,3)&
1537 &     **2+12*gm(2,2)*gm(3,3)))
1538      cm(2,5,13)=cm(2,5,13)-(-3*gm(1,3)**2*gm(2,2)+9*gm(1,2)*gm(1,3)&
1539 &     *gm(2,3)+7.5d0*gm(1,2)**2*gm(3,3)+gm(1,1)*(-3*gm(2,3)**2-1.5d0*gm(2,2)&
1540 &     *gm(3,3)))
1541      cm(2,6,13)=cm(2,6,13)-(12*gm(1,2)*gm(1,3)*gm(2,2)+1.5d0*gm(1,2)&
1542 &     **2*gm(2,3)-4.5d0*gm(1,1)*gm(2,2)*gm(2,3))
1543      cm(2,7,13)=cm(2,7,13)-(3*gm(2,2)**2*gm(2,3))
1544      cm(2,8,13)=cm(2,8,13)-(1.5d0*gm(2,3)**3+7.5d0*gm(2,2)*gm(2,3)&
1545 &     *gm(3,3))
1546      cm(2,9,13)=cm(2,9,13)-(gm(2,2)*(3*gm(2,3)**2+6*gm(2,2)*gm(3,3)))
1547      cm(2,10,13)=cm(2,10,13)-(gm(3,3)*(4.5d0*gm(2,3)**2-1.5d0*gm(2,2)&
1548 &     *gm(3,3)))
1549      cm(2,1,14)=cm(2,1,14)-(2.5d0*gm(1,3)**3-1.5d0*gm(1,1)*gm(1,3)&
1550 &     *gm(3,3))
1551      cm(2,2,14)=cm(2,2,14)-(-3*gm(1,2)*gm(2,3)*gm(3,3)+gm(1,3)*(7.5d0*gm(2,3)&
1552 &     **2-1.5d0*gm(2,2)*gm(3,3)))
1553      cm(2,3,14)=cm(2,3,14)-(3*gm(1,3)*gm(3,3)**2)
1554      cm(2,4,14)=cm(2,4,14)-(gm(3,3)*(9*gm(1,3)*gm(2,3)-3*gm(1,2)*gm(3,3)))
1555      cm(2,5,14)=cm(2,5,14)-(gm(3,3)*(4.5d0*gm(1,3)**2-1.5d0*gm(1,1)&
1556 &     *gm(3,3)))
1557      cm(2,6,14)=cm(2,6,14)-(7.5d0*gm(1,3)**2*gm(2,3)-3*gm(1,2)*gm(1,3)&
1558 &     *gm(3,3)-1.5d0*gm(1,1)*gm(2,3)*gm(3,3))
1559      cm(2,7,14)=cm(2,7,14)-(2.5d0*gm(2,3)**3-1.5d0*gm(2,2)*gm(2,3)&
1560 &     *gm(3,3))
1561      cm(2,8,14)=cm(2,8,14)-(3*gm(2,3)*gm(3,3)**2)
1562      cm(2,9,14)=cm(2,9,14)-(gm(3,3)*(4.5d0*gm(2,3)**2-1.5d0*gm(2,2)&
1563 &     *gm(3,3)))
1564      cm(2,10,14)=cm(2,10,14)-(gm(3,3)**3)
1565      cm(3,1,3)=cm(3,1,3)-(3*gm(1,1)**2*gm(1,3))
1566      cm(3,2,3)=cm(3,2,3)-(1.5d0*gm(1,2)**2*gm(1,3)-4.5d0*gm(1,1)*gm(1,3)&
1567 &     *gm(2,2)+12*gm(1,1)*gm(1,2)*gm(2,3))
1568      cm(3,3,3)=cm(3,3,3)-(1.5d0*gm(1,3)**3+7.5d0*gm(1,1)*gm(1,3)*gm(3,3))
1569      cm(3,4,3)=cm(3,4,3)-(3*gm(1,1)*gm(1,3)*gm(2,3)+gm(1,2)*(3*gm(1,3)&
1570 &     **2+12*gm(1,1)*gm(3,3)))
1571      cm(3,5,3)=cm(3,5,3)-(gm(1,1)*(3*gm(1,3)**2+6*gm(1,1)*gm(3,3)))
1572      cm(3,6,3)=cm(3,6,3)-(gm(1,1)*(3*gm(1,2)*gm(1,3)+6*gm(1,1)*gm(2,3)))
1573      cm(3,7,3)=cm(3,7,3)-(-3*gm(1,2)*gm(1,3)*gm(2,2)+7.5d0*gm(1,2)&
1574 &     **2*gm(2,3)-1.5d0*gm(1,1)*gm(2,2)*gm(2,3))
1575      cm(3,8,3)=cm(3,8,3)-(1.5d0*gm(1,3)**2*gm(2,3)+12*gm(1,2)*gm(1,3)&
1576 &     *gm(3,3)-4.5d0*gm(1,1)*gm(2,3)*gm(3,3))
1577      cm(3,9,3)=cm(3,9,3)-(-3*gm(1,3)**2*gm(2,2)+9*gm(1,2)*gm(1,3)*gm(2,3)&
1578 &     +7.5d0*gm(1,2)**2*gm(3,3)+gm(1,1)*(-3*gm(2,3)**2-1.5d0*gm(2,2)&
1579 &     *gm(3,3)))
1580      cm(3,10,3)=cm(3,10,3)-(gm(3,3)*(4.5d0*gm(1,3)**2-1.5d0*gm(1,1)&
1581 &     *gm(3,3)))
1582      cm(3,1,4)=cm(3,1,4)-(3*gm(1,1)**2*gm(1,2))
1583      cm(3,2,4)=cm(3,2,4)-(1.5d0*gm(1,2)**3+7.5d0*gm(1,1)*gm(1,2)*gm(2,2))
1584      cm(3,3,4)=cm(3,3,4)-(12*gm(1,1)*gm(1,3)*gm(2,3)+gm(1,2)*(1.5d0*gm(1,3)&
1585 &     **2-4.5d0*gm(1,1)*gm(3,3)))
1586      cm(3,4,4)=cm(3,4,4)-(3*gm(1,2)**2*gm(1,3)+12*gm(1,1)*gm(1,3)*gm(2,2)&
1587 &     +3*gm(1,1)*gm(1,2)*gm(2,3))
1588      cm(3,5,4)=cm(3,5,4)-(gm(1,1)*(3*gm(1,2)*gm(1,3)+6*gm(1,1)*gm(2,3)))
1589      cm(3,6,4)=cm(3,6,4)-(gm(1,1)*(3*gm(1,2)**2+6*gm(1,1)*gm(2,2)))
1590      cm(3,7,4)=cm(3,7,4)-(gm(2,2)*(4.5d0*gm(1,2)**2-1.5d0*gm(1,1)*gm(2,2)))
1591      cm(3,8,4)=cm(3,8,4)-(7.5d0*gm(1,3)**2*gm(2,2)+9*gm(1,2)*gm(1,3)&
1592 &     *gm(2,3)-3*gm(1,2)**2*gm(3,3)+gm(1,1)*(-3*gm(2,3)**2-1.5d0*gm(2,2)&
1593 &     *gm(3,3)))
1594      cm(3,9,4)=cm(3,9,4)-(12*gm(1,2)*gm(1,3)*gm(2,2)+1.5d0*gm(1,2)&
1595 &     **2*gm(2,3)-4.5d0*gm(1,1)*gm(2,2)*gm(2,3))
1596      cm(3,10,4)=cm(3,10,4)-(7.5d0*gm(1,3)**2*gm(2,3)-3*gm(1,2)*gm(1,3)&
1597 &     *gm(3,3)-1.5d0*gm(1,1)*gm(2,3)*gm(3,3))
1598      cm(3,1,5)=cm(3,1,5)-(gm(1,1)**3)
1599      cm(3,2,5)=cm(3,2,5)-(gm(1,1)*(4.5d0*gm(1,2)**2-1.5d0*gm(1,1)*gm(2,2)))
1600      cm(3,3,5)=cm(3,3,5)-(gm(1,1)*(4.5d0*gm(1,3)**2-1.5d0*gm(1,1)*gm(3,3)))
1601      cm(3,4,5)=cm(3,4,5)-(gm(1,1)*(9*gm(1,2)*gm(1,3)-3*gm(1,1)*gm(2,3)))
1602      cm(3,5,5)=cm(3,5,5)-(3*gm(1,1)**2*gm(1,3))
1603      cm(3,6,5)=cm(3,6,5)-(3*gm(1,1)**2*gm(1,2))
1604      cm(3,7,5)=cm(3,7,5)-(2.5d0*gm(1,2)**3-1.5d0*gm(1,1)*gm(1,2)*gm(2,2))
1605      cm(3,8,5)=cm(3,8,5)-(-3*gm(1,1)*gm(1,3)*gm(2,3)+gm(1,2)*(7.5d0*gm(1,3)&
1606 &     **2-1.5d0*gm(1,1)*gm(3,3)))
1607      cm(3,9,5)=cm(3,9,5)-(7.5d0*gm(1,2)**2*gm(1,3)-1.5d0*gm(1,1)*gm(1,3)&
1608 &     *gm(2,2)-3*gm(1,1)*gm(1,2)*gm(2,3))
1609      cm(3,10,5)=cm(3,10,5)-(2.5d0*gm(1,3)**3-1.5d0*gm(1,1)*gm(1,3)&
1610 &     *gm(3,3))
1611      cm(3,1,8)=cm(3,1,8)-(gm(1,1)*(9*gm(1,2)*gm(1,3)-3*gm(1,1)*gm(2,3)))
1612      cm(3,2,8)=cm(3,2,8)-(3*gm(1,2)*gm(1,3)*gm(2,2)+3*gm(1,2)**2*gm(2,3)&
1613 &     +12*gm(1,1)*gm(2,2)*gm(2,3))
1614      cm(3,3,8)=cm(3,3,8)-(3*gm(1,3)**2*gm(2,3)+3*gm(1,2)*gm(1,3)*gm(3,3)&
1615 &     +12*gm(1,1)*gm(2,3)*gm(3,3))
1616      cm(3,4,8)=cm(3,4,8)-(9*gm(1,3)**2*gm(2,2)-6*gm(1,2)*gm(1,3)*gm(2,3)&
1617 &     +9*gm(1,2)**2*gm(3,3)+gm(1,1)*(9*gm(2,3)**2+15*gm(2,2)*gm(3,3)))
1618      cm(3,5,8)=cm(3,5,8)-(3*gm(1,1)*gm(1,3)*gm(2,3)+gm(1,2)*(3*gm(1,3)&
1619 &     **2+12*gm(1,1)*gm(3,3)))
1620      cm(3,6,8)=cm(3,6,8)-(3*gm(1,2)**2*gm(1,3)+12*gm(1,1)*gm(1,3)*gm(2,2)&
1621 &     +3*gm(1,1)*gm(1,2)*gm(2,3))
1622      cm(3,7,8)=cm(3,7,8)-(gm(2,2)*(-3*gm(1,3)*gm(2,2)+9*gm(1,2)*gm(2,3)))
1623      cm(3,8,8)=cm(3,8,8)-(3*gm(1,2)*gm(2,3)*gm(3,3)+gm(1,3)*(3*gm(2,3)&
1624 &     **2+12*gm(2,2)*gm(3,3)))
1625      cm(3,9,8)=cm(3,9,8)-(3*gm(1,3)*gm(2,2)*gm(2,3)+gm(1,2)*(3*gm(2,3)&
1626 &     **2+12*gm(2,2)*gm(3,3)))
1627      cm(3,10,8)=cm(3,10,8)-(gm(3,3)*(9*gm(1,3)*gm(2,3)-3*gm(1,2)*gm(3,3)))
1628      cm(3,1,9)=cm(3,1,9)-(gm(1,1)*(4.5d0*gm(1,2)**2-1.5d0*gm(1,1)*gm(2,2)))
1629      cm(3,2,9)=cm(3,2,9)-(gm(2,2)*(3*gm(1,2)**2+6*gm(1,1)*gm(2,2)))
1630      cm(3,3,9)=cm(3,3,9)-(-3*gm(1,3)**2*gm(2,2)+9*gm(1,2)*gm(1,3)*gm(2,3)&
1631 &     -3*gm(1,2)**2*gm(3,3)+gm(1,1)*(7.5d0*gm(2,3)**2-1.5d0*gm(2,2)&
1632 &     *gm(3,3)))
1633      cm(3,4,9)=cm(3,4,9)-(3*gm(1,2)*gm(1,3)*gm(2,2)+3*gm(1,2)**2*gm(2,3)&
1634 &     +12*gm(1,1)*gm(2,2)*gm(2,3))
1635      cm(3,5,9)=cm(3,5,9)-(1.5d0*gm(1,2)**2*gm(1,3)-4.5d0*gm(1,1)*gm(1,3)&
1636 &     *gm(2,2)+12*gm(1,1)*gm(1,2)*gm(2,3))
1637      cm(3,6,9)=cm(3,6,9)-(1.5d0*gm(1,2)**3+7.5d0*gm(1,1)*gm(1,2)*gm(2,2))
1638      cm(3,7,9)=cm(3,7,9)-(3*gm(1,2)*gm(2,2)**2)
1639      cm(3,8,9)=cm(3,8,9)-(12*gm(1,3)*gm(2,2)*gm(2,3)+gm(1,2)*(1.5d0*gm(2,3)&
1640 &     **2-4.5d0*gm(2,2)*gm(3,3)))
1641      cm(3,9,9)=cm(3,9,9)-(gm(2,2)*(6*gm(1,3)*gm(2,2)+3*gm(1,2)*gm(2,3)))
1642      cm(3,10,9)=cm(3,10,9)-(-3*gm(1,2)*gm(2,3)*gm(3,3)+gm(1,3)*(7.5d0*gm(2,3)&
1643 &     **2-1.5d0*gm(2,2)*gm(3,3)))
1644      cm(3,1,10)=cm(3,1,10)-(gm(1,1)*(4.5d0*gm(1,3)**2-1.5d0*gm(1,1)&
1645 &     *gm(3,3)))
1646      cm(3,2,10)=cm(3,2,10)-(-3*gm(1,3)**2*gm(2,2)+9*gm(1,2)*gm(1,3)&
1647 &     *gm(2,3)-3*gm(1,2)**2*gm(3,3)+gm(1,1)*(7.5d0*gm(2,3)**2-1.5d0*gm(2,2)&
1648 &     *gm(3,3)))
1649      cm(3,3,10)=cm(3,3,10)-(gm(3,3)*(3*gm(1,3)**2+6*gm(1,1)*gm(3,3)))
1650      cm(3,4,10)=cm(3,4,10)-(3*gm(1,3)**2*gm(2,3)+3*gm(1,2)*gm(1,3)&
1651 &     *gm(3,3)+12*gm(1,1)*gm(2,3)*gm(3,3))
1652      cm(3,5,10)=cm(3,5,10)-(1.5d0*gm(1,3)**3+7.5d0*gm(1,1)*gm(1,3)&
1653 &     *gm(3,3))
1654      cm(3,6,10)=cm(3,6,10)-(12*gm(1,1)*gm(1,3)*gm(2,3)+gm(1,2)*(1.5d0*gm(1,3)&
1655 &     **2-4.5d0*gm(1,1)*gm(3,3)))
1656      cm(3,7,10)=cm(3,7,10)-(-3*gm(1,3)*gm(2,2)*gm(2,3)+gm(1,2)*(7.5d0*gm(2,3)&
1657 &     **2-1.5d0*gm(2,2)*gm(3,3)))
1658      cm(3,8,10)=cm(3,8,10)-(gm(3,3)*(3*gm(1,3)*gm(2,3)+6*gm(1,2)*gm(3,3)))
1659      cm(3,9,10)=cm(3,9,10)-(12*gm(1,2)*gm(2,3)*gm(3,3)+gm(1,3)*(1.5d0*gm(2,3)&
1660 &     **2-4.5d0*gm(2,2)*gm(3,3)))
1661      cm(3,10,10)=cm(3,10,10)-(3*gm(1,3)*gm(3,3)**2)
1662      cm(3,1,12)=cm(3,1,12)-(7.5d0*gm(1,2)**2*gm(1,3)-1.5d0*gm(1,1)&
1663 &     *gm(1,3)*gm(2,2)-3*gm(1,1)*gm(1,2)*gm(2,3))
1664      cm(3,2,12)=cm(3,2,12)-(gm(2,2)*(6*gm(1,3)*gm(2,2)+3*gm(1,2)*gm(2,3)))
1665      cm(3,3,12)=cm(3,3,12)-(12*gm(1,2)*gm(2,3)*gm(3,3)+gm(1,3)*(1.5d0*gm(2,3)&
1666 &     **2-4.5d0*gm(2,2)*gm(3,3)))
1667      cm(3,4,12)=cm(3,4,12)-(3*gm(1,3)*gm(2,2)*gm(2,3)+gm(1,2)*(3*gm(2,3)&
1668 &     **2+12*gm(2,2)*gm(3,3)))
1669      cm(3,5,12)=cm(3,5,12)-(-3*gm(1,3)**2*gm(2,2)+9*gm(1,2)*gm(1,3)&
1670 &     *gm(2,3)+7.5d0*gm(1,2)**2*gm(3,3)+gm(1,1)*(-3*gm(2,3)**2-1.5d0*gm(2,2)&
1671 &     *gm(3,3)))
1672      cm(3,6,12)=cm(3,6,12)-(12*gm(1,2)*gm(1,3)*gm(2,2)+1.5d0*gm(1,2)&
1673 &     **2*gm(2,3)-4.5d0*gm(1,1)*gm(2,2)*gm(2,3))
1674      cm(3,7,12)=cm(3,7,12)-(3*gm(2,2)**2*gm(2,3))
1675      cm(3,8,12)=cm(3,8,12)-(1.5d0*gm(2,3)**3+7.5d0*gm(2,2)*gm(2,3)&
1676 &     *gm(3,3))
1677      cm(3,9,12)=cm(3,9,12)-(gm(2,2)*(3*gm(2,3)**2+6*gm(2,2)*gm(3,3)))
1678      cm(3,10,12)=cm(3,10,12)-(gm(3,3)*(4.5d0*gm(2,3)**2-1.5d0*gm(2,2)&
1679 &     *gm(3,3)))
1680      cm(3,1,13)=cm(3,1,13)-(2.5d0*gm(1,2)**3-1.5d0*gm(1,1)*gm(1,2)&
1681 &     *gm(2,2))
1682      cm(3,2,13)=cm(3,2,13)-(3*gm(1,2)*gm(2,2)**2)
1683      cm(3,3,13)=cm(3,3,13)-(-3*gm(1,3)*gm(2,2)*gm(2,3)+gm(1,2)*(7.5d0*gm(2,3)&
1684 &     **2-1.5d0*gm(2,2)*gm(3,3)))
1685      cm(3,4,13)=cm(3,4,13)-(gm(2,2)*(-3*gm(1,3)*gm(2,2)+9*gm(1,2)*gm(2,3)))
1686      cm(3,5,13)=cm(3,5,13)-(-3*gm(1,2)*gm(1,3)*gm(2,2)+7.5d0*gm(1,2)&
1687 &     **2*gm(2,3)-1.5d0*gm(1,1)*gm(2,2)*gm(2,3))
1688      cm(3,6,13)=cm(3,6,13)-(gm(2,2)*(4.5d0*gm(1,2)**2-1.5d0*gm(1,1)&
1689 &     *gm(2,2)))
1690      cm(3,7,13)=cm(3,7,13)-(gm(2,2)**3)
1691      cm(3,8,13)=cm(3,8,13)-(gm(2,2)*(4.5d0*gm(2,3)**2-1.5d0*gm(2,2)&
1692 &     *gm(3,3)))
1693      cm(3,9,13)=cm(3,9,13)-(3*gm(2,2)**2*gm(2,3))
1694      cm(3,10,13)=cm(3,10,13)-(2.5d0*gm(2,3)**3-1.5d0*gm(2,2)*gm(2,3)&
1695 &     *gm(3,3))
1696      cm(3,1,14)=cm(3,1,14)-(-3*gm(1,1)*gm(1,3)*gm(2,3)+gm(1,2)*(7.5d0*gm(1,3)&
1697 &     **2-1.5d0*gm(1,1)*gm(3,3)))
1698      cm(3,2,14)=cm(3,2,14)-(12*gm(1,3)*gm(2,2)*gm(2,3)+gm(1,2)*(1.5d0*gm(2,3)&
1699 &     **2-4.5d0*gm(2,2)*gm(3,3)))
1700      cm(3,3,14)=cm(3,3,14)-(gm(3,3)*(3*gm(1,3)*gm(2,3)+6*gm(1,2)*gm(3,3)))
1701      cm(3,4,14)=cm(3,4,14)-(3*gm(1,2)*gm(2,3)*gm(3,3)+gm(1,3)*(3*gm(2,3)&
1702 &     **2+12*gm(2,2)*gm(3,3)))
1703      cm(3,5,14)=cm(3,5,14)-(1.5d0*gm(1,3)**2*gm(2,3)+12*gm(1,2)*gm(1,3)&
1704 &     *gm(3,3)-4.5d0*gm(1,1)*gm(2,3)*gm(3,3))
1705      cm(3,6,14)=cm(3,6,14)-(7.5d0*gm(1,3)**2*gm(2,2)+9*gm(1,2)*gm(1,3)&
1706 &     *gm(2,3)-3*gm(1,2)**2*gm(3,3)+gm(1,1)*(-3*gm(2,3)**2-1.5d0*gm(2,2)&
1707 &     *gm(3,3)))
1708      cm(3,7,14)=cm(3,7,14)-(gm(2,2)*(4.5d0*gm(2,3)**2-1.5d0*gm(2,2)&
1709 &     *gm(3,3)))
1710      cm(3,8,14)=cm(3,8,14)-(gm(3,3)*(3*gm(2,3)**2+6*gm(2,2)*gm(3,3)))
1711      cm(3,9,14)=cm(3,9,14)-(1.5d0*gm(2,3)**3+7.5d0*gm(2,2)*gm(2,3)&
1712 &     *gm(3,3))
1713      cm(3,10,14)=cm(3,10,14)-(3*gm(2,3)*gm(3,3)**2)
1714      cm(3,1,15)=cm(3,1,15)-(2.5d0*gm(1,3)**3-1.5d0*gm(1,1)*gm(1,3)&
1715 &     *gm(3,3))
1716      cm(3,2,15)=cm(3,2,15)-(-3*gm(1,2)*gm(2,3)*gm(3,3)+gm(1,3)*(7.5d0*gm(2,3)&
1717 &     **2-1.5d0*gm(2,2)*gm(3,3)))
1718      cm(3,3,15)=cm(3,3,15)-(3*gm(1,3)*gm(3,3)**2)
1719      cm(3,4,15)=cm(3,4,15)-(gm(3,3)*(9*gm(1,3)*gm(2,3)-3*gm(1,2)*gm(3,3)))
1720      cm(3,5,15)=cm(3,5,15)-(gm(3,3)*(4.5d0*gm(1,3)**2-1.5d0*gm(1,1)&
1721 &     *gm(3,3)))
1722      cm(3,6,15)=cm(3,6,15)-(7.5d0*gm(1,3)**2*gm(2,3)-3*gm(1,2)*gm(1,3)&
1723 &     *gm(3,3)-1.5d0*gm(1,1)*gm(2,3)*gm(3,3))
1724      cm(3,7,15)=cm(3,7,15)-(2.5d0*gm(2,3)**3-1.5d0*gm(2,2)*gm(2,3)&
1725 &     *gm(3,3))
1726      cm(3,8,15)=cm(3,8,15)-(3*gm(2,3)*gm(3,3)**2)
1727      cm(3,9,15)=cm(3,9,15)-(gm(3,3)*(4.5d0*gm(2,3)**2-1.5d0*gm(2,2)&
1728 &     *gm(3,3)))
1729      cm(3,10,15)=cm(3,10,15)-(gm(3,3)**3)
1730    end if
1731  end if
1732 !
1733 !contraction to output 3-vector
1734 !
1735  eisnl(:)=0.d0
1736  do jj=1,((rank+2)*(rank+3))/2
1737    tmp(:,:)=0.d0
1738    do ii=1,((rank+1)*(rank+2))/2
1739      tmp(:,1)=tmp(:,1)+aa(:,ii)*cm(1,ii,jj)
1740      tmp(:,2)=tmp(:,2)+aa(:,ii)*cm(2,ii,jj)
1741      tmp(:,3)=tmp(:,3)+aa(:,ii)*cm(3,ii,jj)
1742    end do
1743    eisnl(:)=eisnl(:)+tmp(1,:)*bb(1,jj)+tmp(2,:)*bb(2,jj)
1744  end do
1745  eisnl(:)=2.0d0*eisnl(:)
1746 !factor of 2 multiplied in to drop call to conjugate contraction
1747 !eisnl(:)=0.5d0*eisnl(:)
1748 
1749  ABI_DEALLOCATE(cm)
1750 
1751 end subroutine contistr01