1 /* Automatically generated functional code: pbex
2    Maxima input:
3     >> pi:3.14159265358979312;
4     >>
5     >> xa:sqrt(grada*grada)/rhoa^(4/3);
6     >> xb:sqrt(gradb*gradb)/rhob^(4/3);
7     >>
8     >> R:0.804;
9     >> d:0.066725;
10     >> mu:d*pi^2/3;
11     >> Sa:xa/(2*(6*pi^2)^(1/3));
12     >> Sb:xb/(2*(6*pi^2)^(1/3));
13     >>
14     >> F(S):=1+R-R/(1+mu*S^2/R);
15     >> Ea(n):=-3/(4*pi)*(3*pi^2)^(1/3)*n^(4/3)*F(Sa);
16     >> Eb(n):=-3/(4*pi)*(3*pi^2)^(1/3)*n^(4/3)*F(Sb);
17     >>
18     >> K(rhoa,grada,rhob,gradb,gradab):=0.5*(Ea(2*rhoa)+Eb(2*rhob));
19 */
20 
21 // add "extern Functional pbexFunctional;" to 'functionals.h'
22 // add "&pbexFunctional," to 'functionals.c'
23 // add "fun-pbex.c" to 'Makefile.in'
24 
25 #include <math.h>
26 #include <stddef.h>
27 #include "general.h"
28 
29 #define __CVERSION__
30 
31 #include "functionals.h"
32 #define LOG log
33 #define ABS fabs
34 #define ASINH asinh
35 #define SQRT sqrt
36 
37 /* INTERFACE PART */
pbex_isgga(void)38 static integer pbex_isgga(void) {return 1;}
39 static integer pbex_read(const char* conf_line);
40 static real pbex_energy(const FunDensProp* dp);
41 static void pbex_first(FunFirstFuncDrv *ds, real factor,
42                            const FunDensProp* dp);
43 static void pbex_second(FunSecondFuncDrv *ds, real factor,
44                             const FunDensProp* dp);
45 static void pbex_third(FunThirdFuncDrv *ds, real factor,
46                            const FunDensProp* dp);
47 
48 static void pbex_fourth(FunFourthFuncDrv *ds, real factor,
49                            const FunDensProp* dp);
50 
51 //static integer fun_true(void) { return 1; }
52 Functional PBExFunctional = {
53   "PBEx",
54   pbex_isgga,
55   3,
56   pbex_read,
57   NULL,
58   pbex_energy,
59   pbex_first,
60   pbex_second,
61   pbex_third,
62   pbex_fourth
63 };
64 
65 /* IMPLEMENTATION PART */
66 static integer
pbex_read(const char * conf_line)67 pbex_read(const char* conf_line)
68 {
69     fun_set_hf_weight(0);
70     return 1;
71 }
72 
73 
74 static real
pbex_energy(const FunDensProp * dp)75 pbex_energy(const FunDensProp* dp)
76 {
77     real t[2],zk;
78     real rhoa = dp->rhoa;
79     real rhob = dp->rhob;
80     real grada = dp->grada;
81     real gradb = dp->gradb;
82     real gradab = dp->gradab;
83 
84     t[1] = pow(2.0,1.333333333333333);
85     zk = 0.5*(-0.73855876638202*t[1]*pow(rhob,1.333333333333333)*(1.804-0.804/(0.0044927994802311*pow(gradb,2.0)/pow(rhob,2.666666666666667)+1.0))-0.73855876638202*t[1]*pow(rhoa,1.333333333333333)*(1.804-0.804/(0.0044927994802311*pow(grada,2.0)/pow(rhoa,2.666666666666667)+1.0)));
86     return zk;
87 }
88 
89 static void
pbex_first(FunFirstFuncDrv * ds,real factor,const FunDensProp * dp)90 pbex_first(FunFirstFuncDrv *ds, real factor, const FunDensProp* dp)
91 {
92     real t[9];
93     real dfdra, dfdrb, dfdga, dfdgb, dfdab;
94     real rhoa = dp->rhoa;
95     real rhob = dp->rhob;
96     real grada = dp->grada;
97     real gradb = dp->gradb;
98     real gradab = dp->gradab;
99 
100     t[1] = pow(2.0,1.333333333333333);
101     t[2] = pow(grada,2.0);
102     t[3] = 0.0044927994802311*t[2]/pow(rhoa,2.666666666666667)+1.0;
103     t[4] = 1/pow(t[3],2.0);
104     t[5] = pow(2.0,3.333333333333334);
105     t[6] = pow(gradb,2.0);
106     t[7] = 0.0044927994802311*t[6]/pow(rhob,2.666666666666667)+1.0;
107     t[8] = 1/pow(t[7],2.0);
108     dfdra = 0.5*(0.0071142131710504*t[1]*t[2]*t[4]/pow(rhoa,2.333333333333334)-0.24618625546067*(1.804-0.804/t[3])*t[5]*pow(rhoa,0.33333333333333));
109     dfdrb = 0.5*(0.0071142131710504*t[1]*t[6]*t[8]/pow(rhob,2.333333333333334)-0.24618625546067*t[5]*(1.804-0.804/t[7])*pow(rhob,0.33333333333333));
110     dfdga = -0.0026678299391439*t[1]*t[4]*grada/pow(rhoa,1.333333333333333);
111     dfdgb = -0.0026678299391439*t[1]*t[8]*gradb/pow(rhob,1.333333333333333);
112     dfdab = 0.0;
113     ds->df1000 += factor*dfdra;
114     ds->df0100 += factor*dfdrb;
115     ds->df0010 += factor*dfdga;
116     ds->df0001 += factor*dfdgb;
117     ds->df00001 += factor*dfdab;
118 }
119 
120 static void
pbex_second(FunSecondFuncDrv * ds,real factor,const FunDensProp * dp)121 pbex_second(FunSecondFuncDrv *ds, real factor, const FunDensProp* dp)
122 {
123     real t[20];
124     real dfdra, dfdrb, dfdga, dfdgb, dfdab;
125     real d2fdraga, d2fdrara, d2fdrarb, d2fdragb, d2fdrbrb;
126     real d2fdrbgb, d2fdgaga, d2fdgbgb, d2fdrbga;
127     real d2fdraab, d2fdrbab;
128     real d2fdgaab, d2fdgbab, d2fdabab, d2fdgagb;
129     real rhoa = dp->rhoa;
130     real rhob = dp->rhob;
131     real grada = dp->grada;
132     real gradb = dp->gradb;
133     real gradab = dp->gradab;
134 
135     t[1] = pow(2.0,1.333333333333333);
136     t[2] = pow(grada,2.0);
137     t[3] = 0.0044927994802311*t[2]/pow(rhoa,2.666666666666667)+1.0;
138     t[4] = 1/pow(t[3],2.0);
139     t[5] = 1/pow(rhoa,2.333333333333334);
140     t[6] = pow(2.0,3.333333333333334);
141     t[7] = 1.804-0.804/t[3];
142     t[8] = pow(gradb,2.0);
143     t[9] = 0.0044927994802311*t[8]/pow(rhob,2.666666666666667)+1.0;
144     t[10] = 1/pow(t[9],2.0);
145     t[11] = 1/pow(rhob,2.333333333333334);
146     t[12] = 1.804-0.804/t[9];
147     t[13] = 1/pow(rhoa,1.333333333333333);
148     t[14] = 1/pow(rhob,1.333333333333333);
149     t[15] = 1/pow(t[3],3.0);
150     t[16] = 1/pow(rhoa,3.333333333333334);
151     t[17] = pow(2.0,2.333333333333334);
152     t[18] = 1/pow(t[9],3.0);
153     t[19] = 1/pow(rhob,3.333333333333334);
154     dfdra = 0.5*(0.0071142131710504*t[1]*t[2]*t[4]*t[5]-0.24618625546067*t[6]*t[7]*pow(rhoa,0.33333333333333));
155     dfdrb = 0.5*(0.0071142131710504*t[1]*t[8]*t[10]*t[11]-0.24618625546067*t[12]*t[6]*pow(rhob,0.33333333333333));
156     dfdga = -0.0026678299391439*t[1]*grada*t[4]*t[13];
157     dfdgb = -0.0026678299391439*t[1]*gradb*t[10]*t[14];
158     dfdab = 0.0;
159     d2fdrara = 0.5*(1.704679105981247E-4*t[1]*t[15]*pow(grada,4.0)/pow(rhoa,6.0)-0.082062085153558*t[6]*t[7]/pow(rhoa,0.66666666666667)+0.0023714043903501*t[6]*t[2]*t[4]*t[16]-0.016599830732451*t[1]*t[2]*t[4]*t[16]);
160     d2fdrarb = 0.0;
161     d2fdraga = 0.5*(-1.2785093294859353E-4*t[1]*t[15]*pow(grada,3.0)/pow(rhoa,5.0)-0.0017785532927626*t[6]*grada*t[4]*t[5]+0.0071142131710504*t[17]*grada*t[4]*t[5]);
162     d2fdragb = 0.0;
163     d2fdraab = 0.0;
164     d2fdrbrb = 0.5*(1.704679105981247E-4*t[1]*t[18]*pow(gradb,4.0)/pow(rhob,6.0)-0.082062085153558*t[12]*t[6]/pow(rhob,0.66666666666667)+0.0023714043903501*t[6]*t[8]*t[10]*t[19]-0.016599830732451*t[1]*t[8]*t[10]*t[19]);
165     d2fdrbga = 0.0;
166     d2fdrbgb = 0.5*(-1.2785093294859353E-4*t[1]*t[18]*pow(gradb,3.0)/pow(rhob,5.0)-0.0017785532927626*t[6]*gradb*t[10]*t[11]+0.0071142131710504*t[17]*gradb*t[10]*t[11]);
167     d2fdrbab = 0.0;
168     d2fdgaga = 4.7944099855722582E-5*t[1]*t[15]*t[2]/pow(rhoa,4.0)-0.0026678299391439*t[1]*t[4]*t[13];
169     d2fdgagb = 0.0;
170     d2fdgaab = 0.0;
171     d2fdgbgb = 4.7944099855722582E-5*t[1]*t[18]*t[8]/pow(rhob,4.0)-0.0026678299391439*t[1]*t[10]*t[14];
172     d2fdgbab = 0.0;
173     d2fdabab = 0.0;
174     ds->df1000 += factor*dfdra;
175     ds->df0100 += factor*dfdrb;
176     ds->df0010 += factor*dfdga;
177     ds->df0001 += factor*dfdgb;
178     ds->df00001 += factor*dfdab;
179     ds->df2000 += factor*d2fdrara;
180     ds->df1100 += factor*d2fdrarb;
181     ds->df1010 += factor*d2fdraga;
182     ds->df1001 += factor*d2fdragb;
183     ds->df10001 += factor*d2fdraab;
184     ds->df0200 += factor*d2fdrbrb;
185     ds->df0110 += factor*d2fdrbga;
186     ds->df0101 += factor*d2fdrbgb;
187     ds->df01001 += factor*d2fdrbab;
188     ds->df0020 += factor*d2fdgaga;
189     ds->df0011 += factor*d2fdgagb;
190     ds->df00101+= factor*d2fdgaab;
191     ds->df0002 += factor*d2fdgbgb;
192     ds->df00011+= factor*d2fdgbab;
193     ds->df00002+= factor*d2fdabab;
194 }
195 
196 static void
pbex_third(FunThirdFuncDrv * ds,real factor,const FunDensProp * dp)197 pbex_third(FunThirdFuncDrv *ds, real factor, const FunDensProp* dp)
198 {
199     real t[37];
200     real dfdra, dfdrb, dfdga, dfdgb, dfdab;
201     real d2fdrara , d2fdrarb , d2fdraga , d2fdragb ;
202     real d2fdraab , d2fdrbrb , d2fdrbga , d2fdrbgb ;
203     real d2fdrbab , d2fdgaga , d2fdgagb , d2fdgaab ;
204     real d2fdgbgb , d2fdgbab , d2fdabab ;
205     real d3fdrarara , d3fdrararb , d3fdraraga , d3fdraragb ;
206     real d3fdraraab , d3fdrarbrb , d3fdrarbga , d3fdrarbgb ;
207     real d3fdrarbab , d3fdragaga , d3fdragagb , d3fdragaab ;
208     real d3fdragbgb , d3fdragbab , d3fdraabab , d3fdrbrbrb ;
209     real d3fdrbrbga , d3fdrbrbgb , d3fdrbrbab , d3fdrbgaga ;
210     real d3fdrbgagb , d3fdrbgaab , d3fdrbgbgb , d3fdrbgbab ;
211     real d3fdrbabab , d3fdgagaga , d3fdgagagb , d3fdgagaab ;
212     real d3fdgagbgb , d3fdgagbab , d3fdgaabab , d3fdgbgbgb ;
213     real d3fdgbgbab , d3fdgbabab , d3fdababab ;
214     real rhoa = dp->rhoa;
215     real rhob = dp->rhob;
216     real grada = dp->grada;
217     real gradb = dp->gradb;
218     real gradab = dp->gradab;
219 
220     t[1] = pow(2.0,1.333333333333333);
221     t[2] = pow(grada,2.0);
222     t[3] = 0.0044927994802311*t[2]/pow(rhoa,2.666666666666667)+1.0;
223     t[4] = 1/pow(t[3],2.0);
224     t[5] = 1/pow(rhoa,2.333333333333334);
225     t[6] = pow(2.0,3.333333333333334);
226     t[7] = 1.804-0.804/t[3];
227     t[8] = pow(gradb,2.0);
228     t[9] = 0.0044927994802311*t[8]/pow(rhob,2.666666666666667)+1.0;
229     t[10] = 1/pow(t[9],2.0);
230     t[11] = 1/pow(rhob,2.333333333333334);
231     t[12] = 1.804-0.804/t[9];
232     t[13] = 1/pow(rhoa,1.333333333333333);
233     t[14] = 1/pow(rhob,1.333333333333333);
234     t[15] = pow(grada,4.0);
235     t[16] = 1/pow(t[3],3.0);
236     t[17] = 1/pow(rhoa,6.0);
237     t[18] = 1/pow(rhoa,3.333333333333334);
238     t[19] = pow(grada,3.0);
239     t[20] = 1/pow(rhoa,5.0);
240     t[21] = pow(2.0,2.333333333333334);
241     t[22] = pow(gradb,4.0);
242     t[23] = 1/pow(t[9],3.0);
243     t[24] = 1/pow(rhob,6.0);
244     t[25] = 1/pow(rhob,3.333333333333334);
245     t[26] = pow(gradb,3.0);
246     t[27] = 1/pow(rhob,5.0);
247     t[28] = 1/pow(rhoa,4.0);
248     t[29] = 1/pow(rhob,4.0);
249     t[30] = 1/pow(t[3],4.0);
250     t[31] = 1/pow(rhoa,7.0);
251     t[32] = 1/pow(rhoa,4.333333333333333);
252     t[33] = pow(2.0,4.333333333333333);
253     t[34] = 1/pow(t[9],4.0);
254     t[35] = 1/pow(rhob,7.0);
255     t[36] = 1/pow(rhob,4.333333333333333);
256     dfdra = 0.5*(0.0071142131710504*t[1]*t[2]*t[4]*t[5]-0.24618625546067*t[6]*t[7]*pow(rhoa,0.33333333333333));
257     dfdrb = 0.5*(0.0071142131710504*t[1]*t[8]*t[10]*t[11]-0.24618625546067*t[12]*t[6]*pow(rhob,0.33333333333333));
258     dfdga = -0.0026678299391439*t[1]*grada*t[4]*t[13];
259     dfdgb = -0.0026678299391439*t[1]*gradb*t[10]*t[14];
260     dfdab = 0.0;
261     d2fdrara = 0.5*(-0.082062085153558*t[6]*t[7]/pow(rhoa,0.66666666666667)+0.0023714043903501*t[6]*t[2]*t[4]*t[18]-0.016599830732451*t[1]*t[2]*t[4]*t[18]+1.704679105981247E-4*t[1]*t[15]*t[16]*t[17]);
262     d2fdrarb = 0.0;
263     d2fdraga = 0.5*(-0.0017785532927626*t[6]*grada*t[4]*t[5]+0.0071142131710504*t[21]*grada*t[4]*t[5]-1.2785093294859353E-4*t[1]*t[19]*t[16]*t[20]);
264     d2fdragb = 0.0;
265     d2fdraab = 0.0;
266     d2fdrbrb = 0.5*(-0.082062085153558*t[12]*t[6]/pow(rhob,0.66666666666667)+0.0023714043903501*t[6]*t[8]*t[10]*t[25]-0.016599830732451*t[1]*t[8]*t[10]*t[25]+1.704679105981247E-4*t[1]*t[22]*t[23]*t[24]);
267     d2fdrbga = 0.0;
268     d2fdrbgb = 0.5*(-0.0017785532927626*t[6]*gradb*t[10]*t[11]+0.0071142131710504*t[21]*gradb*t[10]*t[11]-1.2785093294859353E-4*t[1]*t[26]*t[23]*t[27]);
269     d2fdrbab = 0.0;
270     d2fdgaga = 4.7944099855722582E-5*t[1]*t[2]*t[16]*t[28]-0.0026678299391439*t[1]*t[4]*t[13];
271     d2fdgagb = 0.0;
272     d2fdgaab = 0.0;
273     d2fdgbgb = 4.7944099855722582E-5*t[1]*t[8]*t[23]*t[29]-0.0026678299391439*t[1]*t[10]*t[14];
274     d2fdgbab = 0.0;
275     d2fdabab = 0.0;
276     d3fdrarara = 0.5*(6.1270251210506771E-6*t[1]*t[30]*pow(grada,6.0)/pow(rhoa,9.666666666666666)+0.027354028384519*t[33]*t[7]/pow(rhoa,1.666666666666667)+7.9046813011671027E-4*t[6]*t[2]*t[4]*t[32]-0.0039523406505836*t[33]*t[2]*t[4]*t[32]+0.027666384554085*t[21]*t[2]*t[4]*t[32]+5.6822636866041565E-5*t[6]*t[15]*t[16]*t[31]-5.1140373179437413E-4*t[21]*t[15]*t[16]*t[31]-3.9775845806229101E-4*t[1]*t[15]*t[16]*t[31]);
277     d3fdrararb = 0.0;
278     d3fdraraga = 0.5*(-4.5952688407880086E-6*t[1]*t[30]*pow(grada,5.0)/pow(rhoa,8.666666666666666)-5.9285109758753281E-4*t[6]*grada*t[4]*t[18]+0.0023714043903501*t[33]*grada*t[4]*t[18]-0.016599830732451*t[21]*grada*t[4]*t[18]+1.2785093294859353E-4*t[6]*t[19]*t[16]*t[17]+2.9831884354671826E-4*t[1]*t[19]*t[16]*t[17]);
279     d3fdraragb = 0.0;
280     d3fdraraab = 0.0;
281     d3fdrarbrb = 0.0;
282     d3fdrarbga = 0.0;
283     d3fdrarbgb = 0.0;
284     d3fdrarbab = 0.0;
285     d3fdragaga = 0.5*(3.4464516305910063E-6*t[1]*t[15]*t[30]/pow(rhoa,7.666666666666667)-0.0017785532927626*t[6]*t[4]*t[5]+0.0071142131710504*t[21]*t[4]*t[5]+3.1962733237148383E-5*t[6]*t[2]*t[16]*t[20]-1.2785093294859353E-4*t[21]*t[2]*t[16]*t[20]-3.835527988457806E-4*t[1]*t[2]*t[16]*t[20]);
286     d3fdragagb = 0.0;
287     d3fdragaab = 0.0;
288     d3fdragbgb = 0.0;
289     d3fdragbab = 0.0;
290     d3fdraabab = 0.0;
291     d3fdrbrbrb = 0.5*(6.1270251210506771E-6*t[1]*t[34]*pow(gradb,6.0)/pow(rhob,9.666666666666666)+0.027354028384519*t[12]*t[33]/pow(rhob,1.666666666666667)+7.9046813011671027E-4*t[6]*t[8]*t[10]*t[36]-0.0039523406505836*t[33]*t[8]*t[10]*t[36]+0.027666384554085*t[21]*t[8]*t[10]*t[36]+5.6822636866041565E-5*t[6]*t[22]*t[23]*t[35]-5.1140373179437413E-4*t[21]*t[22]*t[23]*t[35]-3.9775845806229101E-4*t[1]*t[22]*t[23]*t[35]);
292     d3fdrbrbga = 0.0;
293     d3fdrbrbgb = 0.5*(-4.5952688407880086E-6*t[1]*t[34]*pow(gradb,5.0)/pow(rhob,8.666666666666666)-5.9285109758753281E-4*t[6]*gradb*t[10]*t[25]+0.0023714043903501*t[33]*gradb*t[10]*t[25]-0.016599830732451*t[21]*gradb*t[10]*t[25]+1.2785093294859353E-4*t[6]*t[26]*t[23]*t[24]+2.9831884354671826E-4*t[1]*t[26]*t[23]*t[24]);
294     d3fdrbrbab = 0.0;
295     d3fdrbgaga = 0.0;
296     d3fdrbgagb = 0.0;
297     d3fdrbgaab = 0.0;
298     d3fdrbgbgb = 0.5*(3.4464516305910063E-6*t[1]*t[22]*t[34]/pow(rhob,7.666666666666667)+3.1962733237148383E-5*t[6]*t[8]*t[23]*t[27]-1.2785093294859353E-4*t[21]*t[8]*t[23]*t[27]-3.835527988457806E-4*t[1]*t[8]*t[23]*t[27]-0.0017785532927626*t[6]*t[10]*t[11]+0.0071142131710504*t[21]*t[10]*t[11]);
299     d3fdrbgbab = 0.0;
300     d3fdrbabab = 0.0;
301     d3fdgagaga = -1.2924193614716277E-6*t[1]*t[19]*t[30]/pow(rhoa,6.666666666666667)+4.7944099855722582E-5*t[21]*grada*t[16]*t[28]+4.7944099855722582E-5*t[1]*grada*t[16]*t[28];
302     d3fdgagagb = 0.0;
303     d3fdgagaab = 0.0;
304     d3fdgagbgb = 0.0;
305     d3fdgagbab = 0.0;
306     d3fdgaabab = 0.0;
307     d3fdgbgbgb = -1.2924193614716277E-6*t[1]*t[26]*t[34]/pow(rhob,6.666666666666667)+4.7944099855722582E-5*t[21]*gradb*t[23]*t[29]+4.7944099855722582E-5*t[1]*gradb*t[23]*t[29];
308     d3fdgbgbab = 0.0;
309     d3fdgbabab = 0.0;
310     d3fdababab = 0.0;
311     ds->df1000 += factor*dfdra;
312     ds->df0100 += factor*dfdrb;
313     ds->df0010 += factor*dfdga;
314     ds->df0001 += factor*dfdgb;
315     ds->df00001 += factor*dfdab;
316     ds->df2000 += factor*d2fdrara;
317     ds->df1100 += factor*d2fdrarb;
318     ds->df1010 += factor*d2fdraga;
319     ds->df1001 += factor*d2fdragb;
320     ds->df10001 += factor*d2fdraab;
321     ds->df0200 += factor*d2fdrbrb;
322     ds->df0110 += factor*d2fdrbga;
323     ds->df0101 += factor*d2fdrbgb;
324     ds->df01001 += factor*d2fdrbab;
325     ds->df0020 += factor*d2fdgaga;
326     ds->df0011 += factor*d2fdgagb;
327     ds->df00101+= factor*d2fdgaab;
328     ds->df0002 += factor*d2fdgbgb;
329     ds->df00011+= factor*d2fdgbab;
330     ds->df00002+= factor*d2fdabab;
331     ds->df3000 += factor*d3fdrarara;
332     ds->df2100 += factor*d3fdrararb;
333     ds->df2010 += factor*d3fdraraga;
334     ds->df2001 += factor*d3fdraragb;
335     ds->df20001 += factor*d3fdraraab;
336     ds->df1200 += factor*d3fdrarbrb;
337     ds->df1110 += factor*d3fdrarbga;
338     ds->df1101 += factor*d3fdrarbgb;
339     ds->df11001 += factor*d3fdrarbab;
340     ds->df1020 += factor*d3fdragaga;
341     ds->df1011 += factor*d3fdragagb;
342     ds->df10101+= factor*d3fdragaab;
343     ds->df1002 += factor*d3fdragbgb;
344     ds->df10011+= factor*d3fdragbab;
345     ds->df10002+= factor*d3fdraabab;
346     ds->df0300 += factor*d3fdrbrbrb;
347     ds->df0210 += factor*d3fdrbrbga;
348     ds->df0201 += factor*d3fdrbrbgb;
349     ds->df02001 += factor*d3fdrbrbab;
350     ds->df0120 += factor*d3fdrbgaga;
351     ds->df0111 += factor*d3fdrbgagb;
352     ds->df01101+= factor*d3fdrbgaab;
353     ds->df0102 += factor*d3fdrbgbgb;
354     ds->df01011+= factor*d3fdrbgbab;
355     ds->df01002+= factor*d3fdrbgbab;
356     ds->df0030 += factor*d3fdgagaga;
357     ds->df0021 += factor*d3fdgagagb;
358     ds->df00201+= factor*d3fdgagaab;
359     ds->df0012 += factor*d3fdgagbgb;
360     ds->df00111+= factor*d3fdgagbab;
361     ds->df00102+= factor*d3fdgaabab;
362     ds->df0003 += factor*d3fdgbgbgb;
363     ds->df00021+= factor*d3fdgbgbab;
364     ds->df00012+= factor*d3fdgbabab;
365     ds->df00003+= factor*d3fdababab;
366 }
367 
368 static void
pbex_fourth(FunFourthFuncDrv * ds,real factor,const FunDensProp * dp)369 pbex_fourth(FunFourthFuncDrv *ds, real factor, const FunDensProp* dp)
370 {
371     real t[60];
372     real dfdra, dfdrb, dfdga, dfdgb, dfdab;
373     real d2fdrara , d2fdrarb , d2fdraga , d2fdragb ;
374     real d2fdraab , d2fdrbrb , d2fdrbga , d2fdrbgb ;
375     real d2fdrbab , d2fdgaga , d2fdgagb , d2fdgaab ;
376     real d2fdgbgb , d2fdgbab , d2fdabab ;
377     real d3fdrarara , d3fdrararb , d3fdraraga , d3fdraragb ;
378     real d3fdraraab , d3fdrarbrb , d3fdrarbga , d3fdrarbgb ;
379     real d3fdrarbab , d3fdragaga , d3fdragagb , d3fdragaab ;
380     real d3fdragbgb , d3fdragbab , d3fdraabab , d3fdrbrbrb ;
381     real d3fdrbrbga , d3fdrbrbgb , d3fdrbrbab , d3fdrbgaga ;
382     real d3fdrbgagb , d3fdrbgaab , d3fdrbgbgb , d3fdrbgbab ;
383     real d3fdrbabab , d3fdgagaga , d3fdgagagb , d3fdgagaab ;
384     real d3fdgagbgb , d3fdgagbab , d3fdgaabab , d3fdgbgbgb ;
385     real d3fdgbgbab , d3fdgbabab , d3fdababab ;
386     real d4fdrararara , d4fdrarararb , d4fdrararaga , d4fdrararagb ;
387     real d4fdrararaab , d4fdrararbrb , d4fdrararbga , d4fdrararbgb ;
388     real d4fdrararbab , d4fdraragaga , d4fdraragagb , d4fdraragaab ;
389     real d4fdraragbgb , d4fdraragbab , d4fdraraabab , d4fdrarbrbrb ;
390     real d4fdrarbrbga , d4fdrarbrbgb , d4fdrarbrbab , d4fdrarbgaga ;
391     real d4fdrarbgagb , d4fdrarbgaab , d4fdrarbgbgb , d4fdrarbgbab ;
392     real d4fdrarbabab , d4fdragagaga , d4fdragagagb , d4fdragagaab ;
393     real d4fdragagbgb , d4fdragagbab , d4fdragaabab , d4fdragbgbgb ;
394     real d4fdragbgbab , d4fdragbabab , d4fdraababab , d4fdrbrbrbrb ;
395     real d4fdrbrbrbga , d4fdrbrbrbgb , d4fdrbrbrbab , d4fdrbrbgaga ;
396     real d4fdrbrbgagb , d4fdrbrbgaab , d4fdrbrbgbgb , d4fdrbrbgbab ;
397     real d4fdrbrbabab , d4fdrbgagaga , d4fdrbgagagb , d4fdrbgagaab ;
398     real d4fdrbgagbgb , d4fdrbgagbab , d4fdrbgaabab , d4fdrbgbgbgb ;
399     real d4fdrbgbgbab , d4fdrbgbabab , d4fdrbababab , d4fdgagagaga ;
400     real d4fdgagagagb , d4fdgagagaab , d4fdgagagbgb , d4fdgagagbab ;
401     real d4fdgagaabab , d4fdgagbgbgb , d4fdgagbgbab , d4fdgagbabab ;
402     real d4fdgaababab , d4fdgbgbgbgb , d4fdgbgbgbab , d4fdgbgbabab ;
403     real d4fdgbababab , d4fdabababab ;
404     real rhoa = dp->rhoa;
405     real rhob = dp->rhob;
406     real grada = dp->grada;
407     real gradb = dp->gradb;
408     real gradab = dp->gradab;
409 
410     t[1] = pow(2.0,1.333333333333333);
411     t[2] = pow(grada,2.0);
412     t[3] = 1/pow(rhoa,2.666666666666667);
413     t[4] = 0.0044927994802311*t[2]*t[3]+1.0;
414     t[5] = 1/pow(t[4],2.0);
415     t[6] = 1/pow(rhoa,2.333333333333334);
416     t[7] = pow(2.0,3.333333333333334);
417     t[8] = 1.804-0.804/t[4];
418     t[9] = pow(gradb,2.0);
419     t[10] = 1/pow(rhob,2.666666666666667);
420     t[11] = 0.0044927994802311*t[9]*t[10]+1.0;
421     t[12] = 1/pow(t[11],2.0);
422     t[13] = 1/pow(rhob,2.333333333333334);
423     t[14] = 1.804-0.804/t[11];
424     t[15] = 1/pow(rhoa,1.333333333333333);
425     t[16] = 1/pow(rhob,1.333333333333333);
426     t[17] = pow(grada,4.0);
427     t[18] = 1/pow(t[4],3.0);
428     t[19] = 1/pow(rhoa,6.0);
429     t[20] = 1/pow(rhoa,3.333333333333334);
430     t[21] = pow(grada,3.0);
431     t[22] = 1/pow(rhoa,5.0);
432     t[23] = pow(2.0,2.333333333333334);
433     t[24] = pow(gradb,4.0);
434     t[25] = 1/pow(t[11],3.0);
435     t[26] = 1/pow(rhob,6.0);
436     t[27] = 1/pow(rhob,3.333333333333334);
437     t[28] = pow(gradb,3.0);
438     t[29] = 1/pow(rhob,5.0);
439     t[30] = 1/pow(rhoa,4.0);
440     t[31] = 1/pow(rhob,4.0);
441     t[32] = pow(grada,6.0);
442     t[33] = 1/pow(t[4],4.0);
443     t[34] = 1/pow(rhoa,9.666666666666666);
444     t[35] = 1/pow(rhoa,7.0);
445     t[36] = 1/pow(rhoa,4.333333333333333);
446     t[37] = pow(2.0,4.333333333333333);
447     t[38] = pow(grada,5.0);
448     t[39] = 1/pow(rhoa,8.666666666666666);
449     t[40] = 1/pow(rhoa,7.666666666666667);
450     t[41] = pow(gradb,6.0);
451     t[42] = 1/pow(t[11],4.0);
452     t[43] = 1/pow(rhob,9.666666666666666);
453     t[44] = 1/pow(rhob,7.0);
454     t[45] = 1/pow(rhob,4.333333333333333);
455     t[46] = pow(gradb,5.0);
456     t[47] = 1/pow(rhob,8.666666666666666);
457     t[48] = 1/pow(rhob,7.666666666666667);
458     t[49] = 1/pow(rhoa,6.666666666666667);
459     t[50] = 1/pow(rhob,6.666666666666667);
460     t[51] = 1/pow(t[4],5.0);
461     t[52] = 1/pow(rhoa,10.66666666666667);
462     t[53] = 1/pow(rhoa,8.0);
463     t[54] = 1/pow(rhoa,5.333333333333333);
464     t[55] = pow(2.0,5.333333333333333);
465     t[56] = 1/pow(t[11],5.0);
466     t[57] = 1/pow(rhob,10.66666666666667);
467     t[58] = 1/pow(rhob,8.0);
468     t[59] = 1/pow(rhob,5.333333333333333);
469     dfdra = 0.5*(0.0071142131710504*t[1]*t[2]*t[5]*t[6]-0.24618625546067*t[7]*t[8]*pow(rhoa,0.33333333333333));
470     dfdrb = 0.5*(0.0071142131710504*t[1]*t[9]*t[12]*t[13]-0.24618625546067*t[14]*t[7]*pow(rhob,0.33333333333333));
471     dfdga = -0.0026678299391439*t[1]*grada*t[5]*t[15];
472     dfdgb = -0.0026678299391439*t[1]*gradb*t[12]*t[16];
473     dfdab = 0.0;
474     d2fdrara = 0.5*(-0.082062085153558*t[7]*t[8]/pow(rhoa,0.66666666666667)+0.0023714043903501*t[7]*t[2]*t[5]*t[20]-0.016599830732451*t[1]*t[2]*t[5]*t[20]+1.704679105981247E-4*t[1]*t[17]*t[18]*t[19]);
475     d2fdrarb = 0.0;
476     d2fdraga = 0.5*(-0.0017785532927626*t[7]*grada*t[5]*t[6]+0.0071142131710504*t[23]*grada*t[5]*t[6]-1.2785093294859353E-4*t[1]*t[21]*t[18]*t[22]);
477     d2fdragb = 0.0;
478     d2fdraab = 0.0;
479     d2fdrbrb = 0.5*(-0.082062085153558*t[14]*t[7]/pow(rhob,0.66666666666667)+0.0023714043903501*t[7]*t[9]*t[12]*t[27]-0.016599830732451*t[1]*t[9]*t[12]*t[27]+1.704679105981247E-4*t[1]*t[24]*t[25]*t[26]);
480     d2fdrbga = 0.0;
481     d2fdrbgb = 0.5*(-0.0017785532927626*t[7]*gradb*t[12]*t[13]+0.0071142131710504*t[23]*gradb*t[12]*t[13]-1.2785093294859353E-4*t[1]*t[28]*t[25]*t[29]);
482     d2fdrbab = 0.0;
483     d2fdgaga = 4.7944099855722582E-5*t[1]*t[2]*t[18]*t[30]-0.0026678299391439*t[1]*t[5]*t[15];
484     d2fdgagb = 0.0;
485     d2fdgaab = 0.0;
486     d2fdgbgb = 4.7944099855722582E-5*t[1]*t[9]*t[25]*t[31]-0.0026678299391439*t[1]*t[12]*t[16];
487     d2fdgbab = 0.0;
488     d2fdabab = 0.0;
489     d3fdrarara= 0.5*(0.027354028384519*t[37]*t[8]/pow(rhoa,1.666666666666667)+7.9046813011671027E-4*t[7]*t[2]*t[5]*t[36]-0.0039523406505836*t[37]*t[2]*t[5]*t[36]+0.027666384554085*t[23]*t[2]*t[5]*t[36]+5.6822636866041565E-5*t[7]*t[17]*t[18]*t[35]-5.1140373179437413E-4*t[23]*t[17]*t[18]*t[35]-3.9775845806229101E-4*t[1]*t[17]*t[18]*t[35]+6.1270251210506771E-6*t[1]*t[32]*t[33]*t[34]);
490     d3fdrararb = 0.0;
491     d3fdraraga = 0.5*(0.0023714043903501*t[37]*grada*t[5]*t[20]-5.9285109758753281E-4*t[7]*grada*t[5]*t[20]-0.016599830732451*t[23]*grada*t[5]*t[20]+1.2785093294859353E-4*t[7]*t[21]*t[18]*t[19]+2.9831884354671826E-4*t[1]*t[21]*t[18]*t[19]-4.5952688407880086E-6*t[1]*t[38]*t[33]*t[39]);
492     d3fdraragb = 0.0;
493     d3fdraraab = 0.0;
494     d3fdrarbrb = 0.0;
495     d3fdrarbga = 0.0;
496     d3fdrarbgb = 0.0;
497     d3fdrarbab = 0.0;
498     d3fdragaga = 0.5*(-0.0017785532927626*t[7]*t[5]*t[6]+0.0071142131710504*t[23]*t[5]*t[6]+3.1962733237148383E-5*t[7]*t[2]*t[18]*t[22]-1.2785093294859353E-4*t[23]*t[2]*t[18]*t[22]-3.835527988457806E-4*t[1]*t[2]*t[18]*t[22]+3.4464516305910063E-6*t[1]*t[17]*t[33]*t[40]);
499     d3fdragagb = 0.0;
500     d3fdragaab = 0.0;
501     d3fdragbgb = 0.0;
502     d3fdragbab = 0.0;
503     d3fdraabab = 0.0;
504     d3fdrbrbrb = 0.5*(0.027354028384519*t[14]*t[37]/pow(rhob,1.666666666666667)+7.9046813011671027E-4*t[7]*t[9]*t[12]*t[45]-0.0039523406505836*t[37]*t[9]*t[12]*t[45]+0.027666384554085*t[23]*t[9]*t[12]*t[45]+5.6822636866041565E-5*t[7]*t[24]*t[25]*t[44]-5.1140373179437413E-4*t[23]*t[24]*t[25]*t[44]-3.9775845806229101E-4*t[1]*t[24]*t[25]*t[44]+6.1270251210506771E-6*t[1]*t[41]*t[42]*t[43]);
505     d3fdrbrbga = 0.0;
506     d3fdrbrbgb = 0.5*(0.0023714043903501*t[37]*gradb*t[12]*t[27]-5.9285109758753281E-4*t[7]*gradb*t[12]*t[27]-0.016599830732451*t[23]*gradb*t[12]*t[27]+1.2785093294859353E-4*t[7]*t[28]*t[25]*t[26]+2.9831884354671826E-4*t[1]*t[28]*t[25]*t[26]-4.5952688407880086E-6*t[1]*t[46]*t[42]*t[47]);
507     d3fdrbrbab = 0.0;
508     d3fdrbgaga = 0.0;
509     d3fdrbgagb = 0.0;
510     d3fdrbgaab = 0.0;
511     d3fdrbgbgb= 0.5*(-0.0017785532927626*t[7]*t[12]*t[13]+0.0071142131710504*t[23]*t[12]*t[13]+3.1962733237148383E-5*t[7]*t[9]*t[25]*t[29]-1.2785093294859353E-4*t[23]*t[9]*t[25]*t[29]-3.835527988457806E-4*t[1]*t[9]*t[25]*t[29]+3.4464516305910063E-6*t[1]*t[24]*t[42]*t[48]);
512     d3fdrbgbab = 0.0;
513     d3fdrbabab = 0.0;
514     d3fdgagaga = 4.7944099855722582E-5*t[23]*grada*t[18]*t[30]+4.7944099855722582E-5*t[1]*grada*t[18]*t[30]-1.2924193614716277E-6*t[1]*t[21]*t[33]*t[49];
515     d3fdgagagb = 0.0;
516     d3fdgagaab = 0.0;
517     d3fdgagbgb = 0.0;
518     d3fdgagbab = 0.0;
519     d3fdgaabab = 0.0;
520     d3fdgbgbgb = 4.7944099855722582E-5*t[23]*gradb*t[25]*t[31]+4.7944099855722582E-5*t[1]*gradb*t[25]*t[31]-1.2924193614716277E-6*t[1]*t[28]*t[42]*t[50];
521     d3fdgbgbab = 0.0;
522     d3fdgbabab = 0.0;
523     d3fdababab = 0.0;
524     d4fdrararara = 0.5*(2.9362661631167268E-7*t[1]*t[51]*pow(grada,8.0)/pow(rhoa,13.33333333333333)-0.0034253618971724*t[7]*t[2]*t[5]*t[54]+0.016863320109156*t[37]*t[2]*t[5]*t[54]-0.11988766640103*t[23]*t[2]*t[5]*t[54]-3.7881757910694376E-4*t[7]*t[17]*t[18]*t[53]-9.4704394776735953E-5*t[37]*t[17]*t[18]*t[53]+0.0042427568859978*t[23]*t[17]*t[18]*t[53]+0.002784309206436*t[1]*t[17]*t[18]*t[53]+2.0423417070168925E-6*t[7]*t[32]*t[33]*t[52]-1.8381075363152035E-5*t[23]*t[32]*t[33]*t[52]-7.3524301452608125E-5*t[1]*t[32]*t[33]*t[52]-0.045590047307532*t[37]*t[8]*t[3]);
525     d4fdrarararb = 0.0;
526     d4fdrararaga = 0.5*(-2.2021996223375453E-7*t[1]*t[51]*pow(grada,7.0)/pow(rhoa,12.33333333333333)+0.027666384554085*t[7]*grada*t[5]*t[36]-0.0039523406505836*t[55]*grada*t[5]*t[36]+9.8808516264588774E-4*t[37]*grada*t[5]*t[36]-4.1196411727880143E-4*t[7]*t[21]*t[18]*t[35]+5.6822636866041565E-5*t[55]*t[21]*t[18]*t[35]-4.4037543571182216E-4*t[37]*t[21]*t[18]*t[35]-4.9719807257786388E-4*t[23]*t[21]*t[18]*t[35]-1.5317562802626695E-6*t[7]*t[38]*t[33]*t[34]+3.216688188551606E-5*t[23]*t[38]*t[33]*t[34]+1.0722293961838687E-5*t[1]*t[38]*t[33]*t[34]);
527     d4fdrararagb = 0.0;
528     d4fdrararaab = 0.0;
529     d4fdrararbrb = 0.0;
530     d4fdrararbga = 0.0;
531     d4fdrararbgb = 0.0;
532     d4fdrararbab = 0.0;
533     d4fdraragaga = 0.5*(1.6516497167531593E-7*t[1]*t[32]*t[51]/pow(rhoa,11.33333333333333)-3.4464516305910063E-6*t[7]*t[17]*t[33]*t[39]-3.1018064675319062E-5*t[1]*t[17]*t[33]*t[39]-5.9285109758753281E-4*t[7]*t[5]*t[20]+0.0023714043903501*t[37]*t[5]*t[20]-0.016599830732451*t[23]*t[5]*t[20]+3.9420704325816337E-4*t[7]*t[2]*t[18]*t[19]-4.2616977649531175E-5*t[37]*t[2]*t[18]*t[19]+2.9831884354671826E-4*t[23]*t[2]*t[18]*t[19]+8.9495653064015478E-4*t[1]*t[2]*t[18]*t[19]);
534     d4fdraragagb = 0.0;
535     d4fdraragaab = 0.0;
536     d4fdraragbgb = 0.0;
537     d4fdraragbab = 0.0;
538     d4fdraraabab = 0.0;
539     d4fdrarbrbrb = 0.0;
540     d4fdrarbrbga = 0.0;
541     d4fdrarbrbgb = 0.0;
542     d4fdrarbrbab = 0.0;
543     d4fdrarbgaga = 0.0;
544     d4fdrarbgagb = 0.0;
545     d4fdrarbgaab = 0.0;
546     d4fdrarbgbgb = 0.0;
547     d4fdrarbgbab = 0.0;
548     d4fdrarbabab = 0.0;
549     d4fdragagaga = 0.5*(-1.2387372875648693E-7*t[1]*t[38]*t[51]/pow(rhoa,10.33333333333333)+2.5848387229432549E-6*t[7]*t[21]*t[33]*t[40]+3.4464516305910063E-6*t[23]*t[21]*t[33]*t[40]+1.033935489177302E-5*t[1]*t[21]*t[33]*t[40]-9.588819971144515E-5*t[7]*grada*t[18]*t[22]+3.1962733237148383E-5*t[37]*grada*t[18]*t[22]-5.1140373179437413E-4*t[23]*grada*t[18]*t[22]);
550     d4fdragagagb = 0.0;
551     d4fdragagaab = 0.0;
552     d4fdragagbgb = 0.0;
553     d4fdragagbab = 0.0;
554     d4fdragaabab = 0.0;
555     d4fdragbgbgb = 0.0;
556     d4fdragbgbab = 0.0;
557     d4fdragbabab = 0.0;
558     d4fdraababab = 0.0;
559     d4fdrbrbrbrb= 0.5*(2.9362661631167268E-7*t[1]*t[56]*pow(gradb,8.0)/pow(rhob,13.33333333333333)-0.0034253618971724*t[7]*t[9]*t[12]*t[59]+0.016863320109156*t[37]*t[9]*t[12]*t[59]-0.11988766640103*t[23]*t[9]*t[12]*t[59]-3.7881757910694376E-4*t[7]*t[24]*t[25]*t[58]-9.4704394776735953E-5*t[37]*t[24]*t[25]*t[58]+0.0042427568859978*t[23]*t[24]*t[25]*t[58]+0.002784309206436*t[1]*t[24]*t[25]*t[58]+2.0423417070168925E-6*t[7]*t[41]*t[42]*t[57]-1.8381075363152035E-5*t[23]*t[41]*t[42]*t[57]-7.3524301452608125E-5*t[1]*t[41]*t[42]*t[57]-0.045590047307532*t[37]*t[14]*t[10]);
560     d4fdrbrbrbga = 0.0;
561     d4fdrbrbrbgb = 0.5*(-2.2021996223375453E-7*t[1]*t[56]*pow(gradb,7.0)/pow(rhob,12.33333333333333)+0.027666384554085*t[7]*gradb*t[12]*t[45]-0.0039523406505836*t[55]*gradb*t[12]*t[45]+9.8808516264588774E-4*t[37]*gradb*t[12]*t[45]-4.1196411727880143E-4*t[7]*t[28]*t[25]*t[44]+5.6822636866041565E-5*t[55]*t[28]*t[25]*t[44]-4.4037543571182216E-4*t[37]*t[28]*t[25]*t[44]-4.9719807257786388E-4*t[23]*t[28]*t[25]*t[44]-1.5317562802626695E-6*t[7]*t[46]*t[42]*t[43]+3.216688188551606E-5*t[23]*t[46]*t[42]*t[43]+1.0722293961838687E-5*t[1]*t[46]*t[42]*t[43]);
562     d4fdrbrbrbab = 0.0;
563     d4fdrbrbgaga = 0.0;
564     d4fdrbrbgagb = 0.0;
565     d4fdrbrbgaab = 0.0;
566     d4fdrbrbgbgb = 0.5*(1.6516497167531593E-7*t[1]*t[41]*t[56]/pow(rhob,11.33333333333333)-3.4464516305910063E-6*t[7]*t[24]*t[42]*t[47]-3.1018064675319062E-5*t[1]*t[24]*t[42]*t[47]-5.9285109758753281E-4*t[7]*t[12]*t[27]+0.0023714043903501*t[37]*t[12]*t[27]-0.016599830732451*t[23]*t[12]*t[27]+3.9420704325816337E-4*t[7]*t[9]*t[25]*t[26]-4.2616977649531175E-5*t[37]*t[9]*t[25]*t[26]+2.9831884354671826E-4*t[23]*t[9]*t[25]*t[26]+8.9495653064015478E-4*t[1]*t[9]*t[25]*t[26]);
567     d4fdrbrbgbab = 0.0;
568     d4fdrbrbabab = 0.0;
569     d4fdrbgagaga = 0.0;
570     d4fdrbgagagb = 0.0;
571     d4fdrbgagaab = 0.0;
572     d4fdrbgagbgb = 0.0;
573     d4fdrbgagbab = 0.0;
574     d4fdrbgaabab = 0.0;
575     d4fdrbgbgbgb = 0.5*(-1.2387372875648693E-7*t[1]*t[46]*t[56]/pow(rhob,10.33333333333333)+2.5848387229432549E-6*t[7]*t[28]*t[42]*t[48]+3.4464516305910063E-6*t[23]*t[28]*t[42]*t[48]+1.033935489177302E-5*t[1]*t[28]*t[42]*t[48]-9.588819971144515E-5*t[7]*gradb*t[25]*t[29]+3.1962733237148383E-5*t[37]*gradb*t[25]*t[29]-5.1140373179437413E-4*t[23]*gradb*t[25]*t[29]);
576     d4fdrbgbgbab = 0.0;
577     d4fdrbgbabab = 0.0;
578     d4fdrbababab = 0.0;
579     d4fdgagagaga = 4.6452648283682616E-8*t[1]*t[17]*t[51]/pow(rhoa,9.333333333333334)-1.2924193614716277E-6*t[23]*t[2]*t[33]*t[49]-5.1696774458865107E-6*t[1]*t[2]*t[33]*t[49]+4.7944099855722582E-5*t[23]*t[18]*t[30]+4.7944099855722582E-5*t[1]*t[18]*t[30];
580     d4fdgagagagb = 0.0;
581     d4fdgagagaab = 0.0;
582     d4fdgagagbgb = 0.0;
583     d4fdgagagbab = 0.0;
584     d4fdgagaabab = 0.0;
585     d4fdgagbgbgb = 0.0;
586     d4fdgagbgbab = 0.0;
587     d4fdgagbabab = 0.0;
588     d4fdgaababab = 0.0;
589     d4fdgbgbgbgb = 4.6452648283682616E-8*t[1]*t[24]*t[56]/pow(rhob,9.333333333333334)-1.2924193614716277E-6*t[23]*t[9]*t[42]*t[50]-5.1696774458865107E-6*t[1]*t[9]*t[42]*t[50]+4.7944099855722582E-5*t[23]*t[25]*t[31]+4.7944099855722582E-5*t[1]*t[25]*t[31];
590     d4fdgbgbgbab = 0.0;
591     d4fdgbgbabab = 0.0;
592     d4fdgbababab = 0.0;
593     d4fdabababab = 0.0;
594     ds->df1000 += factor*dfdra;
595     ds->df0100 += factor*dfdrb;
596     ds->df0010 += factor*dfdga;
597     ds->df0001 += factor*dfdgb;
598     ds->df00001 += factor*dfdab;
599     ds->df2000 += factor*d2fdrara;
600     ds->df1100 += factor*d2fdrarb;
601     ds->df1010 += factor*d2fdraga;
602     ds->df1001 += factor*d2fdragb;
603     ds->df10001 += factor*d2fdraab;
604     ds->df0200 += factor*d2fdrbrb;
605     ds->df0110 += factor*d2fdrbga;
606     ds->df0101 += factor*d2fdrbgb;
607     ds->df01001 += factor*d2fdrbab;
608     ds->df0020 += factor*d2fdgaga;
609     ds->df0011 += factor*d2fdgagb;
610     ds->df00101 += factor*d2fdgaab;
611     ds->df0002 += factor*d2fdgbgb;
612     ds->df00011 += factor*d2fdgbab;
613     ds->df00002 += factor*d2fdabab;
614     ds->df3000 += factor*d3fdrarara;
615     ds->df2100 += factor*d3fdrararb;
616     ds->df2010 += factor*d3fdraraga;
617     ds->df2001 += factor*d3fdraragb;
618     ds->df20001 += factor*d3fdraraab;
619     ds->df1200 += factor*d3fdrarbrb;
620     ds->df1110 += factor*d3fdrarbga;
621     ds->df1101 += factor*d3fdrarbgb;
622     ds->df11001 += factor*d3fdrarbab;
623     ds->df1020 += factor*d3fdragaga;
624     ds->df1011 += factor*d3fdragagb;
625     ds->df10101 += factor*d3fdragaab;
626     ds->df1002 += factor*d3fdragbgb;
627     ds->df10011 += factor*d3fdragbab;
628     ds->df10002 += factor*d3fdraabab;
629     ds->df0300 += factor*d3fdrbrbrb;
630     ds->df0210 += factor*d3fdrbrbga;
631     ds->df0201 += factor*d3fdrbrbgb;
632     ds->df02001 += factor*d3fdrbrbab;
633     ds->df0120 += factor*d3fdrbgaga;
634     ds->df0111 += factor*d3fdrbgagb;
635     ds->df01101 += factor*d3fdrbgaab;
636     ds->df0102 += factor*d3fdrbgbgb;
637     ds->df01011 += factor*d3fdrbgbab;
638     ds->df01002 += factor*d3fdrbabab;
639     ds->df0030 += factor*d3fdgagaga;
640     ds->df0021 += factor*d3fdgagagb;
641     ds->df00201 += factor*d3fdgagaab;
642     ds->df0012 += factor*d3fdgagbgb;
643     ds->df00111 += factor*d3fdgagbab;
644     ds->df00102 += factor*d3fdgaabab;
645     ds->df0003 += factor*d3fdgbgbgb;
646     ds->df00021 += factor*d3fdgbgbab;
647     ds->df00012 += factor*d3fdgbabab;
648     ds->df00003 += factor*d3fdababab;
649     ds->df4000 += factor*d4fdrararara;
650     ds->df3100 += factor*d4fdrarararb;
651     ds->df3010 += factor*d4fdrararaga;
652     ds->df3001 += factor*d4fdrararagb;
653     ds->df30001 += factor*d4fdrararaab;
654     ds->df2200 += factor*d4fdrararbrb;
655     ds->df2110 += factor*d4fdrararbga;
656     ds->df2101 += factor*d4fdrararbgb;
657     ds->df21001 += factor*d4fdrararbab;
658     ds->df2020 += factor*d4fdraragaga;
659     ds->df2011 += factor*d4fdraragagb;
660     ds->df20101 += factor*d4fdraragaab;
661     ds->df2002 += factor*d4fdraragbgb;
662     ds->df20011 += factor*d4fdraragbab;
663     ds->df20002 += factor*d4fdraraabab;
664     ds->df1300 += factor*d4fdrarbrbrb;
665     ds->df1210 += factor*d4fdrarbrbga;
666     ds->df1201 += factor*d4fdrarbrbgb;
667     ds->df12001 += factor*d4fdrarbrbab;
668     ds->df1120 += factor*d4fdrarbgaga;
669     ds->df1111 += factor*d4fdrarbgagb;
670     ds->df11101 += factor*d4fdrarbgaab;
671     ds->df1102 += factor*d4fdrarbgbgb;
672     ds->df11011 += factor*d4fdrarbgbab;
673     ds->df11002 += factor*d4fdrarbabab;
674     ds->df1030 += factor*d4fdragagaga;
675     ds->df1021 += factor*d4fdragagagb;
676     ds->df10201 += factor*d4fdragagaab;
677     ds->df1012 += factor*d4fdragagbgb;
678     ds->df10111 += factor*d4fdragagbab;
679     ds->df10102 += factor*d4fdragaabab;
680     ds->df1003 += factor*d4fdragbgbgb;
681     ds->df10021 += factor*d4fdragbgbab;
682     ds->df10012 += factor*d4fdragbabab;
683     ds->df10003 += factor*d4fdraababab;
684     ds->df0400 += factor*d4fdrbrbrbrb;
685     ds->df0310 += factor*d4fdrbrbrbga;
686     ds->df0301 += factor*d4fdrbrbrbgb;
687     ds->df03001 += factor*d4fdrbrbrbab;
688     ds->df0220 += factor*d4fdrbrbgaga;
689     ds->df0211 += factor*d4fdrbrbgagb;
690     ds->df02101 += factor*d4fdrbrbgaab;
691     ds->df0202 += factor*d4fdrbrbgbgb;
692     ds->df02011 += factor*d4fdrbrbgbab;
693     ds->df02002 += factor*d4fdrbrbabab;
694     ds->df0130 += factor*d4fdrbgagaga;
695     ds->df0121 += factor*d4fdrbgagagb;
696     ds->df01201 += factor*d4fdrbgagaab;
697     ds->df0112 += factor*d4fdrbgagbgb;
698     ds->df01111 += factor*d4fdrbgagbab;
699     ds->df01102 += factor*d4fdrbgaabab;
700     ds->df0103 += factor*d4fdrbgbgbgb;
701     ds->df01021 += factor*d4fdrbgbgbab;
702     ds->df01012 += factor*d4fdrbgbabab;
703     ds->df01003 += factor*d4fdrbababab;
704     ds->df0040 += factor*d4fdgagagaga;
705     ds->df0031 += factor*d4fdgagagagb;
706     ds->df00301 += factor*d4fdgagagaab;
707     ds->df0022 += factor*d4fdgagagbgb;
708     ds->df00211 += factor*d4fdgagagbab;
709     ds->df00202 += factor*d4fdgagaabab;
710     ds->df0013 += factor*d4fdgagbgbgb;
711     ds->df00121 += factor*d4fdgagbgbab;
712     ds->df00112 += factor*d4fdgagbabab;
713     ds->df00103 += factor*d4fdgaababab;
714     ds->df0004 += factor*d4fdgbgbgbgb;
715     ds->df00031 += factor*d4fdgbgbgbab;
716     ds->df00022 += factor*d4fdgbgbabab;
717     ds->df00013 += factor*d4fdgbababab;
718     ds->df00004 += factor*d4fdabababab;
719 }
720