TABLE OF CONTENTS


ABINIT/contstr26 [ Functions ]

[ Top ] [ Functions ]

NAME

 contstr26

FUNCTION

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

COPYRIGHT

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

INPUTS

  istr1=1,...6 specifies cartesian strain component 11,22,33,32,31,21
  istr2=seccond strain component
  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

  e2nl=contraction for nonlocal 2nd-order strain derivative energy

NOTES

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

PARENTS

      nonlop_pl

CHILDREN

SOURCE

 48 #if defined HAVE_CONFIG_H
 49 #include "config.h"
 50 #endif
 51 
 52 #include "abi_common.h"
 53 
 54 
 55 module m_contstr26
 56 contains
 57 !!**
 58 
 59 
 60 
 61 subroutine contstr26(istr1,istr2,rank,gm,gprimd,e2nl,aa,bb)
 62 
 63  use defs_basis
 64  use m_abicore
 65 
 66 !This section has been created automatically by the script Abilint (TD).
 67 !Do not modify the following lines by hand.
 68 #undef ABI_FUNC
 69 #define ABI_FUNC 'contstr26'
 70 !End of the abilint section
 71 
 72  implicit none
 73 
 74 !Arguments ------------------------------------
 75 !scalars
 76  integer,intent(in) :: istr1,istr2,rank
 77  real(dp),intent(out) :: e2nl
 78 !arrays
 79  real(dp),intent(in) :: aa(2,((rank+1)*(rank+2))/2),bb(2,((rank+1)*(rank+2))/2)
 80  real(dp),intent(in) :: gm(3,3),gprimd(3,3)
 81 
 82 !Local variables-------------------------------
 83 !scalars
 84  integer,parameter :: mrank=3
 85  integer :: ii,jj,ka,kb,kd,kg
 86 !arrays
 87  integer,save :: idx(12)=(/1,1,2,2,3,3,3,2,3,1,2,1/)
 88  real(dp) :: d2gm(3,3),dgm01(3,3),dgm10(3,3),tmp(2)
 89  real(dp),allocatable :: cm(:,:)
 90 
 91 ! *************************************************************************
 92  ABI_ALLOCATE(cm,(((mrank+1)*(mrank+2))/2,((mrank+1)*(mrank+2))/2))
 93 
 94  ka=idx(2*istr1-1);kb=idx(2*istr1);kg=idx(2*istr2-1);kd=idx(2*istr2)
 95 
 96  do ii = 1,3
 97    dgm01(:,ii)=-(gprimd(ka,:)*gprimd(kb,ii)+gprimd(kb,:)*gprimd(ka,ii))
 98    dgm10(:,ii)=-(gprimd(kg,:)*gprimd(kd,ii)+gprimd(kd,:)*gprimd(kg,ii))
 99  end do
