1 /* Ergo, version 3.8, a program for linear scaling electronic structure
2 * calculations.
3 * Copyright (C) 2019 Elias Rudberg, Emanuel H. Rubensson, Pawel Salek,
4 * and Anastasia Kruchinina.
5 *
6 * This program is free software: you can redistribute it and/or modify
7 * it under the terms of the GNU General Public License as published by
8 * the Free Software Foundation, either version 3 of the License, or
9 * (at your option) any later version.
10 *
11 * This program is distributed in the hope that it will be useful,
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 * GNU General Public License for more details.
15 *
16 * You should have received a copy of the GNU General Public License
17 * along with this program. If not, see <http://www.gnu.org/licenses/>.
18 *
19 * Primary academic reference:
20 * Ergo: An open-source program for linear-scaling electronic structure
21 * calculations,
22 * Elias Rudberg, Emanuel H. Rubensson, Pawel Salek, and Anastasia
23 * Kruchinina,
24 * SoftwareX 7, 107 (2018),
25 * <http://dx.doi.org/10.1016/j.softx.2018.03.005>
26 *
27 * For further information about Ergo, see <http://www.ergoscf.org>.
28 */
29
30 /*-*-mode: C; c-indentation-style: "bsd"; c-basic-offset: 4; -*-*/
31 /** @file fun-p86c.c P86c implementation.
32
33 Automatically generated code implementing p86c functional and its
34 derivatives. Generated by func-codegen.pl being a part of a
35 "Automatic code generation framework for analytical functional
36 derivative evaluation", Pawel Salek, 2004
37
38 This functional is connected by making following changes:
39 1. add "extern Functional p86cFunctional;" to 'functionals.h'
40 2. add "&p86cFunctional," to 'functionals.c'
41 3. add "fun-p86c}.c" to 'Makefile.am', 'Makefile.in' or 'Makefile'.
42
43 This functional has been generated from following input:
44 ------ cut here -------
45 c0: 1.667e-3;
46 cn1: 2.568e-3;
47 cn2: 2.3266e-2;
48 cn3: 7.389e-6;
49 cd1: 1.0;
50 cd2: 8.723;
51 cd3: 0.472;
52 cd4: 1.0e4*cn3;
53 cinf: c0+cn1/cd1;
54 phi1: (9.0*%PI)^(1/6);
55 phi2: 0.11;
56
57 rho: rhoa + rhob;
58 grad: sqrt(grada*grada + gradb*gradb + 2*gradab);
59 zeta: (rhoa-rhob)/(rhoa+rhob);
60 dzet: 2^(1/3)*sqrt(((1+zeta)/2)^(5/3) + ((1-zeta)/2)^(5/3));
61
62 rs: (3/(4*%PI*rho))^(1/3);
63 crho: c0 + (cn1+cn2*rs+cn3*rs^2)/(1+cd2*rs+cd3*(rs^2)+cd4*(rs^3));
64 phi: phi1*phi2*(cinf/crho)*grad*(rho^(-7.0/6.0));
65
66 K(rhoa,rhob,grada,gradb,gradab):=exp(-phi)*crho*(grad^2)*(rho^(-4.0/3.0))/dzet;
67
68 ------ cut here -------
69 */
70
71 #include <math.h>
72 #include <stddef.h>
73
74 #define __CVERSION__
75
76 #include "functionals.h"
77
78 /* INTERFACE PART */
p86c_isgga(void)79 static int p86c_isgga(void) { return 1; } /* FIXME: detect! */
80 static int p86c_read(const char *conf_line);
81 static real p86c_energy(const FunDensProp *dp);
82 static void p86c_first(FunFirstFuncDrv *ds, real factor,
83 const FunDensProp *dp);
84 static void p86c_second(FunSecondFuncDrv *ds, real factor,
85 const FunDensProp *dp);
86 static void p86c_third(FunThirdFuncDrv *ds, real factor,
87 const FunDensProp *dp);
88
89 Functional P86cFunctional = {
90 "P86c", /* name */
91 p86c_isgga, /* gga-corrected */
92 p86c_read,
93 NULL,
94 p86c_energy,
95 p86c_first,
96 p86c_second,
97 p86c_third
98 };
99
100 /* IMPLEMENTATION PART */
101 static int
p86c_read(const char * conf_line)102 p86c_read(const char *conf_line)
103 {
104 fun_set_hf_weight(0);
105 return 1;
106 }
107
108 static real
p86c_energy(const FunDensProp * dp)109 p86c_energy(const FunDensProp *dp)
110 {
111 real res;
112 real rhoa = dp->rhoa, rhob = dp->rhob;
113 real grada = dp->grada, gradb = dp->gradb, gradab = dp->gradab;
114
115 real t1, t2, t3, t4, t5, t6, t7, t8, t9, t10;
116 real t11, t12, t13, t14;
117
118 t1 = POW(gradb,2.0)+2.0*gradab+POW(grada,2.0);
119 t2 = rhob+rhoa;
120 t3 = 1/POW(2.0,0.66666666666667);
121 t4 = rhoa-1.0*rhob;
122 t5 = 1/t2;
123 t6 = POW(3.0,0.66666666666667);
124 t7 = 1/POW(4.0,0.66666666666667);
125 t8 = 1/POW(3.141592653589793,0.66666666666667);
126 t9 = 1/POW(t2,0.66666666666667);
127 t10 = POW(3.0,0.33333333333333);
128 t11 = 1/POW(4.0,0.33333333333333);
129 t12 = 1/POW(3.141592653589793,0.33333333333333);
130 t13 = 1/POW(t2,0.33333333333333);
131 t14 = (0.023266*t10*t11*t12*t13+7.389E-6*t6*t7*t8*t9+
132 0.002568)/(0.472*t6*t7*t8*t9+0.01763993811759*t5+8.723000000000001*
133 t10*t11*t12*t13+1.0)+0.001667;
134
135 /* code */
136 res = t1*t14/(POW(2.0,0.33333333333333)*POW(2.718281828459045,
137 6.718719623277062E-4*POW(3.141592653589793,0.16666666666667)*
138 SQRT(t1)/(t14*POW(t2,1.166666666666667)))*POW(t2,1.333333333333333)*
139 SQRT(0.5*t3*POW(t4*t5+1.0,1.666666666666667)+0.5*t3*POW(1.0-
140 1.0*t4*t5,1.666666666666667)));
141
142 return res;
143 }
144
145 static void
p86c_first(FunFirstFuncDrv * ds,real factor,const FunDensProp * dp)146 p86c_first(FunFirstFuncDrv *ds, real factor, const FunDensProp *dp)
147 {
148 real dfdra, dfdrb, dfdga, dfdgb, dfdgab;
149 real rhoa = dp->rhoa, rhob = dp->rhob;
150 real grada = dp->grada, gradb = dp->gradb, gradab = dp->gradab;
151
152 real t1, t2, t3, t4, t5, t6, t7, t8, t9, t10;
153 real t11, t12, t13, t14, t15, t16, t17, t18;
154 real t19, t20, t21, t22, t23, t24, t25, t26;
155 real t27, t28, t29, t30, t31, t32, t33, t34;
156 real t35, t36, t37, t38, t39, t40, t41, t42;
157 real t43;
158
159 t1 = 1/POW(2.0,0.33333333333333);
160 t2 = POW(gradb,2.0)+2.0*gradab+POW(grada,2.0);
161 t3 = rhob+rhoa;
162 t4 = 1/POW(t3,1.333333333333333);
163 t5 = 1/POW(2.0,0.66666666666667);
164 t6 = rhoa-1.0*rhob;
165 t7 = 1/t3;
166 t8 = 1.0-1.0*t6*t7;
167 t9 = t6*t7+1.0;
168 t10 = SQRT(0.5*t5*POW(t9,1.666666666666667)+0.5*t5*POW(t8,
169 1.666666666666667));
170 t11 = 1/t10;
171 t12 = 0.31830988618379;
172 t13 = 1/POW(t3,2.0);
173 t14 = POW(3.0,0.66666666666667);
174 t15 = 1/POW(4.0,0.66666666666667);
175 t16 = 1/POW(3.141592653589793,0.66666666666667);
176 t17 = 1/POW(t3,1.666666666666667);
177 t18 = POW(3.0,0.33333333333333);
178 t19 = 1/POW(4.0,0.33333333333333);
179 t20 = 1/POW(3.141592653589793,0.33333333333333);
180 t21 = 1/POW(t3,1.333333333333333);
181 t22 = 1/POW(t3,0.66666666666667);
182 t23 = 1/POW(t3,0.33333333333333);
183 t24 = 0.023266*t18*t19*t20*t23+7.389E-6*t14*t15*t16*t22+
184 0.002568;
185 t25 = 0.0554175*t12*t7+8.723000000000001*t18*t19*t20*
186 t23+0.472*t14*t15*t16*t22+1.0;
187 t26 = 1/t25;
188 t27 = (-0.00775533333333*t18*t19*t20*t21-4.926E-6*t14*
189 t15*t16*t17)*t26-1.0*(-2.907666666666667*t18*t19*t20*t21-0.31466666666667*
190 t14*t15*t16*t17-0.0554175*t12*t13)*t24/POW(t25,2.0);
191 t28 = POW(3.141592653589793,0.16666666666667);
192 t29 = SQRT(t2);
193 t30 = 1/POW(t3,1.166666666666667);
194 t31 = t24*t26+0.001667;
195 t32 = 1/t31;
196 t33 = 1/POW(2.718281828459045,6.718719623277062E-4*t28*
197 t29*t30*t32);
198 t34 = t1*t2*t4*t11*t27*t33;
199 t35 = t6*t13;
200 t36 = -1.0*t7;
201 t37 = POW(t8,0.66666666666667);
202 t38 = -1.0*t13*t6;
203 t39 = POW(t9,0.66666666666667);
204 t40 = 1/POW(t10,3.0);
205 t41 = -1.333333333333333*t1*t11*t2*t31*t33/POW(t3,2.333333333333333);
206 t42 = t1*
207 t11*t2*t31*(7.838506227156574E-4*t28*t29*t32/POW(t3,2.166666666666667)+
208 6.718719623277062E-4*t27*t28*t29*t30/POW(t31,2.0))*t33*t4;
209 t43 = 1/
210 POW(t3,2.5);
211
212 /* code */
213 dfdra = -0.5*t1*t2*t31*t33*t4*t40*(0.83333333333333*t39*
214 t5*(t7+t38)+0.83333333333333*(t36+t35)*t37*t5)+t42+t41+t34;
215 dfdrb = -
216 0.5*t1*t2*t31*t33*t4*t40*(0.83333333333333*t37*t5*(t7+t35)+
217 0.83333333333333*(t36+t38)*t39*t5)+t42+t41+t34;
218 dfdga = 2.0*t1*t11*t31*t33*t4*grada-6.718719623277062E-4*
219 t1*t28*grada*t29*t43*t11*t33;
220 dfdgb = 2.0*t1*t11*t31*t33*t4*gradb-6.718719623277062E-4*
221 t1*t28*gradb*t29*t43*t11*t33;
222 dfdgab = 2.0*t1*t11*t31*t33*t4-6.718719623277062E-4*t1*
223 t28*t29*t43*t11*t33;
224
225
226 /* Final assignment */
227 ds->df1000 += factor*dfdra;
228 ds->df0100 += factor*dfdrb;
229 ds->df0010 += factor*dfdga;
230 ds->df0001 += factor*dfdgb;
231 ds->df00001+= factor*dfdgab;
232 }
233
234 static void
p86c_second(FunSecondFuncDrv * ds,real factor,const FunDensProp * dp)235 p86c_second(FunSecondFuncDrv *ds, real factor, const FunDensProp *dp)
236 {
237 real dfdra, dfdrb, dfdga, dfdgb, dfdgab;
238 real d2fdrara, d2fdrarb, d2fdraga, d2fdragb, d2fdrbrb, d2fdraab,
239 d2fdgaga, d2fdgbgb, d2fdrbga, d2fdrbgb,
240 d2fdrbgab, d2fdgagb, d2fdgagab, d2fdgbgab, d2fdgabgab;
241 real rhoa = dp->rhoa, rhob = dp->rhob;
242 real grada = dp->grada, gradb = dp->gradb, gradab = dp->gradab;
243
244 real t1, t2, t3, t4, t5, t6, t7, t8, t9, t10;
245 real t11, t12, t13, t14, t15, t16, t17, t18;
246 real t19, t20, t21, t22, t23, t24, t25, t26;
247 real t27, t28, t29, t30, t31, t32, t33, t34;
248 real t35, t36, t37, t38, t39, t40, t41, t42;
249 real t43, t44, t45, t46, t47, t48, t49, t50;
250 real t51, t52, t53, t54, t55, t56, t57, t58;
251 real t59, t60, t61, t62, t63, t64, t65, t66;
252 real t67, t68, t69, t70, t71, t72, t73, t74;
253 real t75, t76, t77, t78, t79, t80, t81, t82;
254 real t83, t84, t85, t86, t87, t88, t89, t90;
255 real t91, t92, t93, t94, t95, t96, t97, t98;
256 real t99, t100, t101, t102, t103, t104, t105;
257
258 t1 = 1/POW(2.0,0.33333333333333);
259 t2 = POW(grada,2.0);
260 t3 = POW(gradb,2.0);
261 t4 = 2.0*gradab+t3+t2;
262 t5 = rhob+rhoa;
263 t6 = 1/POW(t5,1.333333333333333);
264 t7 = 1/POW(2.0,0.66666666666667);
265 t8 = rhoa-1.0*rhob;
266 t9 = 1/t5;
267 t10 = 1.0-1.0*t8*t9;
268 t11 = t8*t9+1.0;
269 t12 = SQRT(0.5*POW(t11,1.666666666666667)*t7+0.5*POW(t10,
270 1.666666666666667)*t7);
271 t13 = 1/t12;
272 t14 = 0.31830988618379;
273 t15 = 1/POW(t5,2.0);
274 t16 = POW(3.0,0.66666666666667);
275 t17 = 1/POW(4.0,0.66666666666667);
276 t18 = 1/POW(3.141592653589793,0.66666666666667);
277 t19 = 1/POW(t5,1.666666666666667);
278 t20 = POW(3.0,0.33333333333333);
279 t21 = 1/POW(4.0,0.33333333333333);
280 t22 = POW(3.141592653589793,0.33333333333333);
281 t23 = 1/t22;
282 t24 = 1/POW(t5,1.333333333333333);
283 t25 = -2.907666666666667*t20*t21*t23*t24-0.31466666666667*
284 t16*t17*t18*t19-0.0554175*t14*t15;
285 t26 = 1/POW(t5,0.66666666666667);
286 t27 = 1/POW(t5,0.33333333333333);
287 t28 = 0.023266*t20*t21*t23*t27+7.389E-6*t16*t17*t18*t26+
288 0.002568;
289 t29 = 0.0554175*t14*t9+8.723000000000001*t20*t21*t23*
290 t27+0.472*t16*t17*t18*t26+1.0;
291 t30 = 1/POW(t29,2.0);
292 t31 = -0.00775533333333*t20*t21*t23*t24-4.926E-6*t16*
293 t17*t18*t19;
294 t32 = 1/t29;
295 t33 = t31*t32-1.0*t25*t28*t30;
296 t34 = POW(3.141592653589793,0.16666666666667);
297 t35 = SQRT(t4);
298 t36 = 1/POW(t5,1.166666666666667);
299 t37 = t28*t32+0.001667;
300 t38 = 1/t37;
301 t39 = 1/POW(2.718281828459045,6.718719623277062E-4*t34*
302 t35*t36*t38);
303 t40 = t1*t4*t6*t13*t33*t39;
304 t41 = t8*t15;
305 t42 = -1.0*t9;
306 t43 = t42+t41;
307 t44 = POW(t10,0.66666666666667);
308 t45 = -1.0*t15*t8;
309 t46 = t9+t45;
310 t47 = POW(t11,0.66666666666667);
311 t48 = 0.83333333333333*t46*t47*t7+0.83333333333333*t43*
312 t44*t7;
313 t49 = 1/POW(t12,3.0);
314 t50 = 1/POW(t5,2.333333333333333);
315 t51 = -1.333333333333333*t1*t4*t50*t13*t37*t39;
316 t52 = 1/POW(t37,2.0);
317 t53 = 1/POW(t5,2.166666666666667);
318 t54 = 7.838506227156574E-4*t34*t35*t53*t38+6.718719623277062E-4*
319 t34*t35*t36*t33*t52;
320 t55 = t1*t4*t6*t13*t37*t54*t39;
321 t56 = t9+t41;
322 t57 = t42+t45;
323 t58 = 0.83333333333333*t47*t57*t7+0.83333333333333*t44*
324 t56*t7;
325 t59 = 1/POW(t5,2.5);
326 t60 = -6.718719623277062E-4*t1*t34*t35*t59*t13*t39;
327 t61 = 2.0*
328 t1*t13*t37*t39*t6;
329 t62 = 1/POW(t5,3.0);
330 t63 = 1/POW(t5,2.666666666666667);
331 t64 = 1/POW(t5,2.333333333333334);
332 t65 = -1.0*t28*t30*(3.876888888888889*t20*t21*t23*t64+
333 0.52444444444444*t16*t17*t18*t63+0.110835*t14*t62)+(0.01034044444444*
334 t20*t21*t23*t64+8.210000000000001E-6*t16*t17*t18*t63)*t32-
335 2.0*t25*t30*t31+2.0*POW(t25,2.0)*t28/POW(t29,3.0);
336 t66 = t1*t4*t6*t13*t65*t39;
337 t67 = -2.666666666666667*t1*t4*t50*t13*t33*t39;
338 t68 = 1/POW(t12,5.0);
339 t69 = 1/POW(t10,0.33333333333333);
340 t70 = -2.0*t62*t8;
341 t71 = 2.0*t15;
342 t72 = 1/POW(t11,0.33333333333333);
343 t73 = 2.0*t62*t8;
344 t74 = -2.0*t15;
345 t75 = 3.111111111111111*t1*t13*t37*t39*t4/POW(t5,3.333333333333333);
346 t76 = t1*
347 t13*t37*t39*t4*(6.718719623277062E-4*t34*t35*t36*t65*t52-0.00156770124543*
348 t34*t35*t53*t33*t52-0.00169834301588*t34*t35*t38/POW(t5,3.166666666666667)-
349 0.00134374392466*POW(t33,2.0)*t34*t35*t36/POW(t37,3.0))*t6;
350 t77 = 2.0*
351 t1*t13*t33*t39*t4*t54*t6;
352 t78 = -2.666666666666667*t1*t4*t50*t13*t37*t54*t39;
353 t79 = t1*
354 t13*t37*t39*t4*POW(t54,2.0)*t6;
355 t80 = 1/POW(t5,3.5);
356 t81 = 8.958292831036083E-4*t1*t34*grada*t35*t80*t13*t39;
357 t82 = 2.0*
358 t1*t13*t33*t39*t6*grada;
359 t83 = -6.718719623277062E-4*t1*t34*grada*t35*t59*t13*
360 t33*t38*t39;
361 t84 = -2.666666666666667*t1*grada*t50*t13*t37*t39;
362 t85 = 1/t35;
363 t86 = t1*t4*t6*t13*t37*(7.838506227156574E-4*t34*grada*
364 t85*t53*t38+6.718719623277062E-4*t34*grada*t85*t36*t33*t52)*
365 t39;
366 t87 = -6.718719623277062E-4*t1*t34*grada*t35*t59*t13*
367 t54*t39;
368 t88 = 2.0*t1*t13*t37*t39*t54*t6*grada;
369 t89 = 8.958292831036083E-4*t1*t34*gradb*t35*t80*t13*t39;
370 t90 = 2.0*
371 t1*t13*t33*t39*t6*gradb;
372 t91 = -6.718719623277062E-4*t1*t34*gradb*t35*t59*t13*
373 t33*t38*t39;
374 t92 = -2.666666666666667*t1*gradb*t50*t13*t37*t39;
375 t93 = t1*t4*t6*t13*t37*(7.838506227156574E-4*t34*gradb*
376 t85*t53*t38+6.718719623277062E-4*t34*gradb*t85*t36*t33*t52)*
377 t39;
378 t94 = -6.718719623277062E-4*t1*t34*gradb*t35*t59*t13*
379 t54*t39;
380 t95 = 2.0*t1*t13*t37*t39*t54*t6*gradb;
381 t96 = 8.958292831036083E-4*t1*t34*t35*t80*t13*t39;
382 t97 = 2.0*t1*t13*t33*t39*t6;
383 t98 = -6.718719623277062E-4*t1*t34*t35*t59*t13*t33*t38*
384 t39;
385 t99 = -2.666666666666667*t1*t50*t13*t37*t39;
386 t100 = t1*t4*t6*t13*t37*(7.838506227156574E-4*t34*t85*
387 t53*t38+6.718719623277062E-4*t34*t85*t36*t33*t52)*t39;
388 t101 = -6.718719623277062E-4*t1*t34*t35*t59*t13*t54*t39;
389 t102 = 2.0*
390 t1*t13*t37*t39*t54*t6;
391 t103 = t99+t98+t97+t96-1.0*t1*t37*t39*t48*t49*t6+3.359359811638531E-4*
392 t1*t34*t35*t59*t48*t49*t39+t102+t101+t100;
393 t104 = t99+t98+t97+t96-1.0*t1*t37*t39*t49*t58*t6+3.359359811638531E-4*
394 t1*t34*t35*t59*t58*t49*t39+t102+t101+t100;
395 t105 = 1/POW(t5,3.666666666666667);
396
397 /* code */
398 dfdra = -0.5*t1*t37*t39*t4*t48*t49*t6+t55+t51+t40;
399 dfdrb = -0.5*t1*t37*t39*t4*t49*t58*t6+t55+t51+t40;
400 dfdga = 2.0*t1*t13*t37*t39*t6*grada-6.718719623277062E-4*
401 t1*t34*grada*t35*t59*t13*t39;
402 dfdgb = 2.0*t1*t13*t37*t39*t6*gradb-6.718719623277062E-4*
403 t1*t34*gradb*t35*t59*t13*t39;
404 dfdgab = t61+t60;
405 d2fdrara = t79+t78+t77+t76+t75-0.5*t1*t37*t39*t4*t49*
406 t6*(0.83333333333333*t47*t7*(t74+t73)+0.55555555555556*POW(t46,
407 2.0)*t7*t72+0.83333333333333*t44*t7*(t71+t70)+0.55555555555556*
408 POW(t43,2.0)*t69*t7)+0.75*t1*t37*t39*t4*POW(t48,2.0)*t6*t68+
409 t67+t66-1.0*t1*t37*t39*t4*t48*t49*t54*t6-1.0*t1*t33*t39*t4*
410 t48*t49*t6+1.333333333333333*t1*t4*t50*t48*t49*t37*t39;
411 d2fdrarb = -0.5*t1*t37*t39*t4*t49*t6*(1.666666666666667*
412 t47*t62*t7*t8-1.666666666666667*t44*t62*t7*t8+0.55555555555556*
413 t46*t57*t7*t72+0.55555555555556*t43*t56*t69*t7)+t79+t78+t77+
414 t76+t75+0.75*t1*t37*t39*t4*t48*t58*t6*t68+t67+t66-0.5*t1*t37*
415 t39*t4*t49*t54*t58*t6-0.5*t1*t33*t39*t4*t49*t58*t6-0.5*t1*
416 t37*t39*t4*t48*t49*t54*t6-0.5*t1*t33*t39*t4*t48*t49*t6+0.66666666666667*
417 t1*t4*t50*t58*t49*t37*t39+0.66666666666667*t1*t4*t50*t48*t49*
418 t37*t39;
419 d2fdraga = -1.0*t1*t37*t39*t48*t49*t6*grada+t88+t87+t86+
420 t84+t83+t82+t81+3.359359811638531E-4*t1*t34*grada*t35*t59*
421 t48*t49*t39;
422 d2fdragb = -1.0*t1*t37*t39*t48*t49*t6*gradb+t95+t94+t93+
423 t92+t91+t90+t89+3.359359811638531E-4*t1*t34*gradb*t35*t59*
424 t48*t49*t39;
425 d2fdrbrb = t79+t78+t77+t76+t75-0.5*t1*t37*t39*t4*t49*
426 t6*(0.83333333333333*t44*t7*(t74+t70)+0.55555555555556*POW(t57,
427 2.0)*t7*t72+0.83333333333333*t47*t7*(t71+t73)+0.55555555555556*
428 POW(t56,2.0)*t69*t7)+0.75*t1*t37*t39*t4*POW(t58,2.0)*t6*t68+
429 t67+t66-1.0*t1*t37*t39*t4*t49*t54*t58*t6-1.0*t1*t33*t39*t4*
430 t49*t58*t6+1.333333333333333*t1*t4*t50*t58*t49*t37*t39;
431 d2fdraab = t103;
432 d2fdgaga = t61+4.514119337620827E-7*t1*t22*t2*t105*t13*
433 t38*t39+t60-0.00201561588698*t1*t34*t2*t85*t59*t13*t39;
434 d2fdgbgb = t61+4.514119337620827E-7*t1*t22*t3*t105*t13*
435 t38*t39+t60-0.00201561588698*t1*t34*t3*t85*t59*t13*t39;
436 d2fdrbga = -1.0*t1*t37*t39*t49*t58*t6*grada+t88+t87+t86+
437 t84+t83+t82+t81+3.359359811638531E-4*t1*t34*grada*t35*t59*
438 t58*t49*t39;
439 d2fdrbgb = -1.0*t1*t37*t39*t49*t58*t6*gradb+t95+t94+t93+
440 t92+t91+t90+t89+3.359359811638531E-4*t1*t34*gradb*t35*t59*
441 t58*t49*t39;
442 d2fdgagb = 4.514119337620827E-7*t1*t22*grada*gradb*t105*
443 t13*t38*t39-0.00201561588698*t1*t34*grada*gradb*t85*t59*t13*
444 t39;
445 d2fdrbgab = t104;
446 d2fdgagab = 4.514119337620827E-7*t1*t22*grada*t105*t13*
447 t38*t39-0.00201561588698*t1*t34*grada*t85*t59*t13*t39;
448 d2fdgbgab = 4.514119337620827E-7*t1*t22*gradb*t105*t13*
449 t38*t39-0.00201561588698*t1*t34*gradb*t85*t59*t13*t39;
450 d2fdgabgab = 4.514119337620827E-7*t1*t22*t105*t13*t38*
451 t39-0.00201561588698*t1*t34*t85*t59*t13*t39;
452
453
454 ds->df1000 += factor*dfdra;
455 ds->df0100 += factor*dfdrb;
456 ds->df0010 += factor*dfdga;
457 ds->df0001 += factor*dfdgb;
458 ds->df00001+= factor*dfdgab;
459
460 ds->df2000 += factor*d2fdrara;
461 ds->df1100 += factor*d2fdrarb;
462 ds->df1010 += factor*d2fdraga;
463 ds->df1001 += factor*d2fdragb;
464 ds->df10001+= factor*d2fdraab;
465 ds->df0200 += factor*d2fdrbrb;
466 ds->df0110 += factor*d2fdrbga;
467 ds->df0101 += factor*d2fdrbgb;
468 ds->df01001+= factor*d2fdrbgab;
469 ds->df0020 += factor*d2fdgaga;
470 ds->df0011 += factor*d2fdgagb;
471 ds->df00101+= factor*d2fdgagab;
472 ds->df0002 += factor*d2fdgbgb;
473 ds->df00011+= factor*d2fdgbgab;
474 ds->df00002+= factor*d2fdgabgab;
475 }
476
477 static void
p86c_third(FunThirdFuncDrv * ds,real factor,const FunDensProp * dp)478 p86c_third(FunThirdFuncDrv *ds, real factor, const FunDensProp *dp)
479 {
480 real dfdra, dfdrb, dfdga, dfdgb, dfdgab;
481 real d2fdrara, d2fdrarb, d2fdraga, d2fdragb, d2fdrbrb, d2fdraab,
482 d2fdrbab, d2fdgaga, d2fdgbgb, d2fdrbga, d2fdrbgb,
483 d2fdgagb;
484 real d3fdrararb, d3fdraraga, d3fdraragb, d3fdrbrbab, d3fdraraab,
485 d3fdrarbrb, d3fdrarbga, d3fdrarbgb, d3fdrarbab, d3fdragaga,
486 d3fdragbgb, d3fdrarara, d3fdrbrbrb, d3fdrbrbga, d3fdrbrbgb,
487 d3fdrbgaga, d3fdrbgbgb, d3fdgagaga,
488 d3fdragagb, d3fdrbgagb, d3fdgagagb, d3fdgagbgb, d3fdgbgbgb;
489 real rhoa = dp->rhoa, rhob = dp->rhob;
490 real grada = dp->grada, gradb = dp->gradb, gradab = dp->gradab;
491
492 real t1, t2, t3, t4, t5, t6, t7, t8, t9, t10;
493 real t11, t12, t13, t14, t15, t16, t17, t18;
494 real t19, t20, t21, t22, t23, t24, t25, t26;
495 real t27, t28, t29, t30, t31, t32, t33, t34;
496 real t35, t36, t37, t38, t39, t40, t41, t42;
497 real t43, t44, t45, t46, t47, t48, t49, t50;
498 real t51, t52, t53, t54, t55, t56, t57, t58;
499 real t59, t60, t61, t62, t63, t64, t65, t66;
500 real t67, t68, t69, t70, t71, t72, t73, t74;
501 real t75, t76, t77, t78, t79, t80, t81, t82;
502 real t83, t84, t85, t86, t87, t88, t89, t90;
503 real t91, t92, t93, t94, t95, t96, t97, t98;
504 real t99, t100, t101, t102, t103, t104, t105;
505 real t106, t107, t108, t109, t110, t111, t112;
506 real t113, t114, t115, t116, t117, t118, t119;
507 real t120, t121, t122, t123, t124, t125, t126;
508 real t127, t128, t129, t130, t131, t132, t133;
509 real t134, t135, t136, t137, t138, t139, t140;
510 real t141, t142, t143, t144, t145, t146, t147;
511 real t148, t149, t150, t151, t152, t153, t154;
512 real t155, t156, t157, t158, t159, t160, t161;
513 real t162, t163, t164, t165, t166, t167, t168;
514 real t169, t170, t171, t172, t173, t174, t175;
515 real t176, t177, t178, t179, t180, t181, t182;
516 real t183, t184, t185, t186, t187, t188, t189;
517 real t190, t191, t192, t193, t194, t195, t196;
518 real t197, t198, t199, t200, t201, t202, t203;
519 real t204, t205, t206, t207, t208, t209, t210;
520 real t211, t212, t213, t214, t215, t216, t217;
521 real t218, t219, t220, t221, t222, t223, t224;
522 real t225, t226, t227, t228, t229, t230, t231;
523 real t232, t233, t234, t235, t236, t237, t238;
524 real t239, t240, t241, t242, t243, t244, t245;
525 real t246, t247, t248, t249, t250, t251, t252;
526 real t253, t254, t255, t256, t257, t258, t259;
527 real t260, t261, t262, t263, t264;
528
529 t1 = 1/POW(2.0,0.33333333333333);
530 t2 = POW(grada,2.0);
531 t3 = POW(gradb,2.0);
532 t4 = 2.0*gradab+t3+t2;
533 t5 = rhob+rhoa;
534 t6 = 1/POW(t5,1.333333333333333);
535 t7 = 1/POW(2.0,0.66666666666667);
536 t8 = rhoa-1.0*rhob;
537 t9 = 1/t5;
538 t10 = 1.0-1.0*t8*t9;
539 t11 = t8*t9+1.0;
540 t12 = SQRT(0.5*POW(t11,1.666666666666667)*t7+0.5*POW(t10,
541 1.666666666666667)*t7);
542 t13 = 1/t12;
543 t14 = 0.31830988618379;
544 t15 = 1/POW(t5,2.0);
545 t16 = POW(3.0,0.66666666666667);
546 t17 = 1/POW(4.0,0.66666666666667);
547 t18 = 1/POW(3.141592653589793,0.66666666666667);
548 t19 = 1/POW(t5,1.666666666666667);
549 t20 = POW(3.0,0.33333333333333);
550 t21 = 1/POW(4.0,0.33333333333333);
551 t22 = POW(3.141592653589793,0.33333333333333);
552 t23 = 1/t22;
553 t24 = 1/POW(t5,1.333333333333333);
554 t25 = -2.907666666666667*t20*t21*t23*t24-0.31466666666667*
555 t16*t17*t18*t19-0.0554175*t14*t15;
556 t26 = 1/POW(t5,0.66666666666667);
557 t27 = 1/POW(t5,0.33333333333333);
558 t28 = 0.023266*t20*t21*t23*t27+7.389E-6*t16*t17*t18*t26+
559 0.002568;
560 t29 = 0.0554175*t14*t9+8.723000000000001*t20*t21*t23*
561 t27+0.472*t16*t17*t18*t26+1.0;
562 t30 = 1/POW(t29,2.0);
563 t31 = -0.00775533333333*t20*t21*t23*t24-4.926E-6*t16*
564 t17*t18*t19;
565 t32 = 1/t29;
566 t33 = t31*t32-1.0*t25*t28*t30;
567 t34 = POW(3.141592653589793,0.16666666666667);
568 t35 = SQRT(t4);
569 t36 = 1/POW(t5,1.166666666666667);
570 t37 = t28*t32+0.001667;
571 t38 = 1/t37;
572 t39 = 1/POW(2.718281828459045,6.718719623277062E-4*t34*
573 t35*t36*t38);
574 t40 = t1*t4*t6*t13*t33*t39;
575 t41 = t8*t15;
576 t42 = -1.0*t9;
577 t43 = t42+t41;
578 t44 = POW(t10,0.66666666666667);
579 t45 = -1.0*t15*t8;
580 t46 = t9+t45;
581 t47 = POW(t11,0.66666666666667);
582 t48 = 0.83333333333333*t46*t47*t7+0.83333333333333*t43*
583 t44*t7;
584 t49 = 1/POW(t12,3.0);
585 t50 = 1/POW(t5,2.333333333333333);
586 t51 = -1.333333333333333*t1*t4*t50*t13*t37*t39;
587 t52 = 1/POW(t37,2.0);
588 t53 = 1/POW(t5,2.166666666666667);
589 t54 = 7.838506227156574E-4*t34*t35*t53*t38+6.718719623277062E-4*
590 t34*t35*t36*t33*t52;
591 t55 = t1*t4*t6*t13*t37*t54*t39;
592 t56 = t9+t41;
593 t57 = t42+t45;
594 t58 = 0.83333333333333*t47*t57*t7+0.83333333333333*t44*
595 t56*t7;
596 t59 = 1/POW(t5,2.5);
597 t60 = -6.718719623277062E-4*t1*t34*t35*t59*t13*t39;
598 t61 = 2.0*
599 t1*t13*t37*t39*t6;
600 t62 = POW(t25,2.0);
601 t63 = 1/POW(t29,3.0);
602 t64 = 1/POW(t5,3.0);
603 t65 = 1/POW(t5,2.666666666666667);
604 t66 = 1/POW(t5,2.333333333333334);
605 t67 = 3.876888888888889*t20*t21*t23*t66+0.52444444444444*
606 t16*t17*t18*t65+0.110835*t14*t64;
607 t68 = 0.01034044444444*t20*t21*t23*t66+8.210000000000001E-6*
608 t16*t17*t18*t65;
609 t69 = -1.0*t28*t30*t67+2.0*t28*t62*t63+t68*t32-2.0*t25*
610 t30*t31;
611 t70 = t1*t4*t6*t13*t69*t39;
612 t71 = -2.666666666666667*t1*t4*t50*t13*t33*t39;
613 t72 = POW(t48,2.0);
614 t73 = 1/POW(t12,5.0);
615 t74 = POW(t43,2.0);
616 t75 = 1/POW(t10,0.33333333333333);
617 t76 = -2.0*t64*t8;
618 t77 = 2.0*t15;
619 t78 = t77+t76;
620 t79 = POW(t46,2.0);
621 t80 = 1/POW(t11,0.33333333333333);
622 t81 = 2.0*t64*t8;
623 t82 = -2.0*t15;
624 t83 = t82+t81;
625 t84 = 0.83333333333333*t47*t7*t83+0.55555555555556*t7*
626 t79*t80+0.83333333333333*t44*t7*t78+0.55555555555556*t7*t74*
627 t75;
628 t85 = 1/POW(t5,3.333333333333333);
629 t86 = 3.111111111111111*t1*t4*t85*t13*t37*t39;
630 t87 = POW(t33,2.0);
631 t88 = 1/POW(t37,3.0);
632 t89 = 1/POW(t5,3.166666666666667);
633 t90 = -0.00169834301588*t34*t35*t89*t38-0.00156770124543*
634 t34*t35*t53*t33*t52+6.718719623277062E-4*t34*t35*t36*t69*t52-
635 0.00134374392466*t34*t35*t36*t87*t88;
636 t91 = t1*t4*t6*t13*t37*t90*t39;
637 t92 = 2.0*t1*t13*t33*t39*t4*t54*t6;
638 t93 = -2.666666666666667*t1*t4*t50*t13*t37*t54*t39;
639 t94 = POW(t54,
640 2.0);
641 t95 = t1*t4*t6*t13*t37*t94*t39;
642 t96 = 0.55555555555556*t46*t57*t7*t80+1.666666666666667*
643 t47*t64*t7*t8-1.666666666666667*t44*t64*t7*t8+0.55555555555556*
644 t43*t56*t7*t75;
645 t97 = 1/POW(t5,3.5);
646 t98 = 8.958292831036083E-4*t1*t34*grada*t35*t97*t13*t39;
647 t99 = 2.0*
648 t1*t13*t33*t39*t6*grada;
649 t100 = -6.718719623277062E-4*t1*t34*grada*t35*t59*t13*
650 t33*t38*t39;
651 t101 = -2.666666666666667*t1*grada*t50*t13*t37*t39;
652 t102 = 1/
653 t35;
654 t103 = 7.838506227156574E-4*t34*grada*t102*t53*t38+6.718719623277062E-4*
655 t34*grada*t102*t36*t33*t52;
656 t104 = t1*t4*t6*t13*t37*t103*t39;
657 t105 = -6.718719623277062E-4*t1*t34*grada*t35*t59*t13*
658 t54*t39;
659 t106 = 2.0*t1*t13*t37*t39*t54*t6*grada;
660 t107 = 8.958292831036083E-4*t1*t34*gradb*t35*t97*t13*
661 t39;
662 t108 = 2.0*t1*t13*t33*t39*t6*gradb;
663 t109 = -6.718719623277062E-4*t1*t34*gradb*t35*t59*t13*
664 t33*t38*t39;
665 t110 = -2.666666666666667*t1*gradb*t50*t13*t37*t39;
666 t111 = 7.838506227156574E-4*t34*gradb*t102*t53*t38+6.718719623277062E-4*
667 t34*gradb*t102*t36*t33*t52;
668 t112 = t1*t4*t6*t13*t37*t111*t39;
669 t113 = -6.718719623277062E-4*t1*t34*gradb*t35*t59*t13*
670 t54*t39;
671 t114 = 2.0*t1*t13*t37*t39*t54*t6*gradb;
672 t115 = POW(t58,2.0);
673 t116 = POW(t56,2.0);
674 t117 = t82+t76;
675 t118 = POW(t57,2.0);
676 t119 = t77+t81;
677 t120 = 0.55555555555556*t118*t7*t80+0.55555555555556*
678 t116*t7*t75+0.83333333333333*t119*t47*t7+0.83333333333333*
679 t117*t44*t7;
680 t121 = 3.359359811638531E-4*t1*t34*t35*t59*t48*t49*t39;
681 t122 = 8.958292831036083E-4*t1*t34*t35*t97*t13*t39;
682 t123 = 2.0*
683 t1*t13*t33*t39*t6;
684 t124 = -6.718719623277062E-4*t1*t34*t35*t59*t13*t33*t38*
685 t39;
686 t125 = -1.0*t1*t37*t39*t48*t49*t6;
687 t126 = -2.666666666666667*t1*t50*t13*t37*t39;
688 t127 = 6.718719623277062E-4*t34*t102*t36*t33*t52;
689 t128 = 7.838506227156574E-4*t34*t102*t53*t38;
690 t129 = t128+t127;
691 t130 = t1*t4*t6*t13*t37*t129*t39;
692 t131 = -6.718719623277062E-4*t1*t34*t35*t59*t13*t54*t39;
693 t132 = 2.0*
694 t1*t13*t37*t39*t54*t6;
695 t133 = t132+t131+t130+t126+t125+t124+t123+t122+t121;
696 t134 = 3.359359811638531E-4*t1*t34*t35*t59*t58*t49*t39;
697 t135 = -
698 1.0*t1*t37*t39*t49*t58*t6;
699 t136 = t132+t131+t130+t126+t135+t124+t123+t122+t134;
700 t137 = 1/
701 POW(t5,3.666666666666667);
702 t138 = -0.00201561588698*t1*t34*grada*t102*t59*t13*t39;
703 t139 = 4.514119337620827E-7*t1*t22*grada*t137*t13*t38*
704 t39;
705 t140 = -0.00201561588698*t1*t34*gradb*t102*t59*t13*t39;
706 t141 = 4.514119337620827E-7*t1*t22*gradb*t137*t13*t38*
707 t39;
708 t142 = 1/POW(t5,4.0);
709 t143 = 1/POW(t5,3.666666666666667);
710 t144 = 1/POW(t5,3.333333333333334);
711 t145 = -3.0*t25*t30*t68+6.0*t25*t28*t63*t67-3.0*t30*t31*
712 t67+6.0*t31*t62*t63+(-0.0241277037037*t20*t21*t23*t144-2.189333333333334E-5*
713 t16*t17*t18*t143)*t32-1.0*(-9.046074074074074*t20*t21*t23*
714 t144-1.398518518518519*t16*t17*t18*t143-0.332505*t14*t142)*
715 t28*t30-6.0*POW(t25,3.0)*t28/POW(t29,4.0);
716 t146 = t1*t4*t6*t13*t145*t39;
717 t147 = -4.0*t1*t4*t50*t13*t69*t39;
718 t148 = 1.5*t1*t33*t39*t4*t48*t58*t6*t73;
719 t149 = -1.0*t1*t33*t39*t4*t49*t6*t96;
720 t150 = 9.333333333333332*t1*t4*t85*t13*t33*t39;
721 t151 = 1/POW(t12,7.0);
722 t152 = -2.0*t1*t4*t50*t58*t48*t73*t37*t39;
723 t153 = 1.333333333333333*t1*t4*t50*t96*t49*t37*t39;
724 t154 = 1/
725 POW(t10,1.333333333333333);
726 t155 = 6.0*t142*t8;
727 t156 = 1/POW(t11,1.333333333333333);
728 t157 = -6.0*t142*t8;
729 t158 = -10.37037037037037*t1*t13*t37*t39*t4/POW(t5,4.333333333333333);
730 t159 = t1*
731 t13*t37*t39*t4*t6*(0.00470310373629*t34*t35*t53*t87*t88-0.00403123177397*
732 t34*t35*t36*t69*t33*t88-0.00235155186815*t34*t35*t53*t69*t52+
733 0.00509502904765*t34*t35*t89*t33*t52+6.718719623277062E-4*
734 t34*t35*t36*t145*t52+0.00537808621697*t34*t35*t38/POW(t5,4.166666666666667)+
735 0.00403123177397*POW(t33,3.0)*t34*t35*t36/POW(t37,4.0));
736 t160 = 3.0*
737 t1*t13*t33*t39*t4*t6*t90;
738 t161 = -4.0*t1*t4*t50*t13*t37*t90*t39;
739 t162 = 3.0*t1*t13*t39*t4*t54*t6*t69;
740 t163 = -8.0*t1*t4*t50*t13*t33*t54*t39;
741 t164 = 1.5*t1*t37*t39*t4*t48*t54*t58*t6*t73;
742 t165 = -1.0*t1*t37*t39*t4*t49*t54*t6*t96;
743 t166 = 9.333333333333332*t1*t4*t85*t13*t37*t54*t39;
744 t167 = 3.0*
745 t1*t13*t37*t39*t4*t54*t6*t90;
746 t168 = 3.0*t1*t13*t33*t39*t4*t6*t94;
747 t169 = -4.0*t1*t4*t50*t13*t37*t94*t39;
748 t170 = t1*t13*t37*t39*t4*POW(t54,3.0)*t6;
749 t171 = 1/POW(t5,4.5);
750 t172 = -0.00209026832724*t1*t34*grada*t35*t171*t13*t39;
751 t173 = 2.0*
752 t1*t13*t39*t6*t69*grada;
753 t174 = -5.333333333333333*t1*grada*t50*t13*t33*t39;
754 t175 = -
755 6.718719623277062E-4*t1*t34*grada*t35*t59*t13*t69*t38*t39;
756 t176 = 0.00179165856621*
757 t1*t34*grada*t35*t97*t13*t33*t38*t39;
758 t177 = 6.222222222222221*t1*grada*t85*t13*t37*t39;
759 t178 = t1*t4*t6*t13*t37*(-0.00169834301588*t34*grada*
760 t102*t89*t38-0.00156770124543*t34*grada*t102*t53*t33*t52+6.718719623277062E-4*
761 t34*grada*t102*t36*t69*t52-0.00134374392466*t34*grada*t102*
762 t36*t87*t88)*t39;
763 t179 = -6.718719623277062E-4*t1*t34*grada*t35*t59*t13*
764 t90*t39;
765 t180 = 2.0*t1*t13*t37*t39*t6*t90*grada;
766 t181 = 2.0*t1*t103*t13*t33*t39*t4*t6;
767 t182 = -2.666666666666667*t1*t4*t50*t13*t37*t103*t39;
768 t183 = 0.00179165856621*
769 t1*t34*grada*t35*t97*t13*t54*t39;
770 t184 = 4.0*t1*t13*t33*t39*t54*t6*grada;
771 t185 = -0.00134374392466*t1*t34*grada*t35*t59*t13*t33*
772 t38*t54*t39;
773 t186 = -5.333333333333333*t1*grada*t50*t13*t37*t54*t39;
774 t187 = 2.0*
775 t1*t103*t13*t37*t39*t4*t54*t6;
776 t188 = -6.718719623277062E-4*t1*t34*grada*t35*t59*t13*
777 t94*t39;
778 t189 = 2.0*t1*t13*t37*t39*t6*t94*grada;
779 t190 = -0.00209026832724*t1*t34*gradb*t35*t171*t13*t39;
780 t191 = 2.0*
781 t1*t13*t39*t6*t69*gradb;
782 t192 = -5.333333333333333*t1*gradb*t50*t13*t33*t39;
783 t193 = -
784 6.718719623277062E-4*t1*t34*gradb*t35*t59*t13*t69*t38*t39;
785 t194 = 0.00179165856621*
786 t1*t34*gradb*t35*t97*t13*t33*t38*t39;
787 t195 = 6.222222222222221*t1*gradb*t85*t13*t37*t39;
788 t196 = t1*t4*t6*t13*t37*(-0.00169834301588*t34*gradb*
789 t102*t89*t38-0.00156770124543*t34*gradb*t102*t53*t33*t52+6.718719623277062E-4*
790 t34*gradb*t102*t36*t69*t52-0.00134374392466*t34*gradb*t102*
791 t36*t87*t88)*t39;
792 t197 = -6.718719623277062E-4*t1*t34*gradb*t35*t59*t13*
793 t90*t39;
794 t198 = 2.0*t1*t13*t37*t39*t6*t90*gradb;
795 t199 = 2.0*t1*t111*t13*t33*t39*t4*t6;
796 t200 = -2.666666666666667*t1*t4*t50*t13*t37*t111*t39;
797 t201 = 0.00179165856621*
798 t1*t34*gradb*t35*t97*t13*t54*t39;
799 t202 = 4.0*t1*t13*t33*t39*t54*t6*gradb;
800 t203 = -0.00134374392466*t1*t34*gradb*t35*t59*t13*t33*
801 t38*t54*t39;
802 t204 = -5.333333333333333*t1*gradb*t50*t13*t37*t54*t39;
803 t205 = 2.0*
804 t1*t111*t13*t37*t39*t4*t54*t6;
805 t206 = -6.718719623277062E-4*t1*t34*gradb*t35*t59*t13*
806 t94*t39;
807 t207 = 2.0*t1*t13*t37*t39*t6*t94*gradb;
808 t208 = -0.00209026832724*t1*t34*t35*t171*t13*t39;
809 t209 = 2.0*t1*t13*t39*t6*t69;
810 t210 = -5.333333333333333*t1*t50*t13*t33*t39;
811 t211 = -6.718719623277062E-4*t1*t34*t35*t59*t13*t69*t38*
812 t39;
813 t212 = 0.00179165856621*t1*t34*t35*t97*t13*t33*t38*t39;
814 t213 = 6.222222222222221*
815 t1*t85*t13*t37*t39;
816 t214 = t1*t4*t6*t13*t37*(-0.00169834301588*t34*t102*t89*
817 t38-0.00156770124543*t34*t102*t53*t33*t52+6.718719623277062E-4*
818 t34*t102*t36*t69*t52-0.00134374392466*t34*t102*t36*t87*t88)*
819 t39;
820 t215 = -6.718719623277062E-4*t1*t34*t35*t59*t13*t90*t39;
821 t216 = 2.0*
822 t1*t13*t37*t39*t6*t90;
823 t217 = 2.0*t1*t129*t13*t33*t39*t4*t6;
824 t218 = -2.666666666666667*t1*t4*t50*t13*t37*t129*t39;
825 t219 = 0.00179165856621*
826 t1*t34*t35*t97*t13*t54*t39;
827 t220 = 4.0*t1*t13*t33*t39*t54*t6;
828 t221 = -0.00134374392466*t1*t34*t35*t59*t13*t33*t38*t54*
829 t39;
830 t222 = -5.333333333333333*t1*t50*t13*t37*t54*t39;
831 t223 = 2.0*t1*t129*t13*t37*t39*t4*t54*t6;
832 t224 = -6.718719623277062E-4*t1*t34*t35*t59*t13*t94*t39;
833 t225 = 2.0*
834 t1*t13*t37*t39*t6*t94;
835 t226 = 0.00268748784931*t1*t34*t2*t102*t97*t13*t39;
836 t227 = 4.514119337620827E-7*t1*t22*t2*t137*t13*t33*t52*
837 t39;
838 t228 = 1/POW(t5,4.666666666666667);
839 t229 = -6.018825783494436E-7*t1*t22*t2*t228*t13*t38*t39;
840 t230 = -
841 0.00201561588698*t1*t34*t2*t102*t59*t13*t33*t38*t39;
842 t231 = 1/POW(t35,3.0);
843 t232 = t1*t4*t6*t13*t37*(t128-7.838506227156574E-4*t34*
844 t2*t231*t53*t38+t127-6.718719623277062E-4*t34*t2*t231*t36*
845 t33*t52)*t39;
846 t233 = -0.00134374392466*t1*t34*grada*t35*t59*t13*t103*
847 t39;
848 t234 = 4.0*t1*t103*t13*t37*t39*t6*grada;
849 t235 = -0.00201561588698*t1*t34*t2*t102*t59*t13*t54*t39;
850 t236 = 4.514119337620827E-7*t1*t22*t2*t137*t13*t38*t54*
851 t39;
852 t237 = 0.00268748784931*t1*t34*grada*gradb*t102*t97*t13*
853 t39;
854 t238 = 4.514119337620827E-7*t1*t22*grada*gradb*t137*t13*
855 t33*t52*t39;
856 t239 = -6.018825783494436E-7*t1*t22*grada*gradb*t228*
857 t13*t38*t39;
858 t240 = -0.00201561588698*t1*t34*grada*gradb*t102*t59*
859 t13*t33*t38*t39;
860 t241 = t1*t4*t6*t13*t37*(-7.838506227156574E-4*t34*grada*
861 gradb*t231*t53*t38-6.718719623277062E-4*t34*grada*gradb*t231*
862 t36*t33*t52)*t39;
863 t242 = -6.718719623277062E-4*t1*t34*gradb*t35*t59*t13*
864 t103*t39;
865 t243 = 2.0*t1*t103*t13*t37*t39*t6*gradb;
866 t244 = -6.718719623277062E-4*t1*t34*grada*t35*t59*t13*
867 t111*t39;
868 t245 = 2.0*t1*t111*t13*t37*t39*t6*grada;
869 t246 = -0.00201561588698*t1*t34*grada*gradb*t102*t59*
870 t13*t54*t39;
871 t247 = 4.514119337620827E-7*t1*t22*grada*gradb*t137*t13*
872 t38*t54*t39;
873 t248 = 0.00268748784931*t1*t34*t3*t102*t97*t13*t39;
874 t249 = 4.514119337620827E-7*t1*t22*t3*t137*t13*t33*t52*
875 t39;
876 t250 = -6.018825783494436E-7*t1*t22*t3*t228*t13*t38*t39;
877 t251 = -
878 0.00201561588698*t1*t34*t3*t102*t59*t13*t33*t38*t39;
879 t252 = t1*t4*t6*t13*t37*(t128-7.838506227156574E-4*t34*
880 t3*t231*t53*t38+t127-6.718719623277062E-4*t34*t3*t231*t36*
881 t33*t52)*t39;
882 t253 = -0.00134374392466*t1*t34*gradb*t35*t59*t13*t111*
883 t39;
884 t254 = 4.0*t1*t111*t13*t37*t39*t6*gradb;
885 t255 = -0.00201561588698*t1*t34*t3*t102*t59*t13*t54*t39;
886 t256 = 4.514119337620827E-7*t1*t22*t3*t137*t13*t38*t54*
887 t39;
888 t257 = -6.0*t64;
889 t258 = 6.0*t64;
890 t259 = t132+t236+t131+t235+t234+t233+t232+t126+t135+t124+
891 t230+t229-2.257059668810414E-7*t1*t22*t2*t137*t58*t49*t38*
892 t39+t227+t123+t122+t226+t134+0.00100780794349*t1*t34*t2*t102*
893 t59*t58*t49*t39;
894 t260 = POW(grada,3.0);
895 t261 = 1.772453850905516;
896 t262 = 1/POW(t5,4.833333333333334);
897 t263 = 1/t4;
898 t264 = POW(gradb,3.0);
899
900 /* code */
901 dfdra = -0.5*t1*t37*t39*t4*t48*t49*t6+t55+t51+t40;
902 dfdrb = -0.5*t1*t37*t39*t4*t49*t58*t6+t55+t51+t40;
903 dfdga = 2.0*t1*t13*t37*t39*t6*grada-6.718719623277062E-4*
904 t1*t34*grada*t35*t59*t13*t39;
905 dfdgb = 2.0*t1*t13*t37*t39*t6*gradb-6.718719623277062E-4*
906 t1*t34*gradb*t35*t59*t13*t39;
907 dfdgab = t61+t60;
908 d2fdrara = t95+t93+t92+t91+t86-0.5*t1*t37*t39*t4*t49*
909 t6*t84+0.75*t1*t37*t39*t4*t6*t72*t73+t71+t70-1.0*t1*t37*t39*
910 t4*t48*t49*t54*t6-1.0*t1*t33*t39*t4*t48*t49*t6+1.333333333333333*
911 t1*t4*t50*t48*t49*t37*t39;
912 d2fdrarb = -0.5*t1*t37*t39*t4*t49*t6*t96+t95+t93+t92+
913 t91+t86+0.75*t1*t37*t39*t4*t48*t58*t6*t73+t71+t70-0.5*t1*t37*
914 t39*t4*t49*t54*t58*t6-0.5*t1*t33*t39*t4*t49*t58*t6-0.5*t1*
915 t37*t39*t4*t48*t49*t54*t6-0.5*t1*t33*t39*t4*t48*t49*t6+0.66666666666667*
916 t1*t4*t50*t58*t49*t37*t39+0.66666666666667*t1*t4*t50*t48*t49*
917 t37*t39;
918 d2fdraga = -1.0*t1*t37*t39*t48*t49*t6*grada+t99+t98+3.359359811638531E-4*
919 t1*t34*grada*t35*t59*t48*t49*t39+t106+t105+t104+t101+t100;
920 d2fdragb = -
921 1.0*t1*t37*t39*t48*t49*t6*gradb+3.359359811638531E-4*t1*t34*
922 gradb*t35*t59*t48*t49*t39+t114+t113+t112+t110+t109+t108+t107;
923 d2fdrbrb = t95+
924 t93+t92+t91+t86+0.75*t1*t115*t37*t39*t4*t6*t73+t71+t70-1.0*
925 t1*t37*t39*t4*t49*t54*t58*t6-1.0*t1*t33*t39*t4*t49*t58*t6-
926 0.5*t1*t120*t37*t39*t4*t49*t6+1.333333333333333*t1*t4*t50*
927 t58*t49*t37*t39;
928 d2fdraab = t133;
929 d2fdrbab = t136;
930 d2fdgaga = t61+4.514119337620827E-7*t1*t22*t2*t137*t13*
931 t38*t39+t60-0.00201561588698*t1*t34*t2*t102*t59*t13*t39;
932 d2fdgbgb = t61+
933 4.514119337620827E-7*t1*t22*t3*t137*t13*t38*t39+t60-0.00201561588698*
934 t1*t34*t3*t102*t59*t13*t39;
935 d2fdrbga = -1.0*t1*t37*t39*t49*t58*t6*grada+t99+t98+3.359359811638531E-4*
936 t1*t34*grada*t35*t59*t58*t49*t39+t106+t105+t104+t101+t100;
937 d2fdrbgb = -
938 1.0*t1*t37*t39*t49*t58*t6*gradb+3.359359811638531E-4*t1*t34*
939 gradb*t35*t59*t58*t49*t39+t114+t113+t112+t110+t109+t108+t107;
940 d2fdgagb = 4.514119337620827E-7*t1*t22*grada*gradb*t137*
941 t13*t38*t39-0.00201561588698*t1*t34*grada*gradb*t102*t59*t13*
942 t39;
943 d3fdrararb = 1.5*t1*t37*t39*t4*t48*t6*t73*t96-0.5*t1*
944 t37*t39*t4*t49*t58*t6*t94-1.0*t1*t37*t39*t4*t48*t49*t6*t94-
945 0.5*t1*t37*t39*t4*t49*t58*t6*t90-1.0*t1*t37*t39*t4*t48*t49*
946 t6*t90+0.75*t1*t37*t39*t4*t58*t6*t73*t84-0.5*t1*t37*t39*t4*
947 t49*t54*t6*t84-0.5*t1*t33*t39*t4*t49*t6*t84-0.5*t1*t37*t39*
948 t4*t49*t6*(0.55555555555556*t57*t7*t80*t83+2.222222222222222*
949 t46*t64*t7*t8*t80-2.222222222222222*t43*t64*t7*t75*t8-0.18518518518519*
950 t156*t57*t7*t79+0.55555555555556*t56*t7*t75*t78-0.18518518518519*
951 t154*t56*t7*t74+0.83333333333333*t47*(2.0*t64+t157)*t7+0.83333333333333*
952 t44*(t155-2.0*t64)*t7)+0.75*t1*t37*t39*t4*t54*t6*t72*t73+0.75*
953 t1*t33*t39*t4*t6*t72*t73-1.875*t1*t151*t37*t39*t4*t58*t6*t72-
954 0.5*t1*t39*t4*t49*t58*t6*t69-1.0*t1*t39*t4*t48*t49*t6*t69-
955 1.0*t1*t33*t39*t4*t49*t54*t58*t6-2.0*t1*t33*t39*t4*t48*t49*
956 t54*t6+1.333333333333333*t1*t4*t50*t58*t49*t37*t54*t39+2.666666666666667*
957 t1*t4*t50*t48*t49*t37*t54*t39-1.0*t1*t4*t50*t72*t73*t37*t39+
958 0.66666666666667*t1*t4*t50*t84*t49*t37*t39-1.555555555555555*
959 t1*t4*t85*t58*t49*t37*t39-3.111111111111111*t1*t4*t85*t48*
960 t49*t37*t39+1.333333333333333*t1*t4*t50*t58*t49*t33*t39+2.666666666666667*
961 t1*t4*t50*t48*t49*t33*t39+t170+t169+t168+t167+t166+t165+t164+
962 t163+t162+t161+t160+t159+t158+t153+t152+t150+t149+t148+t147+
963 t146;
964 d3fdraraga = -1.0*t1*t37*t39*t49*t6*t84*grada+1.5*t1*
965 t37*t39*t6*t72*t73*grada-2.0*t1*t37*t39*t48*t49*t54*t6*grada-
966 2.0*t1*t33*t39*t48*t49*t6*grada-1.0*t1*t103*t37*t39*t4*t48*
967 t49*t6-5.039039717457796E-4*t1*t34*grada*t35*t59*t72*t73*t39+
968 6.718719623277062E-4*t1*t34*grada*t35*t59*t48*t49*t54*t39+
969 3.359359811638531E-4*t1*t34*grada*t35*t59*t84*t49*t39-8.958292831036083E-4*
970 t1*t34*grada*t35*t97*t48*t49*t39+6.718719623277062E-4*t1*t34*
971 grada*t35*t59*t48*t49*t33*t38*t39+2.666666666666667*t1*grada*
972 t50*t48*t49*t37*t39+t189+t188+t187+t186+t185+t184+t183+t182+
973 t181+t180+t179+t178+t177+t176+t175+t174+t173+t172;
974 d3fdraragb = -1.0*t1*t37*t39*t49*t6*t84*gradb+1.5*t1*
975 t37*t39*t6*t72*t73*gradb-2.0*t1*t37*t39*t48*t49*t54*t6*gradb-
976 2.0*t1*t33*t39*t48*t49*t6*gradb-1.0*t1*t111*t37*t39*t4*t48*
977 t49*t6-5.039039717457796E-4*t1*t34*gradb*t35*t59*t72*t73*t39+
978 6.718719623277062E-4*t1*t34*gradb*t35*t59*t48*t49*t54*t39+
979 3.359359811638531E-4*t1*t34*gradb*t35*t59*t84*t49*t39-8.958292831036083E-4*
980 t1*t34*gradb*t35*t97*t48*t49*t39+6.718719623277062E-4*t1*t34*
981 gradb*t35*t59*t48*t49*t33*t38*t39+2.666666666666667*t1*gradb*
982 t50*t48*t49*t37*t39+t207+t206+t205+t204+t203+t202+t201+t200+
983 t199+t198+t197+t196+t195+t194+t193+t192+t191+t190;
984 d3fdrbrbab = 1.5*t1*t115*t37*t39*t6*t73-2.0*t1*t37*t39*
985 t49*t54*t58*t6-1.0*t1*t129*t37*t39*t4*t49*t58*t6-2.0*t1*t33*
986 t39*t49*t58*t6-1.0*t1*t120*t37*t39*t49*t6-5.039039717457796E-4*
987 t1*t34*t35*t59*t115*t73*t39+6.718719623277062E-4*t1*t34*t35*
988 t59*t58*t49*t54*t39-8.958292831036083E-4*t1*t34*t35*t97*t58*
989 t49*t39+3.359359811638531E-4*t1*t34*t35*t59*t120*t49*t39+6.718719623277062E-4*
990 t1*t34*t35*t59*t58*t49*t33*t38*t39+2.666666666666667*t1*t50*
991 t58*t49*t37*t39+t225+t224+t223+t222+t221+t220+t219+t218+t217+
992 t216+t215+t214+t213+t212+t211+t210+t209+t208;
993 d3fdraraab = -1.0*t1*t37*t39*t49*t6*t84+1.5*t1*t37*t39*
994 t6*t72*t73-2.0*t1*t37*t39*t48*t49*t54*t6-1.0*t1*t129*t37*t39*
995 t4*t48*t49*t6-2.0*t1*t33*t39*t48*t49*t6-5.039039717457796E-4*
996 t1*t34*t35*t59*t72*t73*t39+6.718719623277062E-4*t1*t34*t35*
997 t59*t48*t49*t54*t39+3.359359811638531E-4*t1*t34*t35*t59*t84*
998 t49*t39-8.958292831036083E-4*t1*t34*t35*t97*t48*t49*t39+6.718719623277062E-4*
999 t1*t34*t35*t59*t48*t49*t33*t38*t39+2.666666666666667*t1*t50*
1000 t48*t49*t37*t39+t225+t224+t223+t222+t221+t220+t219+t218+t217+
1001 t216+t215+t214+t213+t212+t211+t210+t209+t208;
1002 d3fdrarbrb = 1.5*t1*t37*t39*t4*t58*t6*t73*t96-1.0*t1*
1003 t37*t39*t4*t49*t58*t6*t94-0.5*t1*t37*t39*t4*t48*t49*t6*t94-
1004 1.0*t1*t37*t39*t4*t49*t58*t6*t90-0.5*t1*t37*t39*t4*t48*t49*
1005 t6*t90-0.5*t1*t37*t39*t4*t49*t6*(2.222222222222222*t57*t64*
1006 t7*t8*t80+0.55555555555556*t119*t46*t7*t80-2.222222222222222*
1007 t56*t64*t7*t75*t8-5.0*t142*t47*t7*t8+5.0*t142*t44*t7*t8+0.55555555555556*
1008 t117*t43*t7*t75-1.666666666666667*t47*t64*t7+1.666666666666667*
1009 t44*t64*t7-0.18518518518519*t118*t156*t46*t7-0.18518518518519*
1010 t116*t154*t43*t7)+0.75*t1*t115*t37*t39*t4*t54*t6*t73+0.75*
1011 t1*t120*t37*t39*t4*t48*t6*t73+0.75*t1*t115*t33*t39*t4*t6*t73-
1012 1.0*t1*t39*t4*t49*t58*t6*t69-0.5*t1*t39*t4*t48*t49*t6*t69-
1013 2.0*t1*t33*t39*t4*t49*t54*t58*t6-1.0*t1*t33*t39*t4*t48*t49*
1014 t54*t6-0.5*t1*t120*t37*t39*t4*t49*t54*t6-0.5*t1*t120*t33*t39*
1015 t4*t49*t6-1.875*t1*t115*t151*t37*t39*t4*t48*t6+2.666666666666667*
1016 t1*t4*t50*t58*t49*t37*t54*t39+1.333333333333333*t1*t4*t50*
1017 t48*t49*t37*t54*t39-1.0*t1*t4*t50*t115*t73*t37*t39-3.111111111111111*
1018 t1*t4*t85*t58*t49*t37*t39-1.555555555555555*t1*t4*t85*t48*
1019 t49*t37*t39+0.66666666666667*t1*t4*t50*t120*t49*t37*t39+2.666666666666667*
1020 t1*t4*t50*t58*t49*t33*t39+1.333333333333333*t1*t4*t50*t48*
1021 t49*t33*t39+t170+t169+t168+t167+t166+t165+t164+t163+t162+t161+
1022 t160+t159+t158+t153+t152+t150+t149+t148+t147+t146;
1023 d3fdrarbga = -1.0*t1*t37*t39*t49*t6*t96*grada+1.5*t1*
1024 t37*t39*t48*t58*t6*t73*grada-1.0*t1*t37*t39*t49*t54*t58*t6*
1025 grada-1.0*t1*t33*t39*t49*t58*t6*grada-1.0*t1*t37*t39*t48*t49*
1026 t54*t6*grada-1.0*t1*t33*t39*t48*t49*t6*grada-0.5*t1*t103*t37*
1027 t39*t4*t49*t58*t6-0.5*t1*t103*t37*t39*t4*t48*t49*t6-5.039039717457796E-4*
1028 t1*t34*grada*t35*t59*t58*t48*t73*t39+3.359359811638531E-4*
1029 t1*t34*grada*t35*t59*t58*t49*t54*t39+3.359359811638531E-4*
1030 t1*t34*grada*t35*t59*t48*t49*t54*t39+3.359359811638531E-4*
1031 t1*t34*grada*t35*t59*t96*t49*t39-4.479146415518041E-4*t1*t34*
1032 grada*t35*t97*t58*t49*t39-4.479146415518041E-4*t1*t34*grada*
1033 t35*t97*t48*t49*t39+3.359359811638531E-4*t1*t34*grada*t35*
1034 t59*t58*t49*t33*t38*t39+3.359359811638531E-4*t1*t34*grada*
1035 t35*t59*t48*t49*t33*t38*t39+1.333333333333333*t1*grada*t50*
1036 t58*t49*t37*t39+1.333333333333333*t1*grada*t50*t48*t49*t37*
1037 t39+t189+t188+t187+t186+t185+t184+t183+t182+t181+t180+t179+
1038 t178+t177+t176+t175+t174+t173+t172;
1039 d3fdrarbgb = -1.0*t1*t37*t39*t49*t6*t96*gradb+1.5*t1*
1040 t37*t39*t48*t58*t6*t73*gradb-1.0*t1*t37*t39*t49*t54*t58*t6*
1041 gradb-1.0*t1*t33*t39*t49*t58*t6*gradb-1.0*t1*t37*t39*t48*t49*
1042 t54*t6*gradb-1.0*t1*t33*t39*t48*t49*t6*gradb-0.5*t1*t111*t37*
1043 t39*t4*t49*t58*t6-0.5*t1*t111*t37*t39*t4*t48*t49*t6-5.039039717457796E-4*
1044 t1*t34*gradb*t35*t59*t58*t48*t73*t39+3.359359811638531E-4*
1045 t1*t34*gradb*t35*t59*t58*t49*t54*t39+3.359359811638531E-4*
1046 t1*t34*gradb*t35*t59*t48*t49*t54*t39+3.359359811638531E-4*
1047 t1*t34*gradb*t35*t59*t96*t49*t39-4.479146415518041E-4*t1*t34*
1048 gradb*t35*t97*t58*t49*t39-4.479146415518041E-4*t1*t34*gradb*
1049 t35*t97*t48*t49*t39+3.359359811638531E-4*t1*t34*gradb*t35*
1050 t59*t58*t49*t33*t38*t39+3.359359811638531E-4*t1*t34*gradb*
1051 t35*t59*t48*t49*t33*t38*t39+1.333333333333333*t1*gradb*t50*
1052 t58*t49*t37*t39+1.333333333333333*t1*gradb*t50*t48*t49*t37*
1053 t39+t207+t206+t205+t204+t203+t202+t201+t200+t199+t198+t197+
1054 t196+t195+t194+t193+t192+t191+t190;
1055 d3fdrarbab = -1.0*t1*t37*t39*t49*t6*t96+1.5*t1*t37*t39*
1056 t48*t58*t6*t73-1.0*t1*t37*t39*t49*t54*t58*t6-0.5*t1*t129*t37*
1057 t39*t4*t49*t58*t6-1.0*t1*t33*t39*t49*t58*t6-1.0*t1*t37*t39*
1058 t48*t49*t54*t6-0.5*t1*t129*t37*t39*t4*t48*t49*t6-1.0*t1*t33*
1059 t39*t48*t49*t6-5.039039717457796E-4*t1*t34*t35*t59*t58*t48*
1060 t73*t39+3.359359811638531E-4*t1*t34*t35*t59*t58*t49*t54*t39+
1061 3.359359811638531E-4*t1*t34*t35*t59*t48*t49*t54*t39+3.359359811638531E-4*
1062 t1*t34*t35*t59*t96*t49*t39-4.479146415518041E-4*t1*t34*t35*
1063 t97*t58*t49*t39-4.479146415518041E-4*t1*t34*t35*t97*t48*t49*
1064 t39+3.359359811638531E-4*t1*t34*t35*t59*t58*t49*t33*t38*t39+
1065 3.359359811638531E-4*t1*t34*t35*t59*t48*t49*t33*t38*t39+1.333333333333333*
1066 t1*t50*t58*t49*t37*t39+1.333333333333333*t1*t50*t48*t49*t37*
1067 t39+t225+t224+t223+t222+t221+t220+t219+t218+t217+t216+t215+
1068 t214+t213+t212+t211+t210+t209+t208;
1069 d3fdragaga = t132+t236+t131+t235+t234+t233+t232+t126+
1070 t125+t124+t230+t229-2.257059668810414E-7*t1*t22*t2*t137*t48*
1071 t49*t38*t39+t227+t123+t122+t226+t121+0.00100780794349*t1*t34*
1072 t2*t102*t59*t48*t49*t39;
1073 d3fdragagb = t247+t246+t245+t244+t243+t242+t241+t240+
1074 t239-2.257059668810414E-7*t1*t22*grada*gradb*t137*t48*t49*
1075 t38*t39+t238+t237+0.00100780794349*t1*t34*grada*gradb*t102*
1076 t59*t48*t49*t39;
1077 d3fdrbgagb = t247+t246+t245+t244+t243+t242+t241+t240+
1078 t239-2.257059668810414E-7*t1*t22*grada*gradb*t137*t58*t49*
1079 t38*t39+t238+t237+0.00100780794349*t1*t34*grada*gradb*t102*
1080 t59*t58*t49*t39;
1081 d3fdragbgb = t132+t256+t131+t255+t254+t253+t252+t126+
1082 t125+t124+t251+t250-2.257059668810414E-7*t1*t22*t3*t137*t48*
1083 t49*t38*t39+t249+t123+t122+t248+t121+0.00100780794349*t1*t34*
1084 t3*t102*t59*t48*t49*t39;
1085 d3fdrarara = -1.5*t1*t37*t39*t4*t48*t49*t6*t94-1.5*t1*
1086 t37*t39*t4*t48*t49*t6*t90+2.25*t1*t37*t39*t4*t48*t6*t73*t84-
1087 1.5*t1*t37*t39*t4*t49*t54*t6*t84-1.5*t1*t33*t39*t4*t49*t6*
1088 t84-0.5*t1*t37*t39*t4*t49*t6*(1.666666666666667*t46*t7*t80*
1089 t83+1.666666666666667*t43*t7*t75*t78+0.83333333333333*(t258+
1090 t157)*t47*t7-0.18518518518519*t156*POW(t46,3.0)*t7+0.83333333333333*
1091 (t257+t155)*t44*t7-0.18518518518519*t154*POW(t43,3.0)*t7)+
1092 2.25*t1*t37*t39*t4*t54*t6*t72*t73+2.25*t1*t33*t39*t4*t6*t72*
1093 t73-1.5*t1*t39*t4*t48*t49*t6*t69-3.0*t1*t33*t39*t4*t48*t49*
1094 t54*t6-1.875*t1*t151*t37*t39*t4*POW(t48,3.0)*t6+4.0*t1*t4*
1095 t50*t48*t49*t37*t54*t39-3.0*t1*t4*t50*t72*t73*t37*t39+2.0*
1096 t1*t4*t50*t84*t49*t37*t39-4.666666666666666*t1*t4*t85*t48*
1097 t49*t37*t39+4.0*t1*t4*t50*t48*t49*t33*t39+t170+t169+t168+t167+
1098 t166+t163+t162+t161+t160+t159+t158+t150+t147+t146;
1099 d3fdrbrbrb = -1.5*t1*t37*t39*t4*t49*t58*t6*t94-1.5*t1*
1100 t37*t39*t4*t49*t58*t6*t90-0.5*t1*t37*t39*t4*t49*t6*(1.666666666666667*
1101 t119*t57*t7*t80+1.666666666666667*t117*t56*t7*t75-0.18518518518519*
1102 t156*POW(t57,3.0)*t7-0.18518518518519*t154*POW(t56,3.0)*t7+
1103 0.83333333333333*(t257+t157)*t47*t7+0.83333333333333*(t258+
1104 t155)*t44*t7)+2.25*t1*t120*t37*t39*t4*t58*t6*t73+2.25*t1*t115*
1105 t37*t39*t4*t54*t6*t73+2.25*t1*t115*t33*t39*t4*t6*t73-1.5*t1*
1106 t39*t4*t49*t58*t6*t69-1.875*t1*t151*t37*t39*t4*POW(t58,3.0)*
1107 t6-3.0*t1*t33*t39*t4*t49*t54*t58*t6-1.5*t1*t120*t37*t39*t4*
1108 t49*t54*t6-1.5*t1*t120*t33*t39*t4*t49*t6+4.0*t1*t4*t50*t58*
1109 t49*t37*t54*t39-3.0*t1*t4*t50*t115*t73*t37*t39-4.666666666666666*
1110 t1*t4*t85*t58*t49*t37*t39+2.0*t1*t4*t50*t120*t49*t37*t39+4.0*
1111 t1*t4*t50*t58*t49*t33*t39+t170+t169+t168+t167+t166+t163+t162+
1112 t161+t160+t159+t158+t150+t147+t146;
1113 d3fdrbrbga = 1.5*t1*t115*t37*t39*t6*t73*grada-2.0*t1*
1114 t37*t39*t49*t54*t58*t6*grada-2.0*t1*t33*t39*t49*t58*t6*grada-
1115 1.0*t1*t120*t37*t39*t49*t6*grada-1.0*t1*t103*t37*t39*t4*t49*
1116 t58*t6-5.039039717457796E-4*t1*t34*grada*t35*t59*t115*t73*
1117 t39+6.718719623277062E-4*t1*t34*grada*t35*t59*t58*t49*t54*
1118 t39-8.958292831036083E-4*t1*t34*grada*t35*t97*t58*t49*t39+
1119 3.359359811638531E-4*t1*t34*grada*t35*t59*t120*t49*t39+6.718719623277062E-4*
1120 t1*t34*grada*t35*t59*t58*t49*t33*t38*t39+2.666666666666667*
1121 t1*grada*t50*t58*t49*t37*t39+t189+t188+t187+t186+t185+t184+
1122 t183+t182+t181+t180+t179+t178+t177+t176+t175+t174+t173+t172;
1123 d3fdrbrbgb = 1.5*
1124 t1*t115*t37*t39*t6*t73*gradb-2.0*t1*t37*t39*t49*t54*t58*t6*
1125 gradb-2.0*t1*t33*t39*t49*t58*t6*gradb-1.0*t1*t120*t37*t39*
1126 t49*t6*gradb-1.0*t1*t111*t37*t39*t4*t49*t58*t6-5.039039717457796E-4*
1127 t1*t34*gradb*t35*t59*t115*t73*t39+6.718719623277062E-4*t1*
1128 t34*gradb*t35*t59*t58*t49*t54*t39-8.958292831036083E-4*t1*
1129 t34*gradb*t35*t97*t58*t49*t39+3.359359811638531E-4*t1*t34*
1130 gradb*t35*t59*t120*t49*t39+6.718719623277062E-4*t1*t34*gradb*
1131 t35*t59*t58*t49*t33*t38*t39+2.666666666666667*t1*gradb*t50*
1132 t58*t49*t37*t39+t207+t206+t205+t204+t203+t202+t201+t200+t199+
1133 t198+t197+t196+t195+t194+t193+t192+t191+t190;
1134 d3fdrbgaga = t259;
1135 d3fdrbgbgb = t132+t256+t131+t255+t254+t253+t252+t126+
1136 t135+t124+t251+t250-2.257059668810414E-7*t1*t22*t3*t137*t58*
1137 t49*t38*t39+t249+t123+t122+t248+t134+0.00100780794349*t1*t34*
1138 t3*t102*t59*t58*t49*t39;
1139 d3fdgagaga = 1.354235801286248E-6*t1*t22*t260*t263*t137*
1140 t13*t38*t39+1.354235801286248E-6*t1*t22*grada*t137*t13*t38*
1141 t39-3.03291021754875E-10*t1*t261*t260*t102*t262*t13*t52*t39-
1142 0.00604684766095*t1*t34*grada*t102*t59*t13*t39+0.00201561588698*
1143 t1*t34*t260*t231*t59*t13*t39;
1144 d3fdgagagb = 1.354235801286248E-6*t1*t22*t2*gradb*t263*
1145 t137*t13*t38*t39+t141-3.03291021754875E-10*t1*t261*t2*gradb*
1146 t102*t262*t13*t52*t39+t140+0.00201561588698*t1*t34*t2*gradb*
1147 t231*t59*t13*t39;
1148 d3fdgagbgb = 1.354235801286248E-6*t1*t22*grada*t3*t263*
1149 t137*t13*t38*t39+t139-3.03291021754875E-10*t1*t261*grada*t3*
1150 t102*t262*t13*t52*t39+t138+0.00201561588698*t1*t34*grada*t3*
1151 t231*t59*t13*t39;
1152 d3fdgbgbgb = 1.354235801286248E-6*t1*t22*t264*t263*t137*
1153 t13*t38*t39+1.354235801286248E-6*t1*t22*gradb*t137*t13*t38*
1154 t39-3.03291021754875E-10*t1*t261*t264*t102*t262*t13*t52*t39-
1155 0.00604684766095*t1*t34*gradb*t102*t59*t13*t39+0.00201561588698*
1156 t1*t34*t264*t231*t59*t13*t39;
1157
1158
1159 ds->df1000 += factor*dfdra;
1160 ds->df0100 += factor*dfdrb;
1161 ds->df0010 += factor*dfdga;
1162 ds->df0001 += factor*dfdgb;
1163 ds->df00001+= factor*dfdgab;
1164
1165 ds->df2000 += factor*d2fdrara;
1166 ds->df1100 += factor*d2fdrarb;
1167 ds->df1010 += factor*d2fdraga;
1168 ds->df1001 += factor*d2fdragb;
1169 ds->df0200 += factor*d2fdrbrb;
1170 ds->df10001+= factor*d2fdraab;
1171 ds->df01001+= factor*d2fdrbab;
1172 ds->df0020 += factor*d2fdgaga;
1173 ds->df0011 += factor*d2fdgagb;
1174 ds->df0002 += factor*d2fdgbgb;
1175 ds->df0110 += factor*d2fdrbga;
1176 ds->df0101 += factor*d2fdrbgb;
1177
1178 ds->df2010 += factor*d3fdraraga;
1179 ds->df2001 += factor*d3fdraragb;
1180 ds->df1101 += factor*d3fdrarbgb;
1181 ds->df11001 += factor*d3fdrarbab;
1182 ds->df1020 += factor*d3fdragaga;
1183 ds->df1011 += factor*d3fdragagb;
1184 ds->df0111 += factor*d3fdrbgagb;
1185 ds->df1002 += factor*d3fdragbgb;
1186 ds->df3000 += factor*d3fdrarara;
1187 ds->df2100 += factor*d3fdrararb;
1188 ds->df20001 += factor*d3fdraraab;
1189 ds->df02001 += factor*d3fdrbrbab;
1190 ds->df1200 += factor*d3fdrarbrb;
1191 ds->df1110 += factor*d3fdrarbga;
1192 ds->df0300 += factor*d3fdrbrbrb;
1193 ds->df0210 += factor*d3fdrbrbga;
1194 ds->df0201 += factor*d3fdrbrbgb;
1195 ds->df0120 += factor*d3fdrbgaga;
1196 ds->df0102 += factor*d3fdrbgbgb;
1197 ds->df0030 += factor*d3fdgagaga;
1198 ds->df0021 += factor*d3fdgagagb;
1199 ds->df0012 += factor*d3fdgagbgb;
1200 ds->df0003 += factor*d3fdgbgbgb;
1201 }
1202