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