100 
101  d2gm(:,:)=0.d0
102  do ii = 1,3
103    if(ka==kg) d2gm(:,ii)=d2gm(:,ii)&
104 &   +gprimd(kb,:)*gprimd(kd,ii)+gprimd(kd,:)*gprimd(kb,ii)
105    if(ka==kd) d2gm(:,ii)=d2gm(:,ii)&
106 &   +gprimd(kb,:)*gprimd(kg,ii)+gprimd(kg,:)*gprimd(kb,ii)
107    if(kb==kg) d2gm(:,ii)=d2gm(:,ii)&
108 &   +gprimd(ka,:)*gprimd(kd,ii)+gprimd(kd,:)*gprimd(ka,ii)
109    if(kb==kd) d2gm(:,ii)=d2gm(:,ii)&
110 &   +gprimd(ka,:)*gprimd(kg,ii)+gprimd(kg,:)*gprimd(ka,ii)
111  end do
112  d2gm(:,:)=0.5d0*d2gm(:,:)
113 
114 !
115 !The code below was written by a Mathematica program and formatted by
116 !a combination of editing scripts.  It is not intended to be read
117 !by human beings, and certainly not to be modified by one.  Conceivably
118 !it could be shortened somewhat by identifying common subexpressions.
119 !
120  if(rank==0)then
121    cm(1,1)=0.0d0
122  elseif(rank==1)then
123    cm(1,1)=d2gm(1,1)
124    cm(2,1)=d2gm(1,2)
125    cm(3,1)=d2gm(1,3)
126    cm(1,2)=d2gm(1,2)
127    cm(2,2)=d2gm(2,2)
128    cm(3,2)=d2gm(2,3)
129    cm(1,3)=d2gm(1,3)
130    cm(2,3)=d2gm(2,3)
131    cm(3,3)=d2gm(3,3)
132  elseif(rank==2)then
133    cm(1,1)=2*(dgm01(1,1)*dgm10(1,1)+gm(1,1)*d2gm(1,1))
134    cm(2,1)=-0.5d0*dgm01(2,2)*dgm10(1,1)+3*dgm01(1,2)*dgm10(1,2)&
135 &   -0.5d0*dgm01(1,1)*dgm10(2,2)-0.5d0*gm(2,2)*d2gm(1,1)+3*gm(1,2)&
136 &   *d2gm(1,2)-0.5d0*gm(1,1)*d2gm(2,2)
137    cm(3,1)=-0.5d0*dgm01(3,3)*dgm10(1,1)+3*dgm01(1,3)*dgm10(1,3)&
138 &   -0.5d0*dgm01(1,1)*dgm10(3,3)-0.5d0*gm(3,3)*d2gm(1,1)+3*gm(1,3)&
139 &   *d2gm(1,3)-0.5d0*gm(1,1)*d2gm(3,3)
140    cm(4,1)=-dgm01(2,3)*dgm10(1,1)+3*dgm01(1,3)*dgm10(1,2)+3*dgm01(1,2)&
141 &   *dgm10(1,3)-dgm01(1,1)*dgm10(2,3)-gm(2,3)*d2gm(1,1)+3*gm(1,3)&
142 &   *d2gm(1,2)+3*gm(1,2)*d2gm(1,3)-gm(1,1)*d2gm(2,3)
143    cm(5,1)=2*(dgm01(1,3)*dgm10(1,1)+dgm01(1,1)*dgm10(1,3)+gm(1,3)&
144 &   *d2gm(1,1)+gm(1,1)*d2gm(1,3))
145    cm(6,1)=2*(dgm01(1,2)*dgm10(1,1)+dgm01(1,1)*dgm10(1,2)+gm(1,2)&
146 &   *d2gm(1,1)+gm(1,1)*d2gm(1,2))
147    cm(1,2)=-0.5d0*dgm01(2,2)*dgm10(1,1)+3*dgm01(1,2)*dgm10(1,2)&
148 &   -0.5d0*dgm01(1,1)*dgm10(2,2)-0.5d0*gm(2,2)*d2gm(1,1)+3*gm(1,2)&
149 &   *d2gm(1,2)-0.5d0*gm(1,1)*d2gm(2,2)
150    cm(2,2)=2*(dgm01(2,2)*dgm10(2,2)+gm(2,2)*d2gm(2,2))
151    cm(3,2)=-0.5d0*dgm01(3,3)*dgm10(2,2)+3*dgm01(2,3)*dgm10(2,3)&
152 &   -0.5d0*dgm01(2,2)*dgm10(3,3)-0.5d0*gm(3,3)*d2gm(2,2)+3*gm(2,3)&
153 &   *d2gm(2,3)-0.5d0*gm(2,2)*d2gm(3,3)
154    cm(4,2)=2*(dgm01(2,3)*dgm10(2,2)+dgm01(2,2)*dgm10(2,3)+gm(2,3)&
155 &   *d2gm(2,2)+gm(2,2)*d2gm(2,3))
156    cm(5,2)=3*dgm01(2,3)*dgm10(1,2)-dgm01(2,2)*dgm10(1,3)-dgm01(1,3)&
157 &   *dgm10(2,2)+3*dgm01(1,2)*dgm10(2,3)+3*gm(2,3)*d2gm(1,2)-gm(2,2)&
158 &   *d2gm(1,3)-gm(1,3)*d2gm(2,2)+3*gm(1,2)*d2gm(2,3)
159    cm(6,2)=2*(dgm01(2,2)*dgm10(1,2)+dgm01(1,2)*dgm10(2,2)+gm(2,2)&
160 &   *d2gm(1,2)+gm(1,2)*d2gm(2,2))
161    cm(1,3)=-0.5d0*dgm01(3,3)*dgm10(1,1)+3*dgm01(1,3)*dgm10(1,3)&
162 &   -0.5d0*dgm01(1,1)*dgm10(3,3)-0.5d0*gm(3,3)*d2gm(1,1)+3*gm(1,3)&
163 &   *d2gm(1,3)-0.5d0*gm(1,1)*d2gm(3,3)
164    cm(2,3)=-0.5d0*dgm01(3,3)*dgm10(2,2)+3*dgm01(2,3)*dgm10(2,3)&
165 &   -0.5d0*dgm01(2,2)*dgm10(3,3)-0.5d0*gm(3,3)*d2gm(2,2)+3*gm(2,3)&
166 &   *d2gm(2,3)-0.5d0*gm(2,2)*d2gm(3,3)
167    cm(3,3)=2*(dgm01(3,3)*dgm10(3,3)+gm(3,3)*d2gm(3,3))
168    cm(4,3)=2*(dgm01(3,3)*dgm10(2,3)+dgm01(2,3)*dgm10(3,3)+gm(3,3)&
169 &   *d2gm(2,3)+gm(2,3)*d2gm(3,3))
170    cm(5,3)=2*(dgm01(3,3)*dgm10(1,3)+dgm01(1,3)*dgm10(3,3)+gm(3,3)&
171 &   *d2gm(1,3)+gm(1,3)*d2gm(3,3))
172    cm(6,3)=-dgm01(3,3)*dgm10(1,2)+3*dgm01(2,3)*dgm10(1,3)+3*dgm01(1,3)&
173 &   *dgm10(2,3)-dgm01(1,2)*dgm10(3,3)-gm(3,3)*d2gm(1,2)+3*gm(2,3)&
174 &   *d2gm(1,3)+3*gm(1,3)*d2gm(2,3)-gm(1,2)*d2gm(3,3)
175    cm(1,4)=-dgm01(2,3)*dgm10(1,1)+3*dgm01(1,3)*dgm10(1,2)+3*dgm01(1,2)&
176 &   *dgm10(1,3)-dgm01(1,1)*dgm10(2,3)-gm(2,3)*d2gm(1,1)+3*gm(1,3)&
177 &   *d2gm(1,2)+3*gm(1,2)*d2gm(1,3)-gm(1,1)*d2gm(2,3)
178    cm(2,4)=2*(dgm01(2,3)*dgm10(2,2)+dgm01(2,2)*dgm10(2,3)+gm(2,3)&
179 &   *d2gm(2,2)+gm(2,2)*d2gm(2,3))
180    cm(3,4)=2*(dgm01(3,3)*dgm10(2,3)+dgm01(2,3)*dgm10(3,3)+gm(3,3)&
181 &   *d2gm(2,3)+gm(2,3)*d2gm(3,3))
182    cm(4,4)=3*dgm01(3,3)*dgm10(2,2)+2*dgm01(2,3)*dgm10(2,3)+3*dgm01(2,2)&
183 &   *dgm10(3,3)+3*gm(3,3)*d2gm(2,2)+2*gm(2,3)*d2gm(2,3)+3*gm(2,2)&
184 &   *d2gm(3,3)
185    cm(5,4)=3*dgm01(3,3)*dgm10(1,2)+dgm01(2,3)*dgm10(1,3)+dgm01(1,3)&
186 &   *dgm10(2,3)+3*dgm01(1,2)*dgm10(3,3)+3*gm(3,3)*d2gm(1,2)+gm(2,3)&
187 &   *d2gm(1,3)+gm(1,3)*d2gm(2,3)+3*gm(1,2)*d2gm(3,3)
188    cm(6,4)=dgm01(2,3)*dgm10(1,2)+3*dgm01(2,2)*dgm10(1,3)+3*dgm01(1,3)&
189 &   *dgm10(2,2)+dgm01(1,2)*dgm10(2,3)+gm(2,3)*d2gm(1,2)+3*gm(2,2)&
190 &   *d2gm(1,3)+3*gm(1,3)*d2gm(2,2)+gm(1,2)*d2gm(2,3)
191    cm(1,5)=2*(dgm01(1,3)*dgm10(1,1)+dgm01(1,1)*dgm10(1,3)+gm(1,3)&
192 &   *d2gm(1,1)+gm(1,1)*d2gm(1,3))
193    cm(2,5)=3*dgm01(2,3)*dgm10(1,2)-dgm01(2,2)*dgm10(1,3)-dgm01(1,3)&
194 &   *dgm10(2,2)+3*dgm01(1,2)*dgm10(2,3)+3*gm(2,3)*d2gm(1,2)-gm(2,2)&
195 &   *d2gm(1,3)-gm(1,3)*d2gm(2,2)+3*gm(1,2)*d2gm(2,3)
196    cm(3,5)=2*(dgm01(3,3)*dgm10(1,3)+dgm01(1,3)*dgm10(3,3)+gm(3,3)&
197 &   *d2gm(1,3)+gm(1,3)*d2gm(3,3))
198    cm(4,5)=3*dgm01(3,3)*dgm10(1,2)+dgm01(2,3)*dgm10(1,3)+dgm01(1,3)&
199 &   *dgm10(2,3)+3*dgm01(1,2)*dgm10(3,3)+3*gm(3,3)*d2gm(1,2)+gm(2,3)&
200 &   *d2gm(1,3)+gm(1,3)*d2gm(2,3)+3*gm(1,2)*d2gm(3,3)
201    cm(5,5)=3*dgm01(3,3)*dgm10(1,1)+2*dgm01(1,3)*dgm10(1,3)+3*dgm01(1,1)&
202 &   *dgm10(3,3)+3*gm(3,3)*d2gm(1,1)+2*gm(1,3)*d2gm(1,3)+3*gm(1,1)&
203 &   *d2gm(3,3)
204    cm(6,5)=3*dgm01(2,3)*dgm10(1,1)+dgm01(1,3)*dgm10(1,2)+dgm01(1,2)&
205 &   *dgm10(1,3)+3*dgm01(1,1)*dgm10(2,3)+3*gm(2,3)*d2gm(1,1)+gm(1,3)&
206 &   *d2gm(1,2)+gm(1,2)*d2gm(1,3)+3*gm(1,1)*d2gm(2,3)
207    cm(1,6)=2*(dgm01(1,2)*dgm10(1,1)+dgm01(1,1)*dgm10(1,2)+gm(1,2)&
208 &   *d2gm(1,1)+gm(1,1)*d2gm(1,2))
209    cm(2,6)=2*(dgm01(2,2)*dgm10(1,2)+dgm01(1,2)*dgm10(2,2)+gm(2,2)&
210 &   *d2gm(1,2)+gm(1,2)*d2gm(2,2))
211    cm(3,6)=-dgm01(3,3)*dgm10(1,2)+3*dgm01(2,3)*dgm10(1,3)+3*dgm01(1,3)&
212 &   *dgm10(2,3)-dgm01(1,2)*dgm10(3,3)-gm(3,3)*d2gm(1,2)+3*gm(2,3)&
213 &   *d2gm(1,3)+3*gm(1,3)*d2gm(2,3)-gm(1,2)*d2gm(3,3)
214    cm(4,6)=dgm01(2,3)*dgm10(1,2)+3*dgm01(2,2)*dgm10(1,3)+3*dgm01(1,3)&
215 &   *dgm10(2,2)+dgm01(1,2)*dgm10(2,3)+gm(2,3)*d2gm(1,2)+3*gm(2,2)&
216 &   *d2gm(1,3)+3*gm(1,3)*d2gm(2,2)+gm(1,2)*d2gm(2,3)
217    cm(5,6)=3*dgm01(2,3)*dgm10(1,1)+dgm01(1,3)*dgm10(1,2)+dgm01(1,2)&
218 &   *dgm10(1,3)+3*dgm01(1,1)*dgm10(2,3)+3*gm(2,3)*d2gm(1,1)+gm(1,3)&
219 &   *d2gm(1,2)+gm(1,2)*d2gm(1,3)+3*gm(1,1)*d2gm(2,3)
220    cm(6,6)=3*dgm01(2,2)*dgm10(1,1)+2*dgm01(1,2)*dgm10(1,2)+3*dgm01(1,1)&
221 &   *dgm10(2,2)+3*gm(2,2)*d2gm(1,1)+2*gm(1,2)*d2gm(1,2)+3*gm(1,1)&
222 &   *d2gm(2,2)
223  elseif(rank==3)then
224    cm(1,1)=gm(1,1)*(6*dgm01(1,1)*dgm10(1,1)+3*gm(1,1)*d2gm(1,1))
225    cm(2,1)=4.5d0*gm(1,2)**2*d2gm(1,1)-3*gm(2,2)*(dgm01(1,1)*dgm10(1,1)&
226 &   +gm(1,1)*d2gm(1,1))+9*gm(1,2)*(dgm01(1,2)*dgm10(1,1)+dgm01(1,1)&
227 &   *dgm10(1,2)+gm(1,1)*d2gm(1,2))+gm(1,1)*(-3*dgm01(2,2)*dgm10(1,1)&
228 &   +9*dgm01(1,2)*dgm10(1,2)-3*dgm01(1,1)*dgm10(2,2)-1.5d0*gm(1,1)&
229 &   *d2gm(2,2))
230    cm(3,1)=4.5d0*gm(1,3)**2*d2gm(1,1)-3*gm(3,3)*(dgm01(1,1)*dgm10(1,1)&
231 &   +gm(1,1)*d2gm(1,1))+9*gm(1,3)*(dgm01(1,3)*dgm10(1,1)+dgm01(1,1)&
232 &   *dgm10(1,3)+gm(1,1)*d2gm(1,3))+gm(1,1)*(-3*dgm01(3,3)*dgm10(1,1)&
233 &   +9*dgm01(1,3)*dgm10(1,3)-3*dgm01(1,1)*dgm10(3,3)-1.5d0*gm(1,1)&
234 &   *d2gm(3,3))
235    cm(4,1)=9*gm(1,2)*dgm01(1,3)*dgm10(1,1)-6*gm(1,1)*dgm01(2,3)&
236 &   *dgm10(1,1)+9*gm(1,1)*dgm01(1,3)*dgm10(1,2)+9*gm(1,2)*dgm01(1,1)&
237 &   *dgm10(1,3)+9*gm(1,1)*dgm01(1,2)*dgm10(1,3)-6*gm(1,1)*dgm01(1,1)&
238 &   *dgm10(2,3)-6*gm(2,3)*(dgm01(1,1)*dgm10(1,1)+gm(1,1)*d2gm(1,1))&
239 &   +9*gm(1,3)*(dgm01(1,2)*dgm10(1,1)+dgm01(1,1)*dgm10(1,2)+gm(1,2)&
240 &   *d2gm(1,1)+gm(1,1)*d2gm(1,2))+9*gm(1,1)*gm(1,2)*d2gm(1,3)-3*gm(1,1)&
241 &   **2*d2gm(2,3)
242    cm(5,1)=6*gm(1,3)*(dgm01(1,1)*dgm10(1,1)+gm(1,1)*d2gm(1,1))&
243 &   +gm(1,1)*(6*dgm01(1,3)*dgm10(1,1)+6*dgm01(1,1)*dgm10(1,3)+3*gm(1,1)&
244 &   *d2gm(1,3))
245    cm(6,1)=6*gm(1,2)*(dgm01(1,1)*dgm10(1,1)+gm(1,1)*d2gm(1,1))&
246 &   +gm(1,1)*(6*dgm01(1,2)*dgm10(1,1)+6*dgm01(1,1)*dgm10(1,2)+3*gm(1,1)&
247 &   *d2gm(1,2))
248    cm(7,1)=-1.5d0*gm(1,1)*(dgm01(2,2)*dgm10(1,2)+dgm01(1,2)*dgm10(2,2))&
249 &   +7.5d0*gm(1,2)**2*d2gm(1,2)-1.5d0*gm(2,2)*(dgm01(1,2)*dgm10(1,1)&
250 &   +dgm01(1,1)*dgm10(1,2)+gm(1,2)*d2gm(1,1)+gm(1,1)*d2gm(1,2))+gm(1,2)&
251 &   *(-1.5d0*dgm01(2,2)*dgm10(1,1)+15*dgm01(1,2)*dgm10(1,2)-1.5d0*dgm01(1,1)&
252 &   *dgm10(2,2)-1.5d0*gm(1,1)*d2gm(2,2))
253    cm(8,1)=-3*gm(1,3)*dgm01(2,3)*dgm10(1,1)-1.5d0*gm(1,2)*dgm01(3,3)&
254 &   *dgm10(1,1)+15*gm(1,3)*dgm01(1,3)*dgm10(1,2)-1.5d0*gm(1,1)*dgm01(3,3)&
255 &   *dgm10(1,2)+15*gm(1,3)*dgm01(1,2)*dgm10(1,3)+15*gm(1,2)*dgm01(1,3)&
256 &   *dgm10(1,3)-3*gm(1,1)*dgm01(2,3)*dgm10(1,3)-3*gm(1,3)*dgm01(1,1)&
257 &   *dgm10(2,3)-3*gm(1,1)*dgm01(1,3)*dgm10(2,3)-1.5d0*gm(1,2)*dgm01(1,1)&
258 &   *dgm10(3,3)-1.5d0*gm(1,1)*dgm01(1,2)*dgm10(3,3)+7.5d0*gm(1,3)&
259 &   **2*d2gm(1,2)-1.5d0*gm(3,3)*(dgm01(1,2)*dgm10(1,1)+dgm01(1,1)&
260 &   *dgm10(1,2)+gm(1,2)*d2gm(1,1)+gm(1,1)*d2gm(1,2))+15*gm(1,2)*gm(1,3)&
261 &   *d2gm(1,3)-3*gm(2,3)*(dgm01(1,3)*dgm10(1,1)+dgm01(1,1)*dgm10(1,3)&
262 &   +gm(1,3)*d2gm(1,1)+gm(1,1)*d2gm(1,3))-3*gm(1,1)*gm(1,3)*d2gm(2,3)&
263 &   -1.5d0*gm(1,1)*gm(1,2)*d2gm(3,3)
264    cm(9,1)=-1.5d0*gm(1,3)*dgm01(2,2)*dgm10(1,1)-3*gm(1,2)*dgm01(2,3)&
265 &   *dgm10(1,1)+15*gm(1,3)*dgm01(1,2)*dgm10(1,2)+15*gm(1,2)*dgm01(1,3)&
266 &   *dgm10(1,2)-3*gm(1,1)*dgm01(2,3)*dgm10(1,2)+15*gm(1,2)*dgm01(1,2)&
267 &   *dgm10(1,3)-1.5d0*gm(1,1)*dgm01(2,2)*dgm10(1,3)-1.5d0*gm(1,3)&
268 &   *dgm01(1,1)*dgm10(2,2)-1.5d0*gm(1,1)*dgm01(1,3)*dgm10(2,2)-3*gm(1,2)&
269 &   *dgm01(1,1)*dgm10(2,3)-3*gm(1,1)*dgm01(1,2)*dgm10(2,3)+15*gm(1,2)&
270 &   *gm(1,3)*d2gm(1,2)-3*gm(2,3)*(dgm01(1,2)*dgm10(1,1)+dgm01(1,1)&
271 &   *dgm10(1,2)+gm(1,2)*d2gm(1,1)+gm(1,1)*d2gm(1,2))+7.5d0*gm(1,2)&
272 &   **2*d2gm(1,3)-1.5d0*gm(2,2)*(dgm01(1,3)*dgm10(1,1)+dgm01(1,1)&
273 &   *dgm10(1,3)+gm(1,3)*d2gm(1,1)+gm(1,1)*d2gm(1,3))-1.5d0*gm(1,1)&
274 &   *gm(1,3)*d2gm(2,2)-3*gm(1,1)*gm(1,2)*d2gm(2,3)
275    cm(10,1)=-1.5d0*gm(1,1)*(dgm01(3,3)*dgm10(1,3)+dgm01(1,3)*dgm10(3,3))&
276 &   +7.5d0*gm(1,3)**2*d2gm(1,3)-1.5d0*gm(3,3)*(dgm01(1,3)*dgm10(1,1)&
277 &   +dgm01(1,1)*dgm10(1,3)+gm(1,3)*d2gm(1,1)+gm(1,1)*d2gm(1,3))+gm(1,3)&
278 &   *(-1.5d0*dgm01(3,3)*dgm10(1,1)+15*dgm01(1,3)*dgm10(1,3)-1.5d0*dgm01(1,1)&
279 &   *dgm10(3,3)-1.5d0*gm(1,1)*d2gm(3,3))
280    cm(1,2)=4.5d0*gm(1,2)**2*d2gm(1,1)-3*gm(2,2)*(dgm01(1,1)*dgm10(1,1)&
281 &   +gm(1,1)*d2gm(1,1))+9*gm(1,2)*(dgm01(1,2)*dgm10(1,1)+dgm01(1,1)&
282 &   *dgm10(1,2)+gm(1,1)*d2gm(1,2))+gm(1,1)*(-3*dgm01(2,2)*dgm10(1,1)&
283 &   +9*dgm01(1,2)*dgm10(1,2)-3*dgm01(1,1)*dgm10(2,2)-1.5d0*gm(1,1)&
284 &   *d2gm(2,2))
285    cm(2,2)=12*gm(1,1)*dgm01(2,2)*dgm10(2,2)+6*gm(1,2)*(dgm01(2,2)&
286 &   *dgm10(1,2)+dgm01(1,2)*dgm10(2,2))+6*gm(2,2)**2*d2gm(1,1)+3*gm(1,2)&
287 &   **2*d2gm(2,2)+gm(2,2)*(12*dgm01(2,2)*dgm10(1,1)+6*dgm01(1,2)&
288 &   *dgm10(1,2)+12*dgm01(1,1)*dgm10(2,2)+6*gm(1,2)*d2gm(1,2)+12*gm(1,1)&
289 &   *d2gm(2,2))
290    cm(3,2)=-1.5d0*gm(2,2)*dgm01(3,3)*dgm10(1,1)+9*gm(1,3)*dgm01(2,3)&
291 &   *dgm10(1,2)-6*gm(1,2)*dgm01(3,3)*dgm10(1,2)-6*gm(2,2)*dgm01(1,3)&
292 &   *dgm10(1,3)-6*gm(1,3)*dgm01(2,2)*dgm10(1,3)+9*gm(1,2)*dgm01(2,3)&
293 &   *dgm10(1,3)-6*gm(1,3)*dgm01(1,3)*dgm10(2,2)-1.5d0*gm(1,1)*dgm01(3,3)&
294 &   *dgm10(2,2)+9*gm(1,3)*dgm01(1,2)*dgm10(2,3)+9*gm(1,2)*dgm01(1,3)&
295 &   *dgm10(2,3)+15*gm(1,1)*dgm01(2,3)*dgm10(2,3)-1.5d0*gm(2,2)*dgm01(1,1)&
296 &   *dgm10(3,3)-6*gm(1,2)*dgm01(1,2)*dgm10(3,3)-1.5d0*gm(1,1)*dgm01(2,2)&
297 &   *dgm10(3,3)+7.5d0*gm(2,3)**2*d2gm(1,1)-6*gm(1,3)*gm(2,2)*d2gm(1,3)&
298 &   -3*gm(1,3)**2*d2gm(2,2)+gm(3,3)*(-1.5d0*dgm01(2,2)*dgm10(1,1)&
299 &   -6*dgm01(1,2)*dgm10(1,2)-1.5d0*dgm01(1,1)*dgm10(2,2)-1.5d0*gm(2,2)&
300 &   *d2gm(1,1)-6*gm(1,2)*d2gm(1,2)-1.5d0*gm(1,1)*d2gm(2,2))+9*gm(1,2)&
301 &   *gm(1,3)*d2gm(2,3)+gm(2,3)*(15*dgm01(2,3)*dgm10(1,1)+9*dgm01(1,3)&
302 &   *dgm10(1,2)+9*dgm01(1,2)*dgm10(1,3)+15*dgm01(1,1)*dgm10(2,3)&
303 &   +9*gm(1,3)*d2gm(1,2)+9*gm(1,2)*d2gm(1,3)+15*gm(1,1)*d2gm(2,3))&
304 &   -3*gm(1,2)**2*d2gm(3,3)-1.5d0*gm(1,1)*gm(2,2)*d2gm(3,3)
305    cm(4,2)=3*gm(1,3)*dgm01(2,2)*dgm10(1,2)+6*gm(1,2)*dgm01(2,3)&
306 &   *dgm10(1,2)+3*gm(1,2)*dgm01(2,2)*dgm10(1,3)+3*gm(1,3)*dgm01(1,2)&
307 &   *dgm10(2,2)+3*gm(1,2)*dgm01(1,3)*dgm10(2,2)+12*gm(1,1)*dgm01(2,3)&
308 &   *dgm10(2,2)+6*gm(1,2)*dgm01(1,2)*dgm10(2,3)+12*gm(1,1)*dgm01(2,2)&
309 &   *dgm10(2,3)+3*gm(1,2)*gm(1,3)*d2gm(2,2)+gm(2,3)*(12*dgm01(2,2)&
310 &   *dgm10(1,1)+6*dgm01(1,2)*dgm10(1,2)+12*dgm01(1,1)*dgm10(2,2)&
311 &   +12*gm(2,2)*d2gm(1,1)+6*gm(1,2)*d2gm(1,2)+12*gm(1,1)*d2gm(2,2))&
312 &   +3*gm(1,2)**2*d2gm(2,3)+gm(2,2)*(12*dgm01(2,3)*dgm10(1,1)+3*dgm01(1,3)&
313 &   *dgm10(1,2)+3*dgm01(1,2)*dgm10(1,3)+12*dgm01(1,1)*dgm10(2,3)&
314 &   +3*gm(1,3)*d2gm(1,2)+3*gm(1,2)*d2gm(1,3)+12*gm(1,1)*d2gm(2,3))
315    cm(5,2)=-4.5d0*gm(1,3)*dgm01(2,2)*dgm10(1,1)+12*gm(1,2)*dgm01(2,3)&
316 &   *dgm10(1,1)+3*gm(1,3)*dgm01(1,2)*dgm10(1,2)+3*gm(1,2)*dgm01(1,3)&
317 &   *dgm10(1,2)+12*gm(1,1)*dgm01(2,3)*dgm10(1,2)+3*gm(1,2)*dgm01(1,2)&
318 &   *dgm10(1,3)-4.5d0*gm(1,1)*dgm01(2,2)*dgm10(1,3)-4.5d0*gm(1,3)&
319 &   *dgm01(1,1)*dgm10(2,2)-4.5d0*gm(1,1)*dgm01(1,3)*dgm10(2,2)+12*gm(1,2)&
320 &   *dgm01(1,1)*dgm10(2,3)+12*gm(1,1)*dgm01(1,2)*dgm10(2,3)+3*gm(1,2)&
321 &   *gm(1,3)*d2gm(1,2)+12*gm(2,3)*(dgm01(1,2)*dgm10(1,1)+dgm01(1,1)&
322 &   *dgm10(1,2)+gm(1,2)*d2gm(1,1)+gm(1,1)*d2gm(1,2))+1.5d0*gm(1,2)&
323 &   **2*d2gm(1,3)-4.5d0*gm(2,2)*(dgm01(1,3)*dgm10(1,1)+dgm01(1,1)&
324 &   *dgm10(1,3)+gm(1,3)*d2gm(1,1)+gm(1,1)*d2gm(1,3))-4.5d0*gm(1,1)&
325 &   *gm(1,3)*d2gm(2,2)+12*gm(1,1)*gm(1,2)*d2gm(2,3)
326    cm(6,2)=7.5d0*gm(1,1)*(dgm01(2,2)*dgm10(1,2)+dgm01(1,2)*dgm10(2,2))&
327 &   +4.5d0*gm(1,2)**2*d2gm(1,2)+7.5d0*gm(2,2)*(dgm01(1,2)*dgm10(1,1)&
328 &   +dgm01(1,1)*dgm10(1,2)+gm(1,2)*d2gm(1,1)+gm(1,1)*d2gm(1,2))+gm(1,2)&
329 &   *(7.5d0*dgm01(2,2)*dgm10(1,1)+9*dgm01(1,2)*dgm10(1,2)+7.5d0*dgm01(1,1)&
330 &   *dgm10(2,2)+7.5d0*gm(1,1)*d2gm(2,2))
331    cm(7,2)=6*gm(1,2)*dgm01(2,2)*dgm10(2,2)+3*gm(2,2)**2*d2gm(1,2)&
332 &   +6*gm(2,2)*(dgm01(2,2)*dgm10(1,2)+dgm01(1,2)*dgm10(2,2)+gm(1,2)&
333 &   *d2gm(2,2))
334    cm(8,2)=-4.5d0*gm(2,2)*dgm01(3,3)*dgm10(1,2)+12*gm(2,2)*dgm01(2,3)&
335 &   *dgm10(1,3)+12*gm(1,3)*dgm01(2,3)*dgm10(2,2)-4.5d0*gm(1,2)*dgm01(3,3)&
336 &   *dgm10(2,2)+12*gm(2,2)*dgm01(1,3)*dgm10(2,3)+12*gm(1,3)*dgm01(2,2)&
337 &   *dgm10(2,3)+3*gm(1,2)*dgm01(2,3)*dgm10(2,3)-4.5d0*gm(2,2)*dgm01(1,2)&
338 &   *dgm10(3,3)-4.5d0*gm(1,2)*dgm01(2,2)*dgm10(3,3)+1.5d0*gm(2,3)&
339 &   **2*d2gm(1,2)-4.5d0*gm(3,3)*(dgm01(2,2)*dgm10(1,2)+dgm01(1,2)&
340 &   *dgm10(2,2)+gm(2,2)*d2gm(1,2)+gm(1,2)*d2gm(2,2))+12*gm(1,3)*gm(2,2)&
341 &   *d2gm(2,3)+gm(2,3)*(3*dgm01(2,3)*dgm10(1,2)+12*dgm01(2,2)*dgm10(1,3)&
342 &   +12*dgm01(1,3)*dgm10(2,2)+3*dgm01(1,2)*dgm10(2,3)+12*gm(2,2)&
343 &   *d2gm(1,3)+12*gm(1,3)*d2gm(2,2)+3*gm(1,2)*d2gm(2,3))-4.5d0*gm(1,2)&
344 &   *gm(2,2)*d2gm(3,3)
345    cm(9,2)=12*gm(1,3)*dgm01(2,2)*dgm10(2,2)+3*gm(1,2)*dgm01(2,3)&
346 &   *dgm10(2,2)+3*gm(1,2)*dgm01(2,2)*dgm10(2,3)+6*gm(2,2)**2*d2gm(1,3)&
347 &   +3*gm(2,3)*(dgm01(2,2)*dgm10(1,2)+dgm01(1,2)*dgm10(2,2)+gm(2,2)&
348 &   *d2gm(1,2)+gm(1,2)*d2gm(2,2))+gm(2,2)*(3*dgm01(2,3)*dgm10(1,2)&
349 &   +12*dgm01(2,2)*dgm10(1,3)+12*dgm01(1,3)*dgm10(2,2)+3*dgm01(1,2)&
350 &   *dgm10(2,3)+12*gm(1,3)*d2gm(2,2)+3*gm(1,2)*d2gm(2,3))
351    cm(10,2)=-1.5d0*gm(2,2)*dgm01(3,3)*dgm10(1,3)-1.5d0*gm(1,3)&
352 &   *dgm01(3,3)*dgm10(2,2)+15*gm(1,3)*dgm01(2,3)*dgm10(2,3)-3*gm(1,2)&
353 &   *dgm01(3,3)*dgm10(2,3)-1.5d0*gm(2,2)*dgm01(1,3)*dgm10(3,3)-1.5d0*gm(1,3)&
354 &   *dgm01(2,2)*dgm10(3,3)-3*gm(1,2)*dgm01(2,3)*dgm10(3,3)+7.5d0*gm(2,3)&
355 &   **2*d2gm(1,3)+gm(3,3)*(-3*dgm01(2,3)*dgm10(1,2)-1.5d0*dgm01(2,2)&
356 &   *dgm10(1,3)-1.5d0*dgm01(1,3)*dgm10(2,2)-3*dgm01(1,2)*dgm10(2,3)&
357 &   -3*gm(2,3)*d2gm(1,2)-1.5d0*gm(2,2)*d2gm(1,3)-1.5d0*gm(1,3)*d2gm(2,2)&
358 &   -3*gm(1,2)*d2gm(2,3))-1.5d0*gm(1,3)*gm(2,2)*d2gm(3,3)+gm(2,3)&
359 &   *(-3*dgm01(3,3)*dgm10(1,2)+15*dgm01(2,3)*dgm10(1,3)+15*dgm01(1,3)&
360 &   *dgm10(2,3)-3*dgm01(1,2)*dgm10(3,3)+15*gm(1,3)*d2gm(2,3)-3*gm(1,2)&
361 &   *d2gm(3,3))
362    cm(1,3)=4.5d0*gm(1,3)**2*d2gm(1,1)-3*gm(3,3)*(dgm01(1,1)*dgm10(1,1)&
363 &   +gm(1,1)*d2gm(1,1))+9*gm(1,3)*(dgm01(1,3)*dgm10(1,1)+dgm01(1,1)&
364 &   *dgm10(1,3)+gm(1,1)*d2gm(1,3))+gm(1,1)*(-3*dgm01(3,3)*dgm10(1,1)&
365 &   +9*dgm01(1,3)*dgm10(1,3)-3*dgm01(1,1)*dgm10(3,3)-1.5d0*gm(1,1)&
366 &   *d2gm(3,3))
367    cm(2,3)=-1.5d0*gm(2,2)*dgm01(3,3)*dgm10(1,1)+9*gm(1,3)*dgm01(2,3)&
368 &   *dgm10(1,2)-6*gm(1,2)*dgm01(3,3)*dgm10(1,2)-6*gm(2,2)*dgm01(1,3)&
369 &   *dgm10(1,3)-6*gm(1,3)*dgm01(2,2)*dgm10(1,3)+9*gm(1,2)*dgm01(2,3)&
370 &   *dgm10(1,3)-6*gm(1,3)*dgm01(1,3)*dgm10(2,2)-1.5d0*gm(1,1)*dgm01(3,3)&
371 &   *dgm10(2,2)+9*gm(1,3)*dgm01(1,2)*dgm10(2,3)+9*gm(1,2)*dgm01(1,3)&
372 &   *dgm10(2,3)+15*gm(1,1)*dgm01(2,3)*dgm10(2,3)-1.5d0*gm(2,2)*dgm01(1,1)&
373 &   *dgm10(3,3)-6*gm(1,2)*dgm01(1,2)*dgm10(3,3)-1.5d0*gm(1,1)*dgm01(2,2)&
374 &   *dgm10(3,3)+7.5d0*gm(2,3)**2*d2gm(1,1)-6*gm(1,3)*gm(2,2)*d2gm(1,3)&
375 &   -3*gm(1,3)**2*d2gm(2,2)+gm(3,3)*(-1.5d0*dgm01(2,2)*dgm10(1,1)&
376 &   -6*dgm01(1,2)*dgm10(1,2)-1.5d0*dgm01(1,1)*dgm10(2,2)-1.5d0*gm(2,2)&
377 &   *d2gm(1,1)-6*gm(1,2)*d2gm(1,2)-1.5d0*gm(1,1)*d2gm(2,2))+9*gm(1,2)&
378 &   *gm(1,3)*d2gm(2,3)+gm(2,3)*(15*dgm01(2,3)*dgm10(1,1)+9*dgm01(1,3)&
379 &   *dgm10(1,2)+9*dgm01(1,2)*dgm10(1,3)+15*dgm01(1,1)*dgm10(2,3)&
380 &   +9*gm(1,3)*d2gm(1,2)+9*gm(1,2)*d2gm(1,3)+15*gm(1,1)*d2gm(2,3))&
381 &   -3*gm(1,2)**2*d2gm(3,3)-1.5d0*gm(1,1)*gm(2,2)*d2gm(3,3)
382    cm(3,3)=12*gm(1,1)*dgm01(3,3)*dgm10(3,3)+6*gm(1,3)*(dgm01(3,3)&
383 &   *dgm10(1,3)+dgm01(1,3)*dgm10(3,3))+6*gm(3,3)**2*d2gm(1,1)+3*gm(1,3)&
384 &   **2*d2gm(3,3)+gm(3,3)*(12*dgm01(3,3)*dgm10(1,1)+6*dgm01(1,3)&
385 &   *dgm10(1,3)+12*dgm01(1,1)*dgm10(3,3)+6*gm(1,3)*d2gm(1,3)+12*gm(1,1)&
386 &   *d2gm(3,3))
387    cm(4,3)=3*gm(1,3)*dgm01(3,3)*dgm10(1,2)+6*gm(1,3)*dgm01(2,3)&
388 &   *dgm10(1,3)+3*gm(1,2)*dgm01(3,3)*dgm10(1,3)+6*gm(1,3)*dgm01(1,3)&
389 &   *dgm10(2,3)+12*gm(1,1)*dgm01(3,3)*dgm10(2,3)+3*gm(1,3)*dgm01(1,2)&
390 &   *dgm10(3,3)+3*gm(1,2)*dgm01(1,3)*dgm10(3,3)+12*gm(1,1)*dgm01(2,3)&
391 &   *dgm10(3,3)+3*gm(1,3)**2*d2gm(2,3)+gm(3,3)*(12*dgm01(2,3)*dgm10(1,1)&
392 &   +3*dgm01(1,3)*dgm10(1,2)+3*dgm01(1,2)*dgm10(1,3)+12*dgm01(1,1)&
393 &   *dgm10(2,3)+12*gm(2,3)*d2gm(1,1)+3*gm(1,3)*d2gm(1,2)+3*gm(1,2)&
394 &   *d2gm(1,3)+12*gm(1,1)*d2gm(2,3))+3*gm(1,2)*gm(1,3)*d2gm(3,3)&
395 &   +gm(2,3)*(12*dgm01(3,3)*dgm10(1,1)+6*dgm01(1,3)*dgm10(1,3)+12*dgm01(1,1)&
396 &   *dgm10(3,3)+6*gm(1,3)*d2gm(1,3)+12*gm(1,1)*d2gm(3,3))
397    cm(5,3)=7.5d0*gm(1,1)*(dgm01(3,3)*dgm10(1,3)+dgm01(1,3)*dgm10(3,3))&
398 &   +4.5d0*gm(1,3)**2*d2gm(1,3)+7.5d0*gm(3,3)*(dgm01(1,3)*dgm10(1,1)&
399 &   +dgm01(1,1)*dgm10(1,3)+gm(1,3)*d2gm(1,1)+gm(1,1)*d2gm(1,3))+gm(1,3)&
400 &   *(7.5d0*dgm01(3,3)*dgm10(1,1)+9*dgm01(1,3)*dgm10(1,3)+7.5d0*dgm01(1,1)&
401 &   *dgm10(3,3)+7.5d0*gm(1,1)*d2gm(3,3))
402    cm(6,3)=12*gm(1,3)*dgm01(2,3)*dgm10(1,1)-4.5d0*gm(1,2)*dgm01(3,3)&
403 &   *dgm10(1,1)+3*gm(1,3)*dgm01(1,3)*dgm10(1,2)-4.5d0*gm(1,1)*dgm01(3,3)&
404 &   *dgm10(1,2)+3*gm(1,3)*dgm01(1,2)*dgm10(1,3)+3*gm(1,2)*dgm01(1,3)&
405 &   *dgm10(1,3)+12*gm(1,1)*dgm01(2,3)*dgm10(1,3)+12*gm(1,3)*dgm01(1,1)&
406 &   *dgm10(2,3)+12*gm(1,1)*dgm01(1,3)*dgm10(2,3)-4.5d0*gm(1,2)*dgm01(1,1)&
407 &   *dgm10(3,3)-4.5d0*gm(1,1)*dgm01(1,2)*dgm10(3,3)+1.5d0*gm(1,3)&
408 &   **2*d2gm(1,2)-4.5d0*gm(3,3)*(dgm01(1,2)*dgm10(1,1)+dgm01(1,1)&
409 &   *dgm10(1,2)+gm(1,2)*d2gm(1,1)+gm(1,1)*d2gm(1,2))+3*gm(1,2)*gm(1,3)&
410 &   *d2gm(1,3)+12*gm(2,3)*(dgm01(1,3)*dgm10(1,1)+dgm01(1,1)*dgm10(1,3)&
411 &   +gm(1,3)*d2gm(1,1)+gm(1,1)*d2gm(1,3))+12*gm(1,1)*gm(1,3)*d2gm(2,3)&
412 &   -4.5d0*gm(1,1)*gm(1,2)*d2gm(3,3)
413    cm(7,3)=-1.5d0*gm(2,2)*dgm01(3,3)*dgm10(1,2)-3*gm(2,2)*dgm01(2,3)&
414 &   *dgm10(1,3)-3*gm(1,3)*dgm01(2,3)*dgm10(2,2)-1.5d0*gm(1,2)*dgm01(3,3)&
415 &   *dgm10(2,2)-3*gm(2,2)*dgm01(1,3)*dgm10(2,3)-3*gm(1,3)*dgm01(2,2)&
416 &   *dgm10(2,3)+15*gm(1,2)*dgm01(2,3)*dgm10(2,3)-1.5d0*gm(2,2)*dgm01(1,2)&
417 &   *dgm10(3,3)-1.5d0*gm(1,2)*dgm01(2,2)*dgm10(3,3)+7.5d0*gm(2,3)&
418 &   **2*d2gm(1,2)-1.5d0*gm(3,3)*(dgm01(2,2)*dgm10(1,2)+dgm01(1,2)&
419 &   *dgm10(2,2)+gm(2,2)*d2gm(1,2)+gm(1,2)*d2gm(2,2))-3*gm(1,3)*gm(2,2)&
420 &   *d2gm(2,3)+gm(2,3)*(15*dgm01(2,3)*dgm10(1,2)-3*dgm01(2,2)*dgm10(1,3)&
421 &   -3*dgm01(1,3)*dgm10(2,2)+15*dgm01(1,2)*dgm10(2,3)-3*gm(2,2)*d2gm(1,3)&
422 &   -3*gm(1,3)*d2gm(2,2)+15*gm(1,2)*d2gm(2,3))-1.5d0*gm(1,2)*gm(2,2)&
423 &   *d2gm(3,3)
424    cm(8,3)=3*gm(1,3)*dgm01(3,3)*dgm10(2,3)+3*gm(1,3)*dgm01(2,3)&
425 &   *dgm10(3,3)+12*gm(1,2)*dgm01(3,3)*dgm10(3,3)+6*gm(3,3)**2*d2gm(1,2)&
426 &   +gm(3,3)*(12*dgm01(3,3)*dgm10(1,2)+3*dgm01(2,3)*dgm10(1,3)+3*dgm01(1,3)&
427 &   *dgm10(2,3)+12*dgm01(1,2)*dgm10(3,3)+3*gm(2,3)*d2gm(1,3)+3*gm(1,3)&
428 &   *d2gm(2,3)+12*gm(1,2)*d2gm(3,3))+3*gm(2,3)*(dgm01(3,3)*dgm10(1,3)&
429 &   +dgm01(1,3)*dgm10(3,3)+gm(1,3)*d2gm(3,3))
430    cm(9,3)=-4.5d0*gm(2,2)*dgm01(3,3)*dgm10(1,3)-4.5d0*gm(1,3)*dgm01(3,3)&
431 &   *dgm10(2,2)+3*gm(1,3)*dgm01(2,3)*dgm10(2,3)+12*gm(1,2)*dgm01(3,3)&
432 &   *dgm10(2,3)-4.5d0*gm(2,2)*dgm01(1,3)*dgm10(3,3)-4.5d0*gm(1,3)&
433 &   *dgm01(2,2)*dgm10(3,3)+12*gm(1,2)*dgm01(2,3)*dgm10(3,3)+1.5d0*gm(2,3)&
434 &   **2*d2gm(1,3)+gm(3,3)*(12*dgm01(2,3)*dgm10(1,2)-4.5d0*dgm01(2,2)&
435 &   *dgm10(1,3)-4.5d0*dgm01(1,3)*dgm10(2,2)+12*dgm01(1,2)*dgm10(2,3)&
436 &   +12*gm(2,3)*d2gm(1,2)-4.5d0*gm(2,2)*d2gm(1,3)-4.5d0*gm(1,3)*d2gm(2,2)&
437 &   +12*gm(1,2)*d2gm(2,3))-4.5d0*gm(1,3)*gm(2,2)*d2gm(3,3)+gm(2,3)&
438 &   *(12*dgm01(3,3)*dgm10(1,2)+3*dgm01(2,3)*dgm10(1,3)+3*dgm01(1,3)&
439 &   *dgm10(2,3)+12*dgm01(1,2)*dgm10(3,3)+3*gm(1,3)*d2gm(2,3)+12*gm(1,2)&
440 &   *d2gm(3,3))
441    cm(10,3)=6*gm(1,3)*dgm01(3,3)*dgm10(3,3)+3*gm(3,3)**2*d2gm(1,3)&
442 &   +6*gm(3,3)*(dgm01(3,3)*dgm10(1,3)+dgm01(1,3)*dgm10(3,3)+gm(1,3)&
443 &   *d2gm(3,3))
444    cm(1,4)=9*gm(1,2)*dgm01(1,3)*dgm10(1,1)-6*gm(1,1)*dgm01(2,3)&
445 &   *dgm10(1,1)+9*gm(1,1)*dgm01(1,3)*dgm10(1,2)+9*gm(1,2)*dgm01(1,1)&
446 &   *dgm10(1,3)+9*gm(1,1)*dgm01(1,2)*dgm10(1,3)-6*gm(1,1)*dgm01(1,1)&
447 &   *dgm10(2,3)-6*gm(2,3)*(dgm01(1,1)*dgm10(1,1)+gm(1,1)*d2gm(1,1))&
448 &   +9*gm(1,3)*(dgm01(1,2)*dgm10(1,1)+dgm01(1,1)*dgm10(1,2)+gm(1,2)&
449 &   *d2gm(1,1)+gm(1,1)*d2gm(1,2))+9*gm(1,1)*gm(1,2)*d2gm(1,3)-3*gm(1,1)&
450 &   **2*d2gm(2,3)
451    cm(2,4)=3*gm(1,3)*dgm01(2,2)*dgm10(1,2)+6*gm(1,2)*dgm01(2,3)&
452 &   *dgm10(1,2)+3*gm(1,2)*dgm01(2,2)*dgm10(1,3)+3*gm(1,3)*dgm01(1,2)&
453 &   *dgm10(2,2)+3*gm(1,2)*dgm01(1,3)*dgm10(2,2)+12*gm(1,1)*dgm01(2,3)&
454 &   *dgm10(2,2)+6*gm(1,2)*dgm01(1,2)*dgm10(2,3)+12*gm(1,1)*dgm01(2,2)&
455 &   *dgm10(2,3)+3*gm(1,2)*gm(1,3)*d2gm(2,2)+gm(2,3)*(12*dgm01(2,2)&
456 &   *dgm10(1,1)+6*dgm01(1,2)*dgm10(1,2)+12*dgm01(1,1)*dgm10(2,2)&
457 &   +12*gm(2,2)*d2gm(1,1)+6*gm(1,2)*d2gm(1,2)+12*gm(1,1)*d2gm(2,2))&
458 &   +3*gm(1,2)**2*d2gm(2,3)+gm(2,2)*(12*dgm01(2,3)*dgm10(1,1)+3*dgm01(1,3)&
459 &   *dgm10(1,2)+3*dgm01(1,2)*dgm10(1,3)+12*dgm01(1,1)*dgm10(2,3)&
460 &   +3*gm(1,3)*d2gm(1,2)+3*gm(1,2)*d2gm(1,3)+12*gm(1,1)*d2gm(2,3))
461    cm(3,4)=3*gm(1,3)*dgm01(3,3)*dgm10(1,2)+6*gm(1,3)*dgm01(2,3)&
462 &   *dgm10(1,3)+3*gm(1,2)*dgm01(3,3)*dgm10(1,3)+6*gm(1,3)*dgm01(1,3)&
463 &   *dgm10(2,3)+12*gm(1,1)*dgm01(3,3)*dgm10(2,3)+3*gm(1,3)*dgm01(1,2)&
464 &   *dgm10(3,3)+3*gm(1,2)*dgm01(1,3)*dgm10(3,3)+12*gm(1,1)*dgm01(2,3)&
465 &   *dgm10(3,3)+3*gm(1,3)**2*d2gm(2,3)+gm(3,3)*(12*dgm01(2,3)*dgm10(1,1)&
466 &   +3*dgm01(1,3)*dgm10(1,2)+3*dgm01(1,2)*dgm10(1,3)+12*dgm01(1,1)&
467 &   *dgm10(2,3)+12*gm(2,3)*d2gm(1,1)+3*gm(1,3)*d2gm(1,2)+3*gm(1,2)&
468 &   *d2gm(1,3)+12*gm(1,1)*d2gm(2,3))+3*gm(1,2)*gm(1,3)*d2gm(3,3)&
469 &   +gm(2,3)*(12*dgm01(3,3)*dgm10(1,1)+6*dgm01(1,3)*dgm10(1,3)+12*dgm01(1,1)&
470 &   *dgm10(3,3)+6*gm(1,3)*d2gm(1,3)+12*gm(1,1)*d2gm(3,3))
471    cm(4,4)=15*gm(2,2)*dgm01(3,3)*dgm10(1,1)-6*gm(1,3)*dgm01(2,3)&
472 &   *dgm10(1,2)+18*gm(1,2)*dgm01(3,3)*dgm10(1,2)+18*gm(2,2)*dgm01(1,3)&
473 &   *dgm10(1,3)+18*gm(1,3)*dgm01(2,2)*dgm10(1,3)-6*gm(1,2)*dgm01(2,3)&
474 &   *dgm10(1,3)+18*gm(1,3)*dgm01(1,3)*dgm10(2,2)+15*gm(1,1)*dgm01(3,3)&
475 &   *dgm10(2,2)-6*gm(1,3)*dgm01(1,2)*dgm10(2,3)-6*gm(1,2)*dgm01(1,3)&
476 &   *dgm10(2,3)+18*gm(1,1)*dgm01(2,3)*dgm10(2,3)+15*gm(2,2)*dgm01(1,1)&
477 &   *dgm10(3,3)+18*gm(1,2)*dgm01(1,2)*dgm10(3,3)+15*gm(1,1)*dgm01(2,2)&
478 &   *dgm10(3,3)+9*gm(2,3)**2*d2gm(1,1)+18*gm(1,3)*gm(2,2)*d2gm(1,3)&
479 &   +9*gm(1,3)**2*d2gm(2,2)+gm(3,3)*(15*dgm01(2,2)*dgm10(1,1)+18*dgm01(1,2)&
480 &   *dgm10(1,2)+15*dgm01(1,1)*dgm10(2,2)+15*gm(2,2)*d2gm(1,1)+18*gm(1,2)&
481 &   *d2gm(1,2)+15*gm(1,1)*d2gm(2,2))-6*gm(1,2)*gm(1,3)*d2gm(2,3)&
482 &   +gm(2,3)*(18*dgm01(2,3)*dgm10(1,1)-6*dgm01(1,3)*dgm10(1,2)-6*dgm01(1,2)&
483 &   *dgm10(1,3)+18*dgm01(1,1)*dgm10(2,3)-6*gm(1,3)*d2gm(1,2)-6*gm(1,2)&
484 &   *d2gm(1,3)+18*gm(1,1)*d2gm(2,3))+9*gm(1,2)**2*d2gm(3,3)+15*gm(1,1)&
485 &   *gm(2,2)*d2gm(3,3)
486    cm(5,4)=3*gm(1,3)*dgm01(2,3)*dgm10(1,1)+12*gm(1,2)*dgm01(3,3)&
487 &   *dgm10(1,1)+6*gm(1,3)*dgm01(1,3)*dgm10(1,2)+12*gm(1,1)*dgm01(3,3)&
488 &   *dgm10(1,2)+6*gm(1,3)*dgm01(1,2)*dgm10(1,3)+6*gm(1,2)*dgm01(1,3)&
489 &   *dgm10(1,3)+3*gm(1,1)*dgm01(2,3)*dgm10(1,3)+3*gm(1,3)*dgm01(1,1)&
490 &   *dgm10(2,3)+3*gm(1,1)*dgm01(1,3)*dgm10(2,3)+12*gm(1,2)*dgm01(1,1)&
491 &   *dgm10(3,3)+12*gm(1,1)*dgm01(1,2)*dgm10(3,3)+3*gm(1,3)**2*d2gm(1,2)&
492 &   +12*gm(3,3)*(dgm01(1,2)*dgm10(1,1)+dgm01(1,1)*dgm10(1,2)+gm(1,2)&
493 &   *d2gm(1,1)+gm(1,1)*d2gm(1,2))+6*gm(1,2)*gm(1,3)*d2gm(1,3)+3*gm(2,3)&
494 &   *(dgm01(1,3)*dgm10(1,1)+dgm01(1,1)*dgm10(1,3)+gm(1,3)*d2gm(1,1)&
495 &   +gm(1,1)*d2gm(1,3))+3*gm(1,1)*gm(1,3)*d2gm(2,3)+12*gm(1,1)*gm(1,2)&
496 &   *d2gm(3,3)
497    cm(6,4)=12*gm(1,3)*dgm01(2,2)*dgm10(1,1)+3*gm(1,2)*dgm01(2,3)&
498 &   *dgm10(1,1)+6*gm(1,3)*dgm01(1,2)*dgm10(1,2)+6*gm(1,2)*dgm01(1,3)&
499 &   *dgm10(1,2)+3*gm(1,1)*dgm01(2,3)*dgm10(1,2)+6*gm(1,2)*dgm01(1,2)&
500 &   *dgm10(1,3)+12*gm(1,1)*dgm01(2,2)*dgm10(1,3)+12*gm(1,3)*dgm01(1,1)&
501 &   *dgm10(2,2)+12*gm(1,1)*dgm01(1,3)*dgm10(2,2)+3*gm(1,2)*dgm01(1,1)&
502 &   *dgm10(2,3)+3*gm(1,1)*dgm01(1,2)*dgm10(2,3)+6*gm(1,2)*gm(1,3)&
503 &   *d2gm(1,2)+3*gm(2,3)*(dgm01(1,2)*dgm10(1,1)+dgm01(1,1)*dgm10(1,2)&
504 &   +gm(1,2)*d2gm(1,1)+gm(1,1)*d2gm(1,2))+3*gm(1,2)**2*d2gm(1,3)&
505 &   +12*gm(2,2)*(dgm01(1,3)*dgm10(1,1)+dgm01(1,1)*dgm10(1,3)+gm(1,3)&
506 &   *d2gm(1,1)+gm(1,1)*d2gm(1,3))+12*gm(1,1)*gm(1,3)*d2gm(2,2)+3*gm(1,1)&
507 &   *gm(1,2)*d2gm(2,3)
508    cm(7,4)=-6*gm(1,3)*dgm01(2,2)*dgm10(2,2)+9*gm(1,2)*dgm01(2,3)&
509 &   *dgm10(2,2)+9*gm(1,2)*dgm01(2,2)*dgm10(2,3)-3*gm(2,2)**2*d2gm(1,3)&
510 &   +9*gm(2,3)*(dgm01(2,2)*dgm10(1,2)+dgm01(1,2)*dgm10(2,2)+gm(2,2)&
511 &   *d2gm(1,2)+gm(1,2)*d2gm(2,2))+gm(2,2)*(9*dgm01(2,3)*dgm10(1,2)&
512 &   -6*dgm01(2,2)*dgm10(1,3)-6*dgm01(1,3)*dgm10(2,2)+9*dgm01(1,2)&
513 &   *dgm10(2,3)-6*gm(1,3)*d2gm(2,2)+9*gm(1,2)*d2gm(2,3))
514    cm(8,4)=12*gm(2,2)*dgm01(3,3)*dgm10(1,3)+12*gm(1,3)*dgm01(3,3)&
515 &   *dgm10(2,2)+6*gm(1,3)*dgm01(2,3)*dgm10(2,3)+3*gm(1,2)*dgm01(3,3)&
516 &   *dgm10(2,3)+12*gm(2,2)*dgm01(1,3)*dgm10(3,3)+12*gm(1,3)*dgm01(2,2)&
517 &   *dgm10(3,3)+3*gm(1,2)*dgm01(2,3)*dgm10(3,3)+3*gm(2,3)**2*d2gm(1,3)&
518 &   +gm(3,3)*(3*dgm01(2,3)*dgm10(1,2)+12*dgm01(2,2)*dgm10(1,3)+12*dgm01(1,3)&
519 &   *dgm10(2,2)+3*dgm01(1,2)*dgm10(2,3)+3*gm(2,3)*d2gm(1,2)+12*gm(2,2)&
520 &   *d2gm(1,3)+12*gm(1,3)*d2gm(2,2)+3*gm(1,2)*d2gm(2,3))+12*gm(1,3)&
521 &   *gm(2,2)*d2gm(3,3)+gm(2,3)*(3*dgm01(3,3)*dgm10(1,2)+6*dgm01(2,3)&
522 &   *dgm10(1,3)+6*dgm01(1,3)*dgm10(2,3)+3*dgm01(1,2)*dgm10(3,3)+6*gm(1,3)&
523 &   *d2gm(2,3)+3*gm(1,2)*d2gm(3,3))
524    cm(9,4)=12*gm(2,2)*dgm01(3,3)*dgm10(1,2)+3*gm(2,2)*dgm01(2,3)&
525 &   *dgm10(1,3)+3*gm(1,3)*dgm01(2,3)*dgm10(2,2)+12*gm(1,2)*dgm01(3,3)&
526 &   *dgm10(2,2)+3*gm(2,2)*dgm01(1,3)*dgm10(2,3)+3*gm(1,3)*dgm01(2,2)&
527 &   *dgm10(2,3)+6*gm(1,2)*dgm01(2,3)*dgm10(2,3)+12*gm(2,2)*dgm01(1,2)&
528 &   *dgm10(3,3)+12*gm(1,2)*dgm01(2,2)*dgm10(3,3)+3*gm(2,3)**2*d2gm(1,2)&
529 &   +12*gm(3,3)*(dgm01(2,2)*dgm10(1,2)+dgm01(1,2)*dgm10(2,2)+gm(2,2)&
530 &   *d2gm(1,2)+gm(1,2)*d2gm(2,2))+3*gm(1,3)*gm(2,2)*d2gm(2,3)+gm(2,3)&
531 &   *(6*dgm01(2,3)*dgm10(1,2)+3*dgm01(2,2)*dgm10(1,3)+3*dgm01(1,3)&
532 &   *dgm10(2,2)+6*dgm01(1,2)*dgm10(2,3)+3*gm(2,2)*d2gm(1,3)+3*gm(1,3)&
533 &   *d2gm(2,2)+6*gm(1,2)*d2gm(2,3))+12*gm(1,2)*gm(2,2)*d2gm(3,3)
534    cm(10,4)=9*gm(1,3)*dgm01(3,3)*dgm10(2,3)+9*gm(1,3)*dgm01(2,3)&
535 &   *dgm10(3,3)-6*gm(1,2)*dgm01(3,3)*dgm10(3,3)-3*gm(3,3)**2*d2gm(1,2)&
536 &   +gm(3,3)*(-6*dgm01(3,3)*dgm10(1,2)+9*dgm01(2,3)*dgm10(1,3)+9*dgm01(1,3)&
537 &   *dgm10(2,3)-6*dgm01(1,2)*dgm10(3,3)+9*gm(2,3)*d2gm(1,3)+9*gm(1,3)&
538 &   *d2gm(2,3)-6*gm(1,2)*d2gm(3,3))+9*gm(2,3)*(dgm01(3,3)*dgm10(1,3)&
539 &   +dgm01(1,3)*dgm10(3,3)+gm(1,3)*d2gm(3,3))
540    cm(1,5)=6*gm(1,3)*(dgm01(1,1)*dgm10(1,1)+gm(1,1)*d2gm(1,1))&
541 &   +gm(1,1)*(6*dgm01(1,3)*dgm10(1,1)+6*dgm01(1,1)*dgm10(1,3)+3*gm(1,1)&
542 &   *d2gm(1,3))
543    cm(2,5)=-4.5d0*gm(1,3)*dgm01(2,2)*dgm10(1,1)+12*gm(1,2)*dgm01(2,3)&
544 &   *dgm10(1,1)+3*gm(1,3)*dgm01(1,2)*dgm10(1,2)+3*gm(1,2)*dgm01(1,3)&
545 &   *dgm10(1,2)+12*gm(1,1)*dgm01(2,3)*dgm10(1,2)+3*gm(1,2)*dgm01(1,2)&
546 &   *dgm10(1,3)-4.5d0*gm(1,1)*dgm01(2,2)*dgm10(1,3)-4.5d0*gm(1,3)&
547 &   *dgm01(1,1)*dgm10(2,2)-4.5d0*gm(1,1)*dgm01(1,3)*dgm10(2,2)+12*gm(1,2)&
548 &   *dgm01(1,1)*dgm10(2,3)+12*gm(1,1)*dgm01(1,2)*dgm10(2,3)+3*gm(1,2)&
549 &   *gm(1,3)*d2gm(1,2)+12*gm(2,3)*(dgm01(1,2)*dgm10(1,1)+dgm01(1,1)&
550 &   *dgm10(1,2)+gm(1,2)*d2gm(1,1)+gm(1,1)*d2gm(1,2))+1.5d0*gm(1,2)&
551 &   **2*d2gm(1,3)-4.5d0*gm(2,2)*(dgm01(1,3)*dgm10(1,1)+dgm01(1,1)&
552 &   *dgm10(1,3)+gm(1,3)*d2gm(1,1)+gm(1,1)*d2gm(1,3))-4.5d0*gm(1,1)&
553 &   *gm(1,3)*d2gm(2,2)+12*gm(1,1)*gm(1,2)*d2gm(2,3)
554    cm(3,5)=7.5d0*gm(1,1)*(dgm01(3,3)*dgm10(1,3)+dgm01(1,3)*dgm10(3,3))&
555 &   +4.5d0*gm(1,3)**2*d2gm(1,3)+7.5d0*gm(3,3)*(dgm01(1,3)*dgm10(1,1)&
556 &   +dgm01(1,1)*dgm10(1,3)+gm(1,3)*d2gm(1,1)+gm(1,1)*d2gm(1,3))+gm(1,3)&
557 &   *(7.5d0*dgm01(3,3)*dgm10(1,1)+9*dgm01(1,3)*dgm10(1,3)+7.5d0*dgm01(1,1)&
558 &   *dgm10(3,3)+7.5d0*gm(1,1)*d2gm(3,3))
559    cm(4,5)=3*gm(1,3)*dgm01(2,3)*dgm10(1,1)+12*gm(1,2)*dgm01(3,3)&
560 &   *dgm10(1,1)+6*gm(1,3)*dgm01(1,3)*dgm10(1,2)+12*gm(1,1)*dgm01(3,3)&
561 &   *dgm10(1,2)+6*gm(1,3)*dgm01(1,2)*dgm10(1,3)+6*gm(1,2)*dgm01(1,3)&
562 &   *dgm10(1,3)+3*gm(1,1)*dgm01(2,3)*dgm10(1,3)+3*gm(1,3)*dgm01(1,1)&
563 &   *dgm10(2,3)+3*gm(1,1)*dgm01(1,3)*dgm10(2,3)+12*gm(1,2)*dgm01(1,1)&
564 &   *dgm10(3,3)+12*gm(1,1)*dgm01(1,2)*dgm10(3,3)+3*gm(1,3)**2*d2gm(1,2)&
565 &   +12*gm(3,3)*(dgm01(1,2)*dgm10(1,1)+dgm01(1,1)*dgm10(1,2)+gm(1,2)&
566 &   *d2gm(1,1)+gm(1,1)*d2gm(1,2))+6*gm(1,2)*gm(1,3)*d2gm(1,3)+3*gm(2,3)&
567 &   *(dgm01(1,3)*dgm10(1,1)+dgm01(1,1)*dgm10(1,3)+gm(1,3)*d2gm(1,1)&
568 &   +gm(1,1)*d2gm(1,3))+3*gm(1,1)*gm(1,3)*d2gm(2,3)+12*gm(1,1)*gm(1,2)&
569 &   *d2gm(3,3)
570    cm(5,5)=3*gm(1,3)**2*d2gm(1,1)+12*gm(3,3)*(dgm01(1,1)*dgm10(1,1)&
571 &   +gm(1,1)*d2gm(1,1))+6*gm(1,3)*(dgm01(1,3)*dgm10(1,1)+dgm01(1,1)&
572 &   *dgm10(1,3)+gm(1,1)*d2gm(1,3))+gm(1,1)*(12*dgm01(3,3)*dgm10(1,1)&
573 &   +6*dgm01(1,3)*dgm10(1,3)+12*dgm01(1,1)*dgm10(3,3)+6*gm(1,1)*d2gm(3,3))
574    cm(6,5)=3*gm(1,2)*dgm01(1,3)*dgm10(1,1)+12*gm(1,1)*dgm01(2,3)&
575 &   *dgm10(1,1)+3*gm(1,1)*dgm01(1,3)*dgm10(1,2)+3*gm(1,2)*dgm01(1,1)&
576 &   *dgm10(1,3)+3*gm(1,1)*dgm01(1,2)*dgm10(1,3)+12*gm(1,1)*dgm01(1,1)&
577 &   *dgm10(2,3)+12*gm(2,3)*(dgm01(1,1)*dgm10(1,1)+gm(1,1)*d2gm(1,1))&
578 &   +3*gm(1,3)*(dgm01(1,2)*dgm10(1,1)+dgm01(1,1)*dgm10(1,2)+gm(1,2)&
579 &   *d2gm(1,1)+gm(1,1)*d2gm(1,2))+3*gm(1,1)*gm(1,2)*d2gm(1,3)+6*gm(1,1)&
580 &   **2*d2gm(2,3)
581    cm(7,5)=(-36*gm(1,3)*dgm01(2,2)*dgm10(1,2)+180*gm(1,2)*dgm01(2,3)&
582 &   *dgm10(1,2)-36*gm(1,2)*dgm01(2,2)*dgm10(1,3)-36*gm(1,3)*dgm01(1,2)&
583 &   *dgm10(2,2)-36*gm(1,2)*dgm01(1,3)*dgm10(2,2)-18*gm(1,1)*dgm01(2,3)&
584 &   *dgm10(2,2)+180*gm(1,2)*dgm01(1,2)*dgm10(2,3)-18*gm(1,1)*dgm01(2,2)&
585 &   *dgm10(2,3)-36*gm(1,2)*gm(1,3)*d2gm(2,2)+gm(2,3)*(-18*dgm01(2,2)&
586 &   *dgm10(1,1)+180*dgm01(1,2)*dgm10(1,2)-18*dgm01(1,1)*dgm10(2,2)&
587 &   -18*gm(2,2)*d2gm(1,1)+180*gm(1,2)*d2gm(1,2)-18*gm(1,1)*d2gm(2,2))&
588 &   +90*gm(1,2)**2*d2gm(2,3)+gm(2,2)*(-18*dgm01(2,3)*dgm10(1,1)-36*dgm01(1,3)&
589 &   *dgm10(1,2)-36*dgm01(1,2)*dgm10(1,3)-18*dgm01(1,1)*dgm10(2,3)&
590 &   -36*gm(1,3)*d2gm(1,2)-36*gm(1,2)*d2gm(1,3)-18*gm(1,1)*d2gm(2,3)))&
591 &   /12.d0
592    cm(8,5)=12*gm(1,3)*dgm01(3,3)*dgm10(1,2)+3*gm(1,3)*dgm01(2,3)&
593 &   *dgm10(1,3)+12*gm(1,2)*dgm01(3,3)*dgm10(1,3)+3*gm(1,3)*dgm01(1,3)&
594 &   *dgm10(2,3)-4.5d0*gm(1,1)*dgm01(3,3)*dgm10(2,3)+12*gm(1,3)*dgm01(1,2)&
595 &   *dgm10(3,3)+12*gm(1,2)*dgm01(1,3)*dgm10(3,3)-4.5d0*gm(1,1)*dgm01(2,3)&
596 &   *dgm10(3,3)+1.5d0*gm(1,3)**2*d2gm(2,3)+gm(3,3)*(-4.5d0*dgm01(2,3)&
597 &   *dgm10(1,1)+12*dgm01(1,3)*dgm10(1,2)+12*dgm01(1,2)*dgm10(1,3)&
598 &   -4.5d0*dgm01(1,1)*dgm10(2,3)-4.5d0*gm(2,3)*d2gm(1,1)+12*gm(1,3)&
599 &   *d2gm(1,2)+12*gm(1,2)*d2gm(1,3)-4.5d0*gm(1,1)*d2gm(2,3))+12*gm(1,2)&
600 &   *gm(1,3)*d2gm(3,3)+gm(2,3)*(-4.5d0*dgm01(3,3)*dgm10(1,1)+3*dgm01(1,3)&
601 &   *dgm10(1,3)-4.5d0*dgm01(1,1)*dgm10(3,3)+3*gm(1,3)*d2gm(1,3)-4.5d0*gm(1,1)&
602 &   *d2gm(3,3))
603    cm(9,5)=-1.5d0*gm(2,2)*dgm01(3,3)*dgm10(1,1)+9*gm(1,3)*dgm01(2,3)&
604 &   *dgm10(1,2)+15*gm(1,2)*dgm01(3,3)*dgm10(1,2)-6*gm(2,2)*dgm01(1,3)&
605 &   *dgm10(1,3)-6*gm(1,3)*dgm01(2,2)*dgm10(1,3)+9*gm(1,2)*dgm01(2,3)&
606 &   *dgm10(1,3)-6*gm(1,3)*dgm01(1,3)*dgm10(2,2)-1.5d0*gm(1,1)*dgm01(3,3)&
607 &   *dgm10(2,2)+9*gm(1,3)*dgm01(1,2)*dgm10(2,3)+9*gm(1,2)*dgm01(1,3)&
608 &   *dgm10(2,3)-6*gm(1,1)*dgm01(2,3)*dgm10(2,3)-1.5d0*gm(2,2)*dgm01(1,1)&
609 &   *dgm10(3,3)+15*gm(1,2)*dgm01(1,2)*dgm10(3,3)-1.5d0*gm(1,1)*dgm01(2,2)&
610 &   *dgm10(3,3)-3*gm(2,3)**2*d2gm(1,1)-6*gm(1,3)*gm(2,2)*d2gm(1,3)&
611 &   -3*gm(1,3)**2*d2gm(2,2)+gm(3,3)*(-1.5d0*dgm01(2,2)*dgm10(1,1)&
612 &   +15*dgm01(1,2)*dgm10(1,2)-1.5d0*dgm01(1,1)*dgm10(2,2)-1.5d0*gm(2,2)&
613 &   *d2gm(1,1)+15*gm(1,2)*d2gm(1,2)-1.5d0*gm(1,1)*d2gm(2,2))+9*gm(1,2)&
614 &   *gm(1,3)*d2gm(2,3)+gm(2,3)*(-6*dgm01(2,3)*dgm10(1,1)+9*dgm01(1,3)&
615 &   *dgm10(1,2)+9*dgm01(1,2)*dgm10(1,3)-6*dgm01(1,1)*dgm10(2,3)+9*gm(1,3)&
616 &   *d2gm(1,2)+9*gm(1,2)*d2gm(1,3)-6*gm(1,1)*d2gm(2,3))+7.5d0*gm(1,2)&
617 &   **2*d2gm(3,3)-1.5d0*gm(1,1)*gm(2,2)*d2gm(3,3)
618    cm(10,5)=-3*gm(1,1)*dgm01(3,3)*dgm10(3,3)+9*gm(1,3)*(dgm01(3,3)&
619 &   *dgm10(1,3)+dgm01(1,3)*dgm10(3,3))-1.5d0*gm(3,3)**2*d2gm(1,1)&
620 &   +4.5d0*gm(1,3)**2*d2gm(3,3)+gm(3,3)*(-3*dgm01(3,3)*dgm10(1,1)&
621 &   +9*dgm01(1,3)*dgm10(1,3)-3*dgm01(1,1)*dgm10(3,3)+9*gm(1,3)*d2gm(1,3)&
622 &   -3*gm(1,1)*d2gm(3,3))
623    cm(1,6)=6*gm(1,2)*(dgm01(1,1)*dgm10(1,1)+gm(1,1)*d2gm(1,1))&
624 &   +gm(1,1)*(6*dgm01(1,2)*dgm10(1,1)+6*dgm01(1,1)*dgm10(1,2)+3*gm(1,1)&
625 &   *d2gm(1,2))
626    cm(2,6)=7.5d0*gm(1,1)*(dgm01(2,2)*dgm10(1,2)+dgm01(1,2)*dgm10(2,2))&
627 &   +4.5d0*gm(1,2)**2*d2gm(1,2)+7.5d0*gm(2,2)*(dgm01(1,2)*dgm10(1,1)&
628 &   +dgm01(1,1)*dgm10(1,2)+gm(1,2)*d2gm(1,1)+gm(1,1)*d2gm(1,2))+gm(1,2)&
629 &   *(7.5d0*dgm01(2,2)*dgm10(1,1)+9*dgm01(1,2)*dgm10(1,2)+7.5d0*dgm01(1,1)&
630 &   *dgm10(2,2)+7.5d0*gm(1,1)*d2gm(2,2))
631    cm(3,6)=12*gm(1,3)*dgm01(2,3)*dgm10(1,1)-4.5d0*gm(1,2)*dgm01(3,3)&
632 &   *dgm10(1,1)+3*gm(1,3)*dgm01(1,3)*dgm10(1,2)-4.5d0*gm(1,1)*dgm01(3,3)&
633 &   *dgm10(1,2)+3*gm(1,3)*dgm01(1,2)*dgm10(1,3)+3*gm(1,2)*dgm01(1,3)&
634 &   *dgm10(1,3)+12*gm(1,1)*dgm01(2,3)*dgm10(1,3)+12*gm(1,3)*dgm01(1,1)&
635 &   *dgm10(2,3)+12*gm(1,1)*dgm01(1,3)*dgm10(2,3)-4.5d0*gm(1,2)*dgm01(1,1)&
636 &   *dgm10(3,3)-4.5d0*gm(1,1)*dgm01(1,2)*dgm10(3,3)+1.5d0*gm(1,3)&
637 &   **2*d2gm(1,2)-4.5d0*gm(3,3)*(dgm01(1,2)*dgm10(1,1)+dgm01(1,1)&
638 &   *dgm10(1,2)+gm(1,2)*d2gm(1,1)+gm(1,1)*d2gm(1,2))+3*gm(1,2)*gm(1,3)&
639 &   *d2gm(1,3)+12*gm(2,3)*(dgm01(1,3)*dgm10(1,1)+dgm01(1,1)*dgm10(1,3)&
640 &   +gm(1,3)*d2gm(1,1)+gm(1,1)*d2gm(1,3))+12*gm(1,1)*gm(1,3)*d2gm(2,3)&
641 &   -4.5d0*gm(1,1)*gm(1,2)*d2gm(3,3)
642    cm(4,6)=12*gm(1,3)*dgm01(2,2)*dgm10(1,1)+3*gm(1,2)*dgm01(2,3)&
643 &   *dgm10(1,1)+6*gm(1,3)*dgm01(1,2)*dgm10(1,2)+6*gm(1,2)*dgm01(1,3)&
644 &   *dgm10(1,2)+3*gm(1,1)*dgm01(2,3)*dgm10(1,2)+6*gm(1,2)*dgm01(1,2)&
645 &   *dgm10(1,3)+12*gm(1,1)*dgm01(2,2)*dgm10(1,3)+12*gm(1,3)*dgm01(1,1)&
646 &   *dgm10(2,2)+12*gm(1,1)*dgm01(1,3)*dgm10(2,2)+3*gm(1,2)*dgm01(1,1)&
647 &   *dgm10(2,3)+3*gm(1,1)*dgm01(1,2)*dgm10(2,3)+6*gm(1,2)*gm(1,3)&
648 &   *d2gm(1,2)+3*gm(2,3)*(dgm01(1,2)*dgm10(1,1)+dgm01(1,1)*dgm10(1,2)&
649 &   +gm(1,2)*d2gm(1,1)+gm(1,1)*d2gm(1,2))+3*gm(1,2)**2*d2gm(1,3)&
650 &   +12*gm(2,2)*(dgm01(1,3)*dgm10(1,1)+dgm01(1,1)*dgm10(1,3)+gm(1,3)&
651 &   *d2gm(1,1)+gm(1,1)*d2gm(1,3))+12*gm(1,1)*gm(1,3)*d2gm(2,2)+3*gm(1,1)&
652 &   *gm(1,2)*d2gm(2,3)
653    cm(5,6)=3*gm(1,2)*dgm01(1,3)*dgm10(1,1)+12*gm(1,1)*dgm01(2,3)&
654 &   *dgm10(1,1)+3*gm(1,1)*dgm01(1,3)*dgm10(1,2)+3*gm(1,2)*dgm01(1,1)&
655 &   *dgm10(1,3)+3*gm(1,1)*dgm01(1,2)*dgm10(1,3)+12*gm(1,1)*dgm01(1,1)&
656 &   *dgm10(2,3)+12*gm(2,3)*(dgm01(1,1)*dgm10(1,1)+gm(1,1)*d2gm(1,1))&
657 &   +3*gm(1,3)*(dgm01(1,2)*dgm10(1,1)+dgm01(1,1)*dgm10(1,2)+gm(1,2)&
658 &   *d2gm(1,1)+gm(1,1)*d2gm(1,2))+3*gm(1,1)*gm(1,2)*d2gm(1,3)+6*gm(1,1)&
659 &   **2*d2gm(2,3)
660    cm(6,6)=3*gm(1,2)**2*d2gm(1,1)+12*gm(2,2)*(dgm01(1,1)*dgm10(1,1)&
661 &   +gm(1,1)*d2gm(1,1))+6*gm(1,2)*(dgm01(1,2)*dgm10(1,1)+dgm01(1,1)&
662 &   *dgm10(1,2)+gm(1,1)*d2gm(1,2))+gm(1,1)*(12*dgm01(2,2)*dgm10(1,1)&
663 &   +6*dgm01(1,2)*dgm10(1,2)+12*dgm01(1,1)*dgm10(2,2)+6*gm(1,1)*d2gm(2,2))
664    cm(7,6)=-3*gm(1,1)*dgm01(2,2)*dgm10(2,2)+9*gm(1,2)*(dgm01(2,2)&
665 &   *dgm10(1,2)+dgm01(1,2)*dgm10(2,2))-1.5d0*gm(2,2)**2*d2gm(1,1)&
666 &   +4.5d0*gm(1,2)**2*d2gm(2,2)+gm(2,2)*(-3*dgm01(2,2)*dgm10(1,1)&
667 &   +9*dgm01(1,2)*dgm10(1,2)-3*dgm01(1,1)*dgm10(2,2)+9*gm(1,2)*d2gm(1,2)&
668 &   -3*gm(1,1)*d2gm(2,2))
669    cm(8,6)=-1.5d0*gm(2,2)*dgm01(3,3)*dgm10(1,1)+9*gm(1,3)*dgm01(2,3)&
670 &   *dgm10(1,2)-6*gm(1,2)*dgm01(3,3)*dgm10(1,2)+15*gm(2,2)*dgm01(1,3)&
671 &   *dgm10(1,3)+15*gm(1,3)*dgm01(2,2)*dgm10(1,3)+9*gm(1,2)*dgm01(2,3)&
672 &   *dgm10(1,3)+15*gm(1,3)*dgm01(1,3)*dgm10(2,2)-1.5d0*gm(1,1)*dgm01(3,3)&
673 &   *dgm10(2,2)+9*gm(1,3)*dgm01(1,2)*dgm10(2,3)+9*gm(1,2)*dgm01(1,3)&
674 &   *dgm10(2,3)-6*gm(1,1)*dgm01(2,3)*dgm10(2,3)-1.5d0*gm(2,2)*dgm01(1,1)&
675 &   *dgm10(3,3)-6*gm(1,2)*dgm01(1,2)*dgm10(3,3)-1.5d0*gm(1,1)*dgm01(2,2)&
676 &   *dgm10(3,3)-3*gm(2,3)**2*d2gm(1,1)+15*gm(1,3)*gm(2,2)*d2gm(1,3)&
677 &   +7.5d0*gm(1,3)**2*d2gm(2,2)+gm(3,3)*(-1.5d0*dgm01(2,2)*dgm10(1,1)&
678 &   -6*dgm01(1,2)*dgm10(1,2)-1.5d0*dgm01(1,1)*dgm10(2,2)-1.5d0*gm(2,2)&
679 &   *d2gm(1,1)-6*gm(1,2)*d2gm(1,2)-1.5d0*gm(1,1)*d2gm(2,2))+9*gm(1,2)&
680 &   *gm(1,3)*d2gm(2,3)+gm(2,3)*(-6*dgm01(2,3)*dgm10(1,1)+9*dgm01(1,3)&
681 &   *dgm10(1,2)+9*dgm01(1,2)*dgm10(1,3)-6*dgm01(1,1)*dgm10(2,3)+9*gm(1,3)&
682 &   *d2gm(1,2)+9*gm(1,2)*d2gm(1,3)-6*gm(1,1)*d2gm(2,3))-3*gm(1,2)&
683 &   **2*d2gm(3,3)-1.5d0*gm(1,1)*gm(2,2)*d2gm(3,3)
684    cm(9,6)=12*gm(1,3)*dgm01(2,2)*dgm10(1,2)+3*gm(1,2)*dgm01(2,3)&
685 &   *dgm10(1,2)+12*gm(1,2)*dgm01(2,2)*dgm10(1,3)+12*gm(1,3)*dgm01(1,2)&
686 &   *dgm10(2,2)+12*gm(1,2)*dgm01(1,3)*dgm10(2,2)-4.5d0*gm(1,1)*dgm01(2,3)&
687 &   *dgm10(2,2)+3*gm(1,2)*dgm01(1,2)*dgm10(2,3)-4.5d0*gm(1,1)*dgm01(2,2)&
688 &   *dgm10(2,3)+12*gm(1,2)*gm(1,3)*d2gm(2,2)+gm(2,3)*(-4.5d0*dgm01(2,2)&
689 &   *dgm10(1,1)+3*dgm01(1,2)*dgm10(1,2)-4.5d0*dgm01(1,1)*dgm10(2,2)&
690 &   -4.5d0*gm(2,2)*d2gm(1,1)+3*gm(1,2)*d2gm(1,2)-4.5d0*gm(1,1)*d2gm(2,2))&
691 &   +1.5d0*gm(1,2)**2*d2gm(2,3)+gm(2,2)*(-4.5d0*dgm01(2,3)*dgm10(1,1)&
692 &   +12*dgm01(1,3)*dgm10(1,2)+12*dgm01(1,2)*dgm10(1,3)-4.5d0*dgm01(1,1)&
693 &   *dgm10(2,3)+12*gm(1,3)*d2gm(1,2)+12*gm(1,2)*d2gm(1,3)-4.5d0*gm(1,1)&
694 &   *d2gm(2,3))
695    cm(10,6)=(-36*gm(1,3)*dgm01(3,3)*dgm10(1,2)+180*gm(1,3)*dgm01(2,3)&
696 &   *dgm10(1,3)-36*gm(1,2)*dgm01(3,3)*dgm10(1,3)+180*gm(1,3)*dgm01(1,3)&
697 &   *dgm10(2,3)-18*gm(1,1)*dgm01(3,3)*dgm10(2,3)-36*gm(1,3)*dgm01(1,2)&
698 &   *dgm10(3,3)-36*gm(1,2)*dgm01(1,3)*dgm10(3,3)-18*gm(1,1)*dgm01(2,3)&
699 &   *dgm10(3,3)+90*gm(1,3)**2*d2gm(2,3)+gm(3,3)*(-18*dgm01(2,3)*dgm10(1,1)&
700 &   -36*dgm01(1,3)*dgm10(1,2)-36*dgm01(1,2)*dgm10(1,3)-18*dgm01(1,1)&
701 &   *dgm10(2,3)-18*gm(2,3)*d2gm(1,1)-36*gm(1,3)*d2gm(1,2)-36*gm(1,2)&
702 &   *d2gm(1,3)-18*gm(1,1)*d2gm(2,3))-36*gm(1,2)*gm(1,3)*d2gm(3,3)&
703 &   +gm(2,3)*(-18*dgm01(3,3)*dgm10(1,1)+180*dgm01(1,3)*dgm10(1,3)&
704 &   -18*dgm01(1,1)*dgm10(3,3)+180*gm(1,3)*d2gm(1,3)-18*gm(1,1)*d2gm(3,3)))&
705 &   /12.d0
706    cm(1,7)=-1.5d0*gm(1,1)*(dgm01(2,2)*dgm10(1,2)+dgm01(1,2)*dgm10(2,2))&
707 &   +7.5d0*gm(1,2)**2*d2gm(1,2)-1.5d0*gm(2,2)*(dgm01(1,2)*dgm10(1,1)&
708 &   +dgm01(1,1)*dgm10(1,2)+gm(1,2)*d2gm(1,1)+gm(1,1)*d2gm(1,2))+gm(1,2)&
709 &   *(-1.5d0*dgm01(2,2)*dgm10(1,1)+15*dgm01(1,2)*dgm10(1,2)-1.5d0*dgm01(1,1)&
710 &   *dgm10(2,2)-1.5d0*gm(1,1)*d2gm(2,2))
711    cm(2,7)=6*gm(1,2)*dgm01(2,2)*dgm10(2,2)+3*gm(2,2)**2*d2gm(1,2)&
712 &   +6*gm(2,2)*(dgm01(2,2)*dgm10(1,2)+dgm01(1,2)*dgm10(2,2)+gm(1,2)&
713 &   *d2gm(2,2))
714    cm(3,7)=-1.5d0*gm(2,2)*dgm01(3,3)*dgm10(1,2)-3*gm(2,2)*dgm01(2,3)&
715 &   *dgm10(1,3)-3*gm(1,3)*dgm01(2,3)*dgm10(2,2)-1.5d0*gm(1,2)*dgm01(3,3)&
716 &   *dgm10(2,2)-3*gm(2,2)*dgm01(1,3)*dgm10(2,3)-3*gm(1,3)*dgm01(2,2)&
717 &   *dgm10(2,3)+15*gm(1,2)*dgm01(2,3)*dgm10(2,3)-1.5d0*gm(2,2)*dgm01(1,2)&
718 &   *dgm10(3,3)-1.5d0*gm(1,2)*dgm01(2,2)*dgm10(3,3)+7.5d0*gm(2,3)&
719 &   **2*d2gm(1,2)-1.5d0*gm(3,3)*(dgm01(2,2)*dgm10(1,2)+dgm01(1,2)&
720 &   *dgm10(2,2)+gm(2,2)*d2gm(1,2)+gm(1,2)*d2gm(2,2))-3*gm(1,3)*gm(2,2)&
721 &   *d2gm(2,3)+gm(2,3)*(15*dgm01(2,3)*dgm10(1,2)-3*dgm01(2,2)*dgm10(1,3)&
722 &   -3*dgm01(1,3)*dgm10(2,2)+15*dgm01(1,2)*dgm10(2,3)-3*gm(2,2)*d2gm(1,3)&
723 &   -3*gm(1,3)*d2gm(2,2)+15*gm(1,2)*d2gm(2,3))-1.5d0*gm(1,2)*gm(2,2)&
724 &   *d2gm(3,3)
725    cm(4,7)=-6*gm(1,3)*dgm01(2,2)*dgm10(2,2)+9*gm(1,2)*dgm01(2,3)&
726 &   *dgm10(2,2)+9*gm(1,2)*dgm01(2,2)*dgm10(2,3)-3*gm(2,2)**2*d2gm(1,3)&
727 &   +9*gm(2,3)*(dgm01(2,2)*dgm10(1,2)+dgm01(1,2)*dgm10(2,2)+gm(2,2)&
728 &   *d2gm(1,2)+gm(1,2)*d2gm(2,2))+gm(2,2)*(9*dgm01(2,3)*dgm10(1,2)&
729 &   -6*dgm01(2,2)*dgm10(1,3)-6*dgm01(1,3)*dgm10(2,2)+9*dgm01(1,2)&
730 &   *dgm10(2,3)-6*gm(1,3)*d2gm(2,2)+9*gm(1,2)*d2gm(2,3))
731    cm(5,7)=-3*gm(1,3)*dgm01(2,2)*dgm10(1,2)+15*gm(1,2)*dgm01(2,3)&
732 &   *dgm10(1,2)-3*gm(1,2)*dgm01(2,2)*dgm10(1,3)-3*gm(1,3)*dgm01(1,2)&
733 &   *dgm10(2,2)-3*gm(1,2)*dgm01(1,3)*dgm10(2,2)-1.5d0*gm(1,1)*dgm01(2,3)&
734 &   *dgm10(2,2)+15*gm(1,2)*dgm01(1,2)*dgm10(2,3)-1.5d0*gm(1,1)*dgm01(2,2)&
735 &   *dgm10(2,3)-3*gm(1,2)*gm(1,3)*d2gm(2,2)+gm(2,3)*(-1.5d0*dgm01(2,2)&
736 &   *dgm10(1,1)+15*dgm01(1,2)*dgm10(1,2)-1.5d0*dgm01(1,1)*dgm10(2,2)&
737 &   -1.5d0*gm(2,2)*d2gm(1,1)+15*gm(1,2)*d2gm(1,2)-1.5d0*gm(1,1)*d2gm(2,2))&
738 &   +7.5d0*gm(1,2)**2*d2gm(2,3)+gm(2,2)*(-1.5d0*dgm01(2,3)*dgm10(1,1)&
739 &   -3*dgm01(1,3)*dgm10(1,2)-3*dgm01(1,2)*dgm10(1,3)-1.5d0*dgm01(1,1)&
740 &   *dgm10(2,3)-3*gm(1,3)*d2gm(1,2)-3*gm(1,2)*d2gm(1,3)-1.5d0*gm(1,1)&
741 &   *d2gm(2,3))
742    cm(6,7)=-3*gm(1,1)*dgm01(2,2)*dgm10(2,2)+9*gm(1,2)*(dgm01(2,2)&
743 &   *dgm10(1,2)+dgm01(1,2)*dgm10(2,2))-1.5d0*gm(2,2)**2*d2gm(1,1)&
744 &   +4.5d0*gm(1,2)**2*d2gm(2,2)+gm(2,2)*(-3*dgm01(2,2)*dgm10(1,1)&
745 &   +9*dgm01(1,2)*dgm10(1,2)-3*dgm01(1,1)*dgm10(2,2)+9*gm(1,2)*d2gm(1,2)&
746 &   -3*gm(1,1)*d2gm(2,2))
747    cm(7,7)=gm(2,2)*(6*dgm01(2,2)*dgm10(2,2)+3*gm(2,2)*d2gm(2,2))
748    cm(8,7)=4.5d0*gm(2,3)**2*d2gm(2,2)-3*gm(3,3)*(dgm01(2,2)*dgm10(2,2)&
749 &   +gm(2,2)*d2gm(2,2))+9*gm(2,3)*(dgm01(2,3)*dgm10(2,2)+dgm01(2,2)&
750 &   *dgm10(2,3)+gm(2,2)*d2gm(2,3))+gm(2,2)*(-3*dgm01(3,3)*dgm10(2,2)&
751 &   +9*dgm01(2,3)*dgm10(2,3)-3*dgm01(2,2)*dgm10(3,3)-1.5d0*gm(2,2)&
752 &   *d2gm(3,3))
753    cm(9,7)=6*gm(2,3)*(dgm01(2,2)*dgm10(2,2)+gm(2,2)*d2gm(2,2))&
754 &   +gm(2,2)*(6*dgm01(2,3)*dgm10(2,2)+6*dgm01(2,2)*dgm10(2,3)+3*gm(2,2)&
755 &   *d2gm(2,3))
756    cm(10,7)=-1.5d0*gm(2,2)*(dgm01(3,3)*dgm10(2,3)+dgm01(2,3)*dgm10(3,3))&
757 &   +7.5d0*gm(2,3)**2*d2gm(2,3)-1.5d0*gm(3,3)*(dgm01(2,3)*dgm10(2,2)&
758 &   +dgm01(2,2)*dgm10(2,3)+gm(2,3)*d2gm(2,2)+gm(2,2)*d2gm(2,3))+gm(2,3)&
759 &   *(-1.5d0*dgm01(3,3)*dgm10(2,2)+15*dgm01(2,3)*dgm10(2,3)-1.5d0*dgm01(2,2)&
760 &   *dgm10(3,3)-1.5d0*gm(2,2)*d2gm(3,3))
761    cm(1,8)=-3*gm(1,3)*dgm01(2,3)*dgm10(1,1)-1.5d0*gm(1,2)*dgm01(3,3)&
762 &   *dgm10(1,1)+15*gm(1,3)*dgm01(1,3)*dgm10(1,2)-1.5d0*gm(1,1)*dgm01(3,3)&
763 &   *dgm10(1,2)+15*gm(1,3)*dgm01(1,2)*dgm10(1,3)+15*gm(1,2)*dgm01(1,3)&
764 &   *dgm10(1,3)-3*gm(1,1)*dgm01(2,3)*dgm10(1,3)-3*gm(1,3)*dgm01(1,1)&
765 &   *dgm10(2,3)-3*gm(1,1)*dgm01(1,3)*dgm10(2,3)-1.5d0*gm(1,2)*dgm01(1,1)&
766 &   *dgm10(3,3)-1.5d0*gm(1,1)*dgm01(1,2)*dgm10(3,3)+7.5d0*gm(1,3)&
767 &   **2*d2gm(1,2)-1.5d0*gm(3,3)*(dgm01(1,2)*dgm10(1,1)+dgm01(1,1)&
768 &   *dgm10(1,2)+gm(1,2)*d2gm(1,1)+gm(1,1)*d2gm(1,2))+15*gm(1,2)*gm(1,3)&
769 &   *d2gm(1,3)-3*gm(2,3)*(dgm01(1,3)*dgm10(1,1)+dgm01(1,1)*dgm10(1,3)&
770 &   +gm(1,3)*d2gm(1,1)+gm(1,1)*d2gm(1,3))-3*gm(1,1)*gm(1,3)*d2gm(2,3)&
771 &   -1.5d0*gm(1,1)*gm(1,2)*d2gm(3,3)
772    cm(2,8)=-4.5d0*gm(2,2)*dgm01(3,3)*dgm10(1,2)+12*gm(2,2)*dgm01(2,3)&
773 &   *dgm10(1,3)+12*gm(1,3)*dgm01(2,3)*dgm10(2,2)-4.5d0*gm(1,2)*dgm01(3,3)&
774 &   *dgm10(2,2)+12*gm(2,2)*dgm01(1,3)*dgm10(2,3)+12*gm(1,3)*dgm01(2,2)&
775 &   *dgm10(2,3)+3*gm(1,2)*dgm01(2,3)*dgm10(2,3)-4.5d0*gm(2,2)*dgm01(1,2)&
776 &   *dgm10(3,3)-4.5d0*gm(1,2)*dgm01(2,2)*dgm10(3,3)+1.5d0*gm(2,3)&
777 &   **2*d2gm(1,2)-4.5d0*gm(3,3)*(dgm01(2,2)*dgm10(1,2)+dgm01(1,2)&
778 &   *dgm10(2,2)+gm(2,2)*d2gm(1,2)+gm(1,2)*d2gm(2,2))+12*gm(1,3)*gm(2,2)&
779 &   *d2gm(2,3)+gm(2,3)*(3*dgm01(2,3)*dgm10(1,2)+12*dgm01(2,2)*dgm10(1,3)&
780 &   +12*dgm01(1,3)*dgm10(2,2)+3*dgm01(1,2)*dgm10(2,3)+12*gm(2,2)&
781 &   *d2gm(1,3)+12*gm(1,3)*d2gm(2,2)+3*gm(1,2)*d2gm(2,3))-4.5d0*gm(1,2)&
782 &   *gm(2,2)*d2gm(3,3)
783    cm(3,8)=3*gm(1,3)*dgm01(3,3)*dgm10(2,3)+3*gm(1,3)*dgm01(2,3)&
784 &   *dgm10(3,3)+12*gm(1,2)*dgm01(3,3)*dgm10(3,3)+6*gm(3,3)**2*d2gm(1,2)&
785 &   +gm(3,3)*(12*dgm01(3,3)*dgm10(1,2)+3*dgm01(2,3)*dgm10(1,3)+3*dgm01(1,3)&
786 &   *dgm10(2,3)+12*dgm01(1,2)*dgm10(3,3)+3*gm(2,3)*d2gm(1,3)+3*gm(1,3)&
787 &   *d2gm(2,3)+12*gm(1,2)*d2gm(3,3))+3*gm(2,3)*(dgm01(3,3)*dgm10(1,3)&
788 &   +dgm01(1,3)*dgm10(3,3)+gm(1,3)*d2gm(3,3))
789    cm(4,8)=12*gm(2,2)*dgm01(3,3)*dgm10(1,3)+12*gm(1,3)*dgm01(3,3)&
790 &   *dgm10(2,2)+6*gm(1,3)*dgm01(2,3)*dgm10(2,3)+3*gm(1,2)*dgm01(3,3)&
791 &   *dgm10(2,3)+12*gm(2,2)*dgm01(1,3)*dgm10(3,3)+12*gm(1,3)*dgm01(2,2)&
792 &   *dgm10(3,3)+3*gm(1,2)*dgm01(2,3)*dgm10(3,3)+3*gm(2,3)**2*d2gm(1,3)&
793 &   +gm(3,3)*(3*dgm01(2,3)*dgm10(1,2)+12*dgm01(2,2)*dgm10(1,3)+12*dgm01(1,3)&
794 &   *dgm10(2,2)+3*dgm01(1,2)*dgm10(2,3)+3*gm(2,3)*d2gm(1,2)+12*gm(2,2)&
795 &   *d2gm(1,3)+12*gm(1,3)*d2gm(2,2)+3*gm(1,2)*d2gm(2,3))+12*gm(1,3)&
796 &   *gm(2,2)*d2gm(3,3)+gm(2,3)*(3*dgm01(3,3)*dgm10(1,2)+6*dgm01(2,3)&
797 &   *dgm10(1,3)+6*dgm01(1,3)*dgm10(2,3)+3*dgm01(1,2)*dgm10(3,3)+6*gm(1,3)&
798 &   *d2gm(2,3)+3*gm(1,2)*d2gm(3,3))
799    cm(5,8)=12*gm(1,3)*dgm01(3,3)*dgm10(1,2)+3*gm(1,3)*dgm01(2,3)&
800 &   *dgm10(1,3)+12*gm(1,2)*dgm01(3,3)*dgm10(1,3)+3*gm(1,3)*dgm01(1,3)&
801 &   *dgm10(2,3)-4.5d0*gm(1,1)*dgm01(3,3)*dgm10(2,3)+12*gm(1,3)*dgm01(1,2)&
802 &   *dgm10(3,3)+12*gm(1,2)*dgm01(1,3)*dgm10(3,3)-4.5d0*gm(1,1)*dgm01(2,3)&
803 &   *dgm10(3,3)+1.5d0*gm(1,3)**2*d2gm(2,3)+gm(3,3)*(-4.5d0*dgm01(2,3)&
804 &   *dgm10(1,1)+12*dgm01(1,3)*dgm10(1,2)+12*dgm01(1,2)*dgm10(1,3)&
805 &   -4.5d0*dgm01(1,1)*dgm10(2,3)-4.5d0*gm(2,3)*d2gm(1,1)+12*gm(1,3)&
806 &   *d2gm(1,2)+12*gm(1,2)*d2gm(1,3)-4.5d0*gm(1,1)*d2gm(2,3))+12*gm(1,2)&
807 &   *gm(1,3)*d2gm(3,3)+gm(2,3)*(-4.5d0*dgm01(3,3)*dgm10(1,1)+3*dgm01(1,3)&
808 &   *dgm10(1,3)-4.5d0*dgm01(1,1)*dgm10(3,3)+3*gm(1,3)*d2gm(1,3)-4.5d0*gm(1,1)&
809 &   *d2gm(3,3))
810    cm(6,8)=-1.5d0*gm(2,2)*dgm01(3,3)*dgm10(1,1)+9*gm(1,3)*dgm01(2,3)&
811 &   *dgm10(1,2)-6*gm(1,2)*dgm01(3,3)*dgm10(1,2)+15*gm(2,2)*dgm01(1,3)&
812 &   *dgm10(1,3)+15*gm(1,3)*dgm01(2,2)*dgm10(1,3)+9*gm(1,2)*dgm01(2,3)&
813 &   *dgm10(1,3)+15*gm(1,3)*dgm01(1,3)*dgm10(2,2)-1.5d0*gm(1,1)*dgm01(3,3)&
814 &   *dgm10(2,2)+9*gm(1,3)*dgm01(1,2)*dgm10(2,3)+9*gm(1,2)*dgm01(1,3)&
815 &   *dgm10(2,3)-6*gm(1,1)*dgm01(2,3)*dgm10(2,3)-1.5d0*gm(2,2)*dgm01(1,1)&
816 &   *dgm10(3,3)-6*gm(1,2)*dgm01(1,2)*dgm10(3,3)-1.5d0*gm(1,1)*dgm01(2,2)&
817 &   *dgm10(3,3)-3*gm(2,3)**2*d2gm(1,1)+15*gm(1,3)*gm(2,2)*d2gm(1,3)&
818 &   +7.5d0*gm(1,3)**2*d2gm(2,2)+gm(3,3)*(-1.5d0*dgm01(2,2)*dgm10(1,1)&
819 &   -6*dgm01(1,2)*dgm10(1,2)-1.5d0*dgm01(1,1)*dgm10(2,2)-1.5d0*gm(2,2)&
820 &   *d2gm(1,1)-6*gm(1,2)*d2gm(1,2)-1.5d0*gm(1,1)*d2gm(2,2))+9*gm(1,2)&
821 &   *gm(1,3)*d2gm(2,3)+gm(2,3)*(-6*dgm01(2,3)*dgm10(1,1)+9*dgm01(1,3)&
822 &   *dgm10(1,2)+9*dgm01(1,2)*dgm10(1,3)-6*dgm01(1,1)*dgm10(2,3)+9*gm(1,3)&
823 &   *d2gm(1,2)+9*gm(1,2)*d2gm(1,3)-6*gm(1,1)*d2gm(2,3))-3*gm(1,2)&
824 &   **2*d2gm(3,3)-1.5d0*gm(1,1)*gm(2,2)*d2gm(3,3)
825    cm(7,8)=4.5d0*gm(2,3)**2*d2gm(2,2)-3*gm(3,3)*(dgm01(2,2)*dgm10(2,2)&
826 &   +gm(2,2)*d2gm(2,2))+9*gm(2,3)*(dgm01(2,3)*dgm10(2,2)+dgm01(2,2)&
827 &   *dgm10(2,3)+gm(2,2)*d2gm(2,3))+gm(2,2)*(-3*dgm01(3,3)*dgm10(2,2)&
828 &   +9*dgm01(2,3)*dgm10(2,3)-3*dgm01(2,2)*dgm10(3,3)-1.5d0*gm(2,2)&
829 &   *d2gm(3,3))
830    cm(8,8)=12*gm(2,2)*dgm01(3,3)*dgm10(3,3)+6*gm(2,3)*(dgm01(3,3)&
831 &   *dgm10(2,3)+dgm01(2,3)*dgm10(3,3))+6*gm(3,3)**2*d2gm(2,2)+3*gm(2,3)&
832 &   **2*d2gm(3,3)+gm(3,3)*(12*dgm01(3,3)*dgm10(2,2)+6*dgm01(2,3)&
833 &   *dgm10(2,3)+12*dgm01(2,2)*dgm10(3,3)+6*gm(2,3)*d2gm(2,3)+12*gm(2,2)&
834 &   *d2gm(3,3))
835    cm(9,8)=7.5d0*gm(2,2)*(dgm01(3,3)*dgm10(2,3)+dgm01(2,3)*dgm10(3,3))&
836 &   +4.5d0*gm(2,3)**2*d2gm(2,3)+7.5d0*gm(3,3)*(dgm01(2,3)*dgm10(2,2)&
837 &   +dgm01(2,2)*dgm10(2,3)+gm(2,3)*d2gm(2,2)+gm(2,2)*d2gm(2,3))+gm(2,3)&
838 &   *(7.5d0*dgm01(3,3)*dgm10(2,2)+9*dgm01(2,3)*dgm10(2,3)+7.5d0*dgm01(2,2)&
839 &   *dgm10(3,3)+7.5d0*gm(2,2)*d2gm(3,3))
840    cm(10,8)=6*gm(2,3)*dgm01(3,3)*dgm10(3,3)+3*gm(3,3)**2*d2gm(2,3)&
841 &   +6*gm(3,3)*(dgm01(3,3)*dgm10(2,3)+dgm01(2,3)*dgm10(3,3)+gm(2,3)&
842 &   *d2gm(3,3))
843    cm(1,9)=(-18*gm(1,3)*dgm01(2,2)*dgm10(1,1)-36*gm(1,2)*dgm01(2,3)&
844 &   *dgm10(1,1)+180*gm(1,3)*dgm01(1,2)*dgm10(1,2)+180*gm(1,2)*dgm01(1,3)&
845 &   *dgm10(1,2)-36*gm(1,1)*dgm01(2,3)*dgm10(1,2)+180*gm(1,2)*dgm01(1,2)&
846 &   *dgm10(1,3)-18*gm(1,1)*dgm01(2,2)*dgm10(1,3)-18*gm(1,3)*dgm01(1,1)&
847 &   *dgm10(2,2)-18*gm(1,1)*dgm01(1,3)*dgm10(2,2)-36*gm(1,2)*dgm01(1,1)&
848 &   *dgm10(2,3)-36*gm(1,1)*dgm01(1,2)*dgm10(2,3)+180*gm(1,2)*gm(1,3)&
849 &   *d2gm(1,2)-36*gm(2,3)*(dgm01(1,2)*dgm10(1,1)+dgm01(1,1)*dgm10(1,2)&
850 &   +gm(1,2)*d2gm(1,1)+gm(1,1)*d2gm(1,2))+90*gm(1,2)**2*d2gm(1,3)&
851 &   -18*gm(2,2)*(dgm01(1,3)*dgm10(1,1)+dgm01(1,1)*dgm10(1,3)+gm(1,3)&
852 &   *d2gm(1,1)+gm(1,1)*d2gm(1,3))-18*gm(1,1)*gm(1,3)*d2gm(2,2)-36*gm(1,1)&
853 &   *gm(1,2)*d2gm(2,3))/12.d0
854    cm(2,9)=12*gm(1,3)*dgm01(2,2)*dgm10(2,2)+3*gm(1,2)*dgm01(2,3)&
855 &   *dgm10(2,2)+3*gm(1,2)*dgm01(2,2)*dgm10(2,3)+6*gm(2,2)**2*d2gm(1,3)&
856 &   +3*gm(2,3)*(dgm01(2,2)*dgm10(1,2)+dgm01(1,2)*dgm10(2,2)+gm(2,2)&
857 &   *d2gm(1,2)+gm(1,2)*d2gm(2,2))+gm(2,2)*(3*dgm01(2,3)*dgm10(1,2)&
858 &   +12*dgm01(2,2)*dgm10(1,3)+12*dgm01(1,3)*dgm10(2,2)+3*dgm01(1,2)&
859 &   *dgm10(2,3)+12*gm(1,3)*d2gm(2,2)+3*gm(1,2)*d2gm(2,3))
860    cm(3,9)=-4.5d0*gm(2,2)*dgm01(3,3)*dgm10(1,3)-4.5d0*gm(1,3)*dgm01(3,3)&
861 &   *dgm10(2,2)+3*gm(1,3)*dgm01(2,3)*dgm10(2,3)+12*gm(1,2)*dgm01(3,3)&
862 &   *dgm10(2,3)-4.5d0*gm(2,2)*dgm01(1,3)*dgm10(3,3)-4.5d0*gm(1,3)&
863 &   *dgm01(2,2)*dgm10(3,3)+12*gm(1,2)*dgm01(2,3)*dgm10(3,3)+1.5d0*gm(2,3)&
864 &   **2*d2gm(1,3)+gm(3,3)*(12*dgm01(2,3)*dgm10(1,2)-4.5d0*dgm01(2,2)&
865 &   *dgm10(1,3)-4.5d0*dgm01(1,3)*dgm10(2,2)+12*dgm01(1,2)*dgm10(2,3)&
866 &   +12*gm(2,3)*d2gm(1,2)-4.5d0*gm(2,2)*d2gm(1,3)-4.5d0*gm(1,3)*d2gm(2,2)&
867 &   +12*gm(1,2)*d2gm(2,3))-4.5d0*gm(1,3)*gm(2,2)*d2gm(3,3)+gm(2,3)&
868 &   *(12*dgm01(3,3)*dgm10(1,2)+3*dgm01(2,3)*dgm10(1,3)+3*dgm01(1,3)&
869 &   *dgm10(2,3)+12*dgm01(1,2)*dgm10(3,3)+3*gm(1,3)*d2gm(2,3)+12*gm(1,2)&
870 &   *d2gm(3,3))
871    cm(4,9)=12*gm(2,2)*dgm01(3,3)*dgm10(1,2)+3*gm(2,2)*dgm01(2,3)&
872 &   *dgm10(1,3)+3*gm(1,3)*dgm01(2,3)*dgm10(2,2)+12*gm(1,2)*dgm01(3,3)&
873 &   *dgm10(2,2)+3*gm(2,2)*dgm01(1,3)*dgm10(2,3)+3*gm(1,3)*dgm01(2,2)&
874 &   *dgm10(2,3)+6*gm(1,2)*dgm01(2,3)*dgm10(2,3)+12*gm(2,2)*dgm01(1,2)&
875 &   *dgm10(3,3)+12*gm(1,2)*dgm01(2,2)*dgm10(3,3)+3*gm(2,3)**2*d2gm(1,2)&
876 &   +12*gm(3,3)*(dgm01(2,2)*dgm10(1,2)+dgm01(1,2)*dgm10(2,2)+gm(2,2)&
877 &   *d2gm(1,2)+gm(1,2)*d2gm(2,2))+3*gm(1,3)*gm(2,2)*d2gm(2,3)+gm(2,3)&
878 &   *(6*dgm01(2,3)*dgm10(1,2)+3*dgm01(2,2)*dgm10(1,3)+3*dgm01(1,3)&
879 &   *dgm10(2,2)+6*dgm01(1,2)*dgm10(2,3)+3*gm(2,2)*d2gm(1,3)+3*gm(1,3)&
880 &   *d2gm(2,2)+6*gm(1,2)*d2gm(2,3))+12*gm(1,2)*gm(2,2)*d2gm(3,3)
881    cm(5,9)=-1.5d0*gm(2,2)*dgm01(3,3)*dgm10(1,1)+9*gm(1,3)*dgm01(2,3)&
882 &   *dgm10(1,2)+15*gm(1,2)*dgm01(3,3)*dgm10(1,2)-6*gm(2,2)*dgm01(1,3)&
883 &   *dgm10(1,3)-6*gm(1,3)*dgm01(2,2)*dgm10(1,3)+9*gm(1,2)*dgm01(2,3)&
884 &   *dgm10(1,3)-6*gm(1,3)*dgm01(1,3)*dgm10(2,2)-1.5d0*gm(1,1)*dgm01(3,3)&
885 &   *dgm10(2,2)+9*gm(1,3)*dgm01(1,2)*dgm10(2,3)+9*gm(1,2)*dgm01(1,3)&
886 &   *dgm10(2,3)-6*gm(1,1)*dgm01(2,3)*dgm10(2,3)-1.5d0*gm(2,2)*dgm01(1,1)&
887 &   *dgm10(3,3)+15*gm(1,2)*dgm01(1,2)*dgm10(3,3)-1.5d0*gm(1,1)*dgm01(2,2)&
888 &   *dgm10(3,3)-3*gm(2,3)**2*d2gm(1,1)-6*gm(1,3)*gm(2,2)*d2gm(1,3)&
889 &   -3*gm(1,3)**2*d2gm(2,2)+gm(3,3)*(-1.5d0*dgm01(2,2)*dgm10(1,1)&
890 &   +15*dgm01(1,2)*dgm10(1,2)-1.5d0*dgm01(1,1)*dgm10(2,2)-1.5d0*gm(2,2)&
891 &   *d2gm(1,1)+15*gm(1,2)*d2gm(1,2)-1.5d0*gm(1,1)*d2gm(2,2))+9*gm(1,2)&
892 &   *gm(1,3)*d2gm(2,3)+gm(2,3)*(-6*dgm01(2,3)*dgm10(1,1)+9*dgm01(1,3)&
893 &   *dgm10(1,2)+9*dgm01(1,2)*dgm10(1,3)-6*dgm01(1,1)*dgm10(2,3)+9*gm(1,3)&
894 &   *d2gm(1,2)+9*gm(1,2)*d2gm(1,3)-6*gm(1,1)*d2gm(2,3))+7.5d0*gm(1,2)&
895 &   **2*d2gm(3,3)-1.5d0*gm(1,1)*gm(2,2)*d2gm(3,3)
896    cm(6,9)=12*gm(1,3)*dgm01(2,2)*dgm10(1,2)+3*gm(1,2)*dgm01(2,3)&
897 &   *dgm10(1,2)+12*gm(1,2)*dgm01(2,2)*dgm10(1,3)+12*gm(1,3)*dgm01(1,2)&
898 &   *dgm10(2,2)+12*gm(1,2)*dgm01(1,3)*dgm10(2,2)-4.5d0*gm(1,1)*dgm01(2,3)&
899 &   *dgm10(2,2)+3*gm(1,2)*dgm01(1,2)*dgm10(2,3)-4.5d0*gm(1,1)*dgm01(2,2)&
900 &   *dgm10(2,3)+12*gm(1,2)*gm(1,3)*d2gm(2,2)+gm(2,3)*(-4.5d0*dgm01(2,2)&
901 &   *dgm10(1,1)+3*dgm01(1,2)*dgm10(1,2)-4.5d0*dgm01(1,1)*dgm10(2,2)&
902 &   -4.5d0*gm(2,2)*d2gm(1,1)+3*gm(1,2)*d2gm(1,2)-4.5d0*gm(1,1)*d2gm(2,2))&
903 &   +1.5d0*gm(1,2)**2*d2gm(2,3)+gm(2,2)*(-4.5d0*dgm01(2,3)*dgm10(1,1)&
904 &   +12*dgm01(1,3)*dgm10(1,2)+12*dgm01(1,2)*dgm10(1,3)-4.5d0*dgm01(1,1)&
905 &   *dgm10(2,3)+12*gm(1,3)*d2gm(1,2)+12*gm(1,2)*d2gm(1,3)-4.5d0*gm(1,1)&
906 &   *d2gm(2,3))
907    cm(7,9)=6*gm(2,3)*(dgm01(2,2)*dgm10(2,2)+gm(2,2)*d2gm(2,2))&
908 &   +gm(2,2)*(6*dgm01(2,3)*dgm10(2,2)+6*dgm01(2,2)*dgm10(2,3)+3*gm(2,2)&
909 &   *d2gm(2,3))
910    cm(8,9)=7.5d0*gm(2,2)*(dgm01(3,3)*dgm10(2,3)+dgm01(2,3)*dgm10(3,3))&
911 &   +4.5d0*gm(2,3)**2*d2gm(2,3)+7.5d0*gm(3,3)*(dgm01(2,3)*dgm10(2,2)&
912 &   +dgm01(2,2)*dgm10(2,3)+gm(2,3)*d2gm(2,2)+gm(2,2)*d2gm(2,3))+gm(2,3)&
913 &   *(7.5d0*dgm01(3,3)*dgm10(2,2)+9*dgm01(2,3)*dgm10(2,3)+7.5d0*dgm01(2,2)&
914 &   *dgm10(3,3)+7.5d0*gm(2,2)*d2gm(3,3))
915    cm(9,9)=3*gm(2,3)**2*d2gm(2,2)+12*gm(3,3)*(dgm01(2,2)*dgm10(2,2)&
916 &   +gm(2,2)*d2gm(2,2))+6*gm(2,3)*(dgm01(2,3)*dgm10(2,2)+dgm01(2,2)&
917 &   *dgm10(2,3)+gm(2,2)*d2gm(2,3))+gm(2,2)*(12*dgm01(3,3)*dgm10(2,2)&
918 &   +6*dgm01(2,3)*dgm10(2,3)+12*dgm01(2,2)*dgm10(3,3)+6*gm(2,2)*d2gm(3,3))
919    cm(10,9)=-3*gm(2,2)*dgm01(3,3)*dgm10(3,3)+9*gm(2,3)*(dgm01(3,3)&
920 &   *dgm10(2,3)+dgm01(2,3)*dgm10(3,3))-1.5d0*gm(3,3)**2*d2gm(2,2)&
921 &   +4.5d0*gm(2,3)**2*d2gm(3,3)+gm(3,3)*(-3*dgm01(3,3)*dgm10(2,2)&
922 &   +9*dgm01(2,3)*dgm10(2,3)-3*dgm01(2,2)*dgm10(3,3)+9*gm(2,3)*d2gm(2,3)&
923 &   -3*gm(2,2)*d2gm(3,3))
924    cm(1,10)=-1.5d0*gm(1,1)*(dgm01(3,3)*dgm10(1,3)+dgm01(1,3)*dgm10(3,3))&
925 &   +7.5d0*gm(1,3)**2*d2gm(1,3)-1.5d0*gm(3,3)*(dgm01(1,3)*dgm10(1,1)&
926 &   +dgm01(1,1)*dgm10(1,3)+gm(1,3)*d2gm(1,1)+gm(1,1)*d2gm(1,3))+gm(1,3)&
927 &   *(-1.5d0*dgm01(3,3)*dgm10(1,1)+15*dgm01(1,3)*dgm10(1,3)-1.5d0*dgm01(1,1)&
928 &   *dgm10(3,3)-1.5d0*gm(1,1)*d2gm(3,3))
929    cm(2,10)=-1.5d0*gm(2,2)*dgm01(3,3)*dgm10(1,3)-1.5d0*gm(1,3)&
930 &   *dgm01(3,3)*dgm10(2,2)+15*gm(1,3)*dgm01(2,3)*dgm10(2,3)-3*gm(1,2)&
931 &   *dgm01(3,3)*dgm10(2,3)-1.5d0*gm(2,2)*dgm01(1,3)*dgm10(3,3)-1.5d0*gm(1,3)&
932 &   *dgm01(2,2)*dgm10(3,3)-3*gm(1,2)*dgm01(2,3)*dgm10(3,3)+7.5d0*gm(2,3)&
933 &   **2*d2gm(1,3)+gm(3,3)*(-3*dgm01(2,3)*dgm10(1,2)-1.5d0*dgm01(2,2)&
934 &   *dgm10(1,3)-1.5d0*dgm01(1,3)*dgm10(2,2)-3*dgm01(1,2)*dgm10(2,3)&
935 &   -3*gm(2,3)*d2gm(1,2)-1.5d0*gm(2,2)*d2gm(1,3)-1.5d0*gm(1,3)*d2gm(2,2)&
936 &   -3*gm(1,2)*d2gm(2,3))-1.5d0*gm(1,3)*gm(2,2)*d2gm(3,3)+gm(2,3)&
937 &   *(-3*dgm01(3,3)*dgm10(1,2)+15*dgm01(2,3)*dgm10(1,3)+15*dgm01(1,3)&
938 &   *dgm10(2,3)-3*dgm01(1,2)*dgm10(3,3)+15*gm(1,3)*d2gm(2,3)-3*gm(1,2)&
939 &   *d2gm(3,3))
940    cm(3,10)=6*gm(1,3)*dgm01(3,3)*dgm10(3,3)+3*gm(3,3)**2*d2gm(1,3)&
941 &   +6*gm(3,3)*(dgm01(3,3)*dgm10(1,3)+dgm01(1,3)*dgm10(3,3)+gm(1,3)&
942 &   *d2gm(3,3))
943    cm(4,10)=9*gm(1,3)*dgm01(3,3)*dgm10(2,3)+9*gm(1,3)*dgm01(2,3)&
944 &   *dgm10(3,3)-6*gm(1,2)*dgm01(3,3)*dgm10(3,3)-3*gm(3,3)**2*d2gm(1,2)&
945 &   +gm(3,3)*(-6*dgm01(3,3)*dgm10(1,2)+9*dgm01(2,3)*dgm10(1,3)+9*dgm01(1,3)&
946 &   *dgm10(2,3)-6*dgm01(1,2)*dgm10(3,3)+9*gm(2,3)*d2gm(1,3)+9*gm(1,3)&
947 &   *d2gm(2,3)-6*gm(1,2)*d2gm(3,3))+9*gm(2,3)*(dgm01(3,3)*dgm10(1,3)&
948 &   +dgm01(1,3)*dgm10(3,3)+gm(1,3)*d2gm(3,3))
949    cm(5,10)=-3*gm(1,1)*dgm01(3,3)*dgm10(3,3)+9*gm(1,3)*(dgm01(3,3)&
950 &   *dgm10(1,3)+dgm01(1,3)*dgm10(3,3))-1.5d0*gm(3,3)**2*d2gm(1,1)&
951 &   +4.5d0*gm(1,3)**2*d2gm(3,3)+gm(3,3)*(-3*dgm01(3,3)*dgm10(1,1)&
952 &   +9*dgm01(1,3)*dgm10(1,3)-3*dgm01(1,1)*dgm10(3,3)+9*gm(1,3)*d2gm(1,3)&
953 &   -3*gm(1,1)*d2gm(3,3))
954    cm(6,10)=-3*gm(1,3)*dgm01(3,3)*dgm10(1,2)+15*gm(1,3)*dgm01(2,3)&
955 &   *dgm10(1,3)-3*gm(1,2)*dgm01(3,3)*dgm10(1,3)+15*gm(1,3)*dgm01(1,3)&
956 &   *dgm10(2,3)-1.5d0*gm(1,1)*dgm01(3,3)*dgm10(2,3)-3*gm(1,3)*dgm01(1,2)&
957 &   *dgm10(3,3)-3*gm(1,2)*dgm01(1,3)*dgm10(3,3)-1.5d0*gm(1,1)*dgm01(2,3)&
958 &   *dgm10(3,3)+7.5d0*gm(1,3)**2*d2gm(2,3)+gm(3,3)*(-1.5d0*dgm01(2,3)&
959 &   *dgm10(1,1)-3*dgm01(1,3)*dgm10(1,2)-3*dgm01(1,2)*dgm10(1,3)-1.5d0*dgm01(1,1)&
960 &   *dgm10(2,3)-1.5d0*gm(2,3)*d2gm(1,1)-3*gm(1,3)*d2gm(1,2)-3*gm(1,2)&
961 &   *d2gm(1,3)-1.5d0*gm(1,1)*d2gm(2,3))-3*gm(1,2)*gm(1,3)*d2gm(3,3)&
962 &   +gm(2,3)*(-1.5d0*dgm01(3,3)*dgm10(1,1)+15*dgm01(1,3)*dgm10(1,3)&
963 &   -1.5d0*dgm01(1,1)*dgm10(3,3)+15*gm(1,3)*d2gm(1,3)-1.5d0*gm(1,1)&
964 &   *d2gm(3,3))
965    cm(7,10)=-1.5d0*gm(2,2)*(dgm01(3,3)*dgm10(2,3)+dgm01(2,3)*dgm10(3,3))&
966 &   +7.5d0*gm(2,3)**2*d2gm(2,3)-1.5d0*gm(3,3)*(dgm01(2,3)*dgm10(2,2)&
967 &   +dgm01(2,2)*dgm10(2,3)+gm(2,3)*d2gm(2,2)+gm(2,2)*d2gm(2,3))+gm(2,3)&
968 &   *(-1.5d0*dgm01(3,3)*dgm10(2,2)+15*dgm01(2,3)*dgm10(2,3)-1.5d0*dgm01(2,2)&
969 &   *dgm10(3,3)-1.5d0*gm(2,2)*d2gm(3,3))
970    cm(8,10)=6*gm(2,3)*dgm01(3,3)*dgm10(3,3)+3*gm(3,3)**2*d2gm(2,3)&
971 &   +6*gm(3,3)*(dgm01(3,3)*dgm10(2,3)+dgm01(2,3)*dgm10(3,3)+gm(2,3)&
972 &   *d2gm(3,3))
973    cm(9,10)=-3*gm(2,2)*dgm01(3,3)*dgm10(3,3)+9*gm(2,3)*(dgm01(3,3)&
974 &   *dgm10(2,3)+dgm01(2,3)*dgm10(3,3))-1.5d0*gm(3,3)**2*d2gm(2,2)&
975 &   +4.5d0*gm(2,3)**2*d2gm(3,3)+gm(3,3)*(-3*dgm01(3,3)*dgm10(2,2)&
976 &   +9*dgm01(2,3)*dgm10(2,3)-3*dgm01(2,2)*dgm10(3,3)+9*gm(2,3)*d2gm(2,3)&
977 &   -3*gm(2,2)*d2gm(3,3))
978    cm(10,10)=gm(3,3)*(6*dgm01(3,3)*dgm10(3,3)+3*gm(3,3)*d2gm(3,3))
979 
980 
981  end if
982 !
983 !contraction to output scalar
984 !
985  e2nl=0.d0
986  do jj=1,((rank+1)*(rank+2))/2
987    tmp(:)=0.d0
988    do ii=1,((rank+1)*(rank+2))/2
989      tmp(:)=tmp(:)+aa(:,ii)*cm(ii,jj)
990    end do
991    e2nl=e2nl+tmp(1)*bb(1,jj)+tmp(2)*bb(2,jj)
992  end do
993 
994  ABI_DEALLOCATE(cm)
995 
996 end subroutine contstr26