LibreOffice Module sc (master) 1
op_statistical.cxx
Go to the documentation of this file.
1/* -*- Mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- */
2/*
3 * This file is part of the LibreOffice project.
4 *
5 * This Source Code Form is subject to the terms of the Mozilla Public
6 * License, v. 2.0. If a copy of the MPL was not distributed with this
7 * file, You can obtain one at http://mozilla.org/MPL/2.0/.
8 */
9
10#include "op_statistical.hxx"
11
13#include <sstream>
15
16using namespace formula;
17
18#include "op_math_helpers.hxx"
19
20namespace sc::opencl {
21
22void OpZTest::BinInlineFun(std::set<std::string>& decls,
23 std::set<std::string>& funs)
24{
25 decls.insert(phiDecl);
26 funs.insert(phi);
27 decls.insert(taylorDecl);
28 funs.insert(taylor);
29 decls.insert(gaussDecl);
30 funs.insert(gauss);
31}
33 const std::string &sSymName, SubArguments &vSubArguments)
34{
36 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
37 ss << "{\n";
38 ss << " int gid0 = get_global_id(0);\n";
39 ss << " double fSum = 0.0;\n";
40 ss << " double fSumSqr = 0.0;\n";
41 ss << " double mue = 0.0;\n";
42 ss << " double fCount = 0.0;\n";
43 GenerateRangeArg( 0, vSubArguments, ss, SkipEmpty,
44 " fSum += arg;\n"
45 " fSumSqr += arg * arg;\n"
46 " fCount += 1.0;\n"
47 );
48 ss << " if(fCount <= 1.0)\n";
49 ss << " return CreateDoubleError(DivisionByZero);\n";
50 ss << " mue = fSum / fCount;\n";
51 GenerateArg( "mu", 1, vSubArguments, ss );
52 if(vSubArguments.size() == 3)
53 {
54 GenerateArg( "sigma", 2, vSubArguments, ss );
55 ss << " if(sigma <= 0.0)\n";
56 ss << " return CreateDoubleError(IllegalArgument);\n";
57 ss << " return 0.5 - gauss((mue-mu)*sqrt(fCount)/sigma);\n";
58 }
59 else
60 {
61 ss << " double sigma = (fSumSqr-fSum*fSum/fCount)/(fCount-1.0);\n";
62 ss << " if(sigma == 0.0)\n";
63 ss << " return CreateDoubleError(DivisionByZero);\n";
64 ss << " return 0.5 - gauss((mue-mu)/sqrt(sigma/fCount));\n";
65 }
66 ss << "}\n";
67}
68
69void OpTTest::BinInlineFun(std::set<std::string>& decls,
70 std::set<std::string>& funs)
71{
72 decls.insert(fMachEpsDecl);
73 funs.insert("");
74 decls.insert(fMaxGammaArgumentDecl);
75 funs.insert("");
76 decls.insert(lcl_getLanczosSumDecl);
77 funs.insert(lcl_getLanczosSum);
78 decls.insert(GetBetaDecl);
79 funs.insert(GetBeta);
80 decls.insert(GetLogBetaDecl);
81 funs.insert(GetLogBeta);
82 decls.insert(GetBetaDistPDFDecl);
83 funs.insert(GetBetaDistPDF);
85 funs.insert(lcl_GetBetaHelperContFrac);
86 decls.insert(GetBetaDistDecl);
87 funs.insert(GetBetaDist);
88 decls.insert(GetTDistDecl);
89 funs.insert(GetTDist);
90}
91
93 const std::string &sSymName, SubArguments &vSubArguments)
94{
96 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
97 ss << "{\n";
98 ss << " int gid0 = get_global_id(0);\n";
99 ss << " double fSum1 = 0.0;\n";
100 ss << " double fSum2 = 0.0;\n";
101 ss << " double fSumSqr1 = 0.0;\n";
102 ss << " double fSumSqr2 = 0.0;\n";
103 ss << " double fCount1 = 0.0;\n";
104 ss << " double fCount2 = 0.0;\n";
105 ss << " double fT = 0.0;\n";
106 ss << " double fF = 0.0;\n";
107 GenerateArg( "mode", 2, vSubArguments, ss );
108 GenerateArg( "type", 3, vSubArguments, ss );
109 ss << " mode = floor(mode);\n";
110 ss << " type = floor(type);\n";
111 ss << " if(mode != 1.0 && mode != 2.0)\n";
112 ss << " return CreateDoubleError(IllegalArgument);\n";
113 ss << " if(type != 1.0 && type != 2.0 && type != 3.0)\n";
114 ss << " return CreateDoubleError(IllegalArgument);\n";
115
116 ss << " if(type == 1.0)\n";
117 ss << " {\n";
118 GenerateRangeArgPair( 0, 1, vSubArguments, ss, SkipEmpty,
119 " fSum1 += arg1;\n"
120 " fSum2 += arg2;\n"
121 " fSumSqr1 += (arg1 - arg2)*(arg1 - arg2);\n"
122 " fCount1 += 1;\n"
123 );
124 ss << " if(fCount1 < 1.0)\n";
125 ss << " return CreateDoubleError(NoValue);\n";
126 ss << " double divider = sqrt(fCount1 * fSumSqr1 - (fSum1-fSum2)*(fSum1-fSum2));\n";
127 ss << " if(divider == 0)\n";
128 ss << " return CreateDoubleError(DivisionByZero);\n";
129 ss << " fT = sqrt(fCount1-1.0) * fabs(fSum1 - fSum2) / divider;\n";
130 ss << " fF = fCount1 - 1.0;\n";
131 ss << " }\n";
132 ss << " if(type == 2.0 || type == 3.0)\n";
133 ss << " {\n";
134 GenerateRangeArg( 0, vSubArguments, ss, SkipEmpty,
135 " fSum1 += arg;\n"
136 " fSumSqr1 += arg * arg;\n"
137 " fCount1 += 1;\n"
138 );
139 GenerateRangeArg( 1, vSubArguments, ss, SkipEmpty,
140 " fSum2 += arg;\n"
141 " fSumSqr2 += arg * arg;\n"
142 " fCount2 += 1;\n"
143 );
144 ss << " if (fCount1 < 2.0 || fCount2 < 2.0)\n";
145 ss << " return CreateDoubleError(NoValue);\n";
146 ss << " }\n";
147 ss << " if(type == 3.0)\n";
148 ss << " {\n";
149 ss << " double fS1 = (fSumSqr1-fSum1*fSum1/fCount1)\n";
150 ss << " /(fCount1-1.0)/fCount1;\n";
151 ss << " double fS2 = (fSumSqr2-fSum2*fSum2/fCount2)\n";
152 ss << " /(fCount2-1.0)/fCount2;\n";
153 ss << " if (fS1 + fS2 == 0.0)\n";
154 ss << " return CreateDoubleError(NoValue);\n";
155 ss << " fT = fabs(fSum1/fCount1 - fSum2/fCount2)\n";
156 ss << " /sqrt(fS1+fS2);\n";
157 ss << " double c = fS1/(fS1+fS2);\n";
158 ss << " fF = 1.0/(c*c/(fCount1-1.0)+(1.0-c)*(1.0-c)\n";
159 ss << " /(fCount2-1.0));\n";
160 ss << " }\n";
161 ss << " if(type == 2.0)\n";
162 ss << " {\n";
163 ss << " double fS1 = (fSumSqr1 - fSum1*fSum1/fCount1)\n";
164 ss << " /(fCount1 - 1.0);\n";
165 ss << " double fS2 = (fSumSqr2 - fSum2*fSum2/fCount2)\n";
166 ss << " /(fCount2 - 1.0);\n";
167 ss << " fT = fabs( fSum1/fCount1 - fSum2/fCount2 )\n";
168 ss << " /sqrt( (fCount1-1.0)*fS1 + (fCount2-1.0)*fS2 )\n";
169 ss << " *sqrt( fCount1*fCount2*(fCount1+fCount2-2)\n";
170 ss << " /(fCount1+fCount2) );\n";
171 ss << " fF = fCount1 + fCount2 - 2;\n";
172 ss << " }\n";
173 ss << " double tdist=GetTDist(fT, fF);\n";
174 ss << " if (mode==1)\n";
175 ss << " return tdist;\n";
176 ss << " else\n";
177 ss << " return 2.0*tdist;\n";
178 ss << "}\n";
179}
180
181void OpTDist::BinInlineFun(std::set<std::string>& decls,
182 std::set<std::string>& funs)
183{
184 decls.insert(fMachEpsDecl);
185 funs.insert("");
186 decls.insert(fMaxGammaArgumentDecl);
187 funs.insert("");
188 decls.insert(lcl_getLanczosSumDecl);
189 funs.insert(lcl_getLanczosSum);
190 decls.insert(GetBetaDecl);
191 funs.insert(GetBeta);
192 decls.insert(GetLogBetaDecl);
193 funs.insert(GetLogBeta);
194 decls.insert(GetBetaDistPDFDecl);
195 funs.insert(GetBetaDistPDF);
196 decls.insert(lcl_GetBetaHelperContFracDecl);
197 funs.insert(lcl_GetBetaHelperContFrac);
198 decls.insert(GetBetaDistDecl);
199 funs.insert(GetBetaDist);
200 decls.insert(GetTDistDecl);
201 funs.insert(GetTDist);
202}
204 const std::string &sSymName, SubArguments &vSubArguments)
205{
206 CHECK_PARAMETER_COUNT( 3, 3 );
207 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
208 ss << "{\n";
209 ss << " int gid0 = get_global_id(0);\n";
210 GenerateArg( "x", 0, vSubArguments, ss );
211 GenerateArg( "fDF", 1, vSubArguments, ss );
212 GenerateArg( "fFlag", 2, vSubArguments, ss );
213 ss << " fDF = floor( fDF );\n";
214 ss << " fFlag = floor( fFlag );\n";
215 ss << " if(fDF < 1.0 || x < 0.0 || (fFlag != 1.0 && fFlag != 2.0))\n";
216 ss << " return CreateDoubleError(IllegalArgument);\n";
217 ss << " double R = GetTDist(x, fDF);\n";
218 ss << " if (fFlag == 1.0)\n";
219 ss << " return R;\n";
220 ss << " else\n";
221 ss << " return 2.0 * R;\n";
222 ss << "}\n";
223}
224
226 const std::string &sSymName, SubArguments &vSubArguments)
227{
228 CHECK_PARAMETER_COUNT( 3, 3 );
229 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
230 ss << "{\n";
231 ss << " double tmp = 0;\n";
232 ss << " int gid0 = get_global_id(0);\n";
233 GenerateArg( "rx", 0, vSubArguments, ss );
234 GenerateArg( "rlambda", 1, vSubArguments, ss );
235 GenerateArg( "rkum", 2, vSubArguments, ss );
236 ss <<" if(rlambda <= 0.0)\n";
237 ss <<" return CreateDoubleError(IllegalArgument);\n";
238 ss <<" else if(rkum == 0)\n";
239 ss <<" {\n";
240 ss <<" if(rx >= 0)\n";
241 ss <<" tmp = rlambda*exp(-rlambda*rx);\n";
242 ss <<" else\n";
243 ss <<" tmp = 0.0;\n";
244 ss <<" }\n";
245 ss <<" else\n";
246 ss <<" {\n";
247 ss <<" if(rx > 0)\n";
248 ss <<" tmp = 1.0 - exp(-rlambda*rx);\n";
249 ss <<" else\n";
250 ss <<" tmp = 0.0;\n";
251 ss <<" }\n";
252 ss <<" return tmp;\n";
253 ss <<"}";
254}
255void OpFdist::BinInlineFun(std::set<std::string>& decls,
256 std::set<std::string>& funs)
257{
258 decls.insert(GetFDistDecl);decls.insert(GetBetaDistDecl);
259 decls.insert(GetBetaDecl);decls.insert(fMaxGammaArgumentDecl);
260 decls.insert(lcl_GetBetaHelperContFracDecl);
261 decls.insert(GetBetaDistPDFDecl);
262 decls.insert(GetLogBetaDecl);decls.insert(lcl_getLanczosSumDecl);
263 decls.insert(fMachEpsDecl);
264 funs.insert(GetFDist);funs.insert(GetBetaDist);
265 funs.insert(GetBeta);
266 funs.insert(lcl_GetBetaHelperContFrac);funs.insert(GetBetaDistPDF);
267 funs.insert(GetLogBeta);
268 funs.insert(lcl_getLanczosSum);
269}
271 const std::string &sSymName, SubArguments &vSubArguments)
272{
273 CHECK_PARAMETER_COUNT( 3, 3 );
274 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
275 ss << "{\n";
276 ss << " double tmp = 0;\n";
277 ss << " int gid0 = get_global_id(0);\n";
278 GenerateArg( "rX", 0, vSubArguments, ss );
279 GenerateArg( "rF1", 1, vSubArguments, ss );
280 GenerateArg( "rF2", 2, vSubArguments, ss );
281 ss <<" rF1 = floor(rF1);\n";
282 ss <<" rF2 = floor(rF2);\n";
283 ss <<" if (rX < 0.0 || rF1 < 1.0 || rF2 < 1.0 || rF1 >= 1.0E10 ||";
284 ss <<"rF2 >= 1.0E10)\n";
285 ss <<" return CreateDoubleError(IllegalArgument);\n";
286 ss <<" tmp = GetFDist(rX, rF1, rF2);\n";
287 ss <<" return tmp;\n";
288 ss <<"}";
289}
290
292 const std::string &sSymName, SubArguments &vSubArguments)
293{
294 CHECK_PARAMETER_COUNT( 3, 3 );
295 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
296 ss << "{\n";
297 ss << " int gid0 = get_global_id(0);\n";
298 GenerateArg( "x", 0, vSubArguments, ss );
299 GenerateArg( "mu", 1, vSubArguments, ss );
300 GenerateArg( "sigma", 2, vSubArguments, ss );
301 ss << " if(sigma < 0.0)\n";
302 ss << " return CreateDoubleError(IllegalArgument);\n";
303 ss << " else if(sigma == 0.0)\n";
304 ss << " return CreateDoubleError(DivisionByZero);\n";
305 ss << " else\n";
306 ss << " return (x - mu)/sigma;\n";
307 ss << "}";
308}
309
311 const std::string &sSymName, SubArguments &vSubArguments)
312{
313 CHECK_PARAMETER_COUNT( 4, 4 );
314 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
315 ss << "{\n";
316 ss << " int gid0 = get_global_id(0);\n";
317 GenerateArg( "x", 0, vSubArguments, ss );
318 GenerateArg( "alpha", 1, vSubArguments, ss );
319 GenerateArg( "beta", 2, vSubArguments, ss );
320 GenerateArg( "kum", 3, vSubArguments, ss );
321 ss << " if(alpha <= 0.0 || beta <=0.0 || x < 0.0)\n";
322 ss << " return CreateDoubleError(IllegalArgument);\n";
323 ss << " if (kum == 0.0)\n";
324 ss << " return alpha/pow(beta,alpha)*pow(x,alpha-1.0)*\n";
325 ss << " exp(-pow(x/beta,alpha));\n";
326 ss << " else\n";
327 ss << " return 1.0 - exp(-pow(x/beta,alpha));\n";
328 ss << "}\n";
329}
330
331void OpTInv::BinInlineFun(std::set<std::string>& decls,
332 std::set<std::string>& funs)
333{
334 decls.insert(fMachEpsDecl);
335 funs.insert("");
336 decls.insert(fMaxGammaArgumentDecl);
337 funs.insert("");
338 decls.insert(lcl_getLanczosSumDecl);
339 funs.insert(lcl_getLanczosSum);
340 decls.insert(GetBetaDecl);
341 funs.insert(GetBeta);
342 decls.insert(GetLogBetaDecl);
343 funs.insert(GetLogBeta);
344 decls.insert(GetBetaDistPDFDecl);
345 funs.insert(GetBetaDistPDF);
346 decls.insert(lcl_GetBetaHelperContFracDecl);
347 funs.insert(lcl_GetBetaHelperContFrac);
348 decls.insert(GetBetaDistDecl);
349 funs.insert(GetBetaDist);
350 decls.insert(GetTDistDecl);
351 funs.insert(GetTDist);
352 decls.insert(GetValueDecl);
353 funs.insert(GetValue);
354 decls.insert(lcl_HasChangeOfSignDecl);
355 funs.insert(lcl_HasChangeOfSign);
356 decls.insert(lcl_IterateInverseDecl);
357 funs.insert(lcl_IterateInverse);
358}
359
361 const std::string &sSymName, SubArguments &vSubArguments)
362{
363 CHECK_PARAMETER_COUNT( 2, 2 );
364 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
365 ss << "{\n";
366 ss << " int gid0 = get_global_id(0);\n";
367 GenerateArg( "x", 0, vSubArguments, ss );
368 GenerateArg( "fDF", 1, vSubArguments, ss );
369 ss << " fDF = floor(fDF);\n";
370 ss << " if (x > 1.0||fDF < 1.0 || fDF > 1.0E10 || x <= 0.0)\n";
371 ss << " return CreateDoubleError(IllegalArgument);\n";
372 ss << " bool bConvError;\n";
373 ss << " double fVal = lcl_IterateInverse(\n";
374 ss << " fDF*0.5, fDF, &bConvError,x,fDF );\n";
375 ss << " if (bConvError)\n";
376 ss << " return CreateDoubleError(IllegalArgument);\n";
377 ss << " return fVal;\n";
378 ss << "}\n";
379}
380
382 const std::string &sSymName, SubArguments &vSubArguments)
383{
384 CHECK_PARAMETER_COUNT( 1, 1 );
385 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
386 ss << "{\n";
387 ss << " int gid0=get_global_id(0);\n";
388 GenerateArg( 0, vSubArguments, ss );
389 ss << " if (fabs(arg0) >= 1.0)\n";
390 ss << " return CreateDoubleError(IllegalArgument);\n";
391 ss << " double tmp=0.5*log((1+arg0)/(1-arg0));\n";
392 ss << " return tmp;\n";
393 ss << "}\n";
394}
395
397 outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments)
398{
399 CHECK_PARAMETER_COUNT( 1, 1 );
400 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
401 ss << "{\n";
402 ss << " int gid0=get_global_id(0);\n";
403 GenerateArg( 0, vSubArguments, ss );
404 ss << " double tmp=tanh(arg0);\n";
405 ss << " return tmp;\n";
406 ss << "}\n";
407}
408
410 outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments)
411{
412 CHECK_PARAMETER_COUNT( 1, 1 );
413 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
414 ss << "{\n";
415 ss <<" int gid0=get_global_id(0);\n";
416 GenerateArg( 0, vSubArguments, ss );
417 ss << " double tmp=tgamma(arg0);\n";
418 ss << " return tmp;\n";
419 ss << "}\n";
420}
421
423 outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments)
424{
425 CHECK_PARAMETER_COUNT( 3, 3 );
426 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
427 ss << "{\n";
428 ss << " int gid0=get_global_id(0);\n";
429 GenerateArg( "f", 0, vSubArguments, ss );
430 GenerateArg( "s", 1, vSubArguments, ss );
431 GenerateArg( "p", 2, vSubArguments, ss );
432 ss << " f = floor( f );\n";
433 ss << " s = floor( s );\n";
434 ss << " if ((f + s) <= 1.0 || p < 0.0 || p > 1.0)\n";
435 ss << " return CreateDoubleError(IllegalArgument);\n";
436 ss << " double q = 1.0 - p;\n";
437 ss << " double fFactor = pow(p,s);\n";
438 ss << " for(int i=0; i<f; i++)\n";
439 ss << " fFactor *= (i+s)/(i+1.0)*q;\n";
440 ss << " return fFactor;\n";
441 ss << "}\n";
442}
443
445 outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments)
446{
447 CHECK_PARAMETER_COUNT( 1, 1 );
448 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
449 ss << "{\n\t";
450 ss <<"int gid0=get_global_id(0);\n\t";
451 GenerateArg( 0, vSubArguments, ss );
452 ss << "double tmp=lgamma(arg0);\n\t";
453 ss << "return tmp;\n";
454 ss << "}\n";
455}
456void OpGauss::BinInlineFun(std::set<std::string>& decls,
457 std::set<std::string>& funs)
458{
459 decls.insert(taylorDecl);decls.insert(phiDecl);
460 decls.insert(gaussDecl);
461 funs.insert(taylor);funs.insert(phi);
462 funs.insert(gauss);
463}
464
466 outputstream &ss, const std::string &sSymName, SubArguments &
467vSubArguments)
468{
469 CHECK_PARAMETER_COUNT( 1, 1 );
470 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
471 ss << "{\n";
472 ss <<" int gid0=get_global_id(0);\n";
473 GenerateArg( 0, vSubArguments, ss );
474 ss << " double tmp=gauss(arg0);\n";
475 ss << " return tmp;\n";
476 ss << "}\n";
477}
478
480 outputstream &ss, const std::string &sSymName, SubArguments &
481vSubArguments)
482{
483 CHECK_PARAMETER_COUNT( 1, 30 );
484 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
485 ss << "{\n";
486 ss << " int gid0 = get_global_id(0);\n";
487 ss << " double nVal=0.0;\n";
488 ss << " double tmp = 0;\n";
489 ss << " int length;\n";
490 ss << " int totallength=0;\n";
491 GenerateRangeArgs( vSubArguments, ss, SkipEmpty,
492 " if( arg < 0 )\n"
493 " return CreateDoubleError(IllegalArgument);\n"
494 " if( arg == 0 )\n"
495 " return 0;\n"
496 " nVal += log(arg);\n"
497 " ++totallength;\n"
498 );
499 ss << " return exp(nVal/totallength);\n";
500 ss << "}";
501}
502
504 outputstream &ss, const std::string &sSymName, SubArguments &
505vSubArguments)
506{
507 CHECK_PARAMETER_COUNT( 1, 30 );
508 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
509 ss << "{\n";
510 ss << " int gid0 = get_global_id(0);\n";
511 ss << " double nVal=0.0;\n";
512 ss << " double tmp = 0;\n";
513 ss << " int length;\n";
514 ss << " int totallength=0;\n";
515 GenerateRangeArgs( vSubArguments, ss, SkipEmpty,
516 " if( arg <= 0 )\n"
517 " return CreateDoubleError(IllegalArgument);\n"
518 " nVal += (1.0 / arg);\n"
519 " ++totallength;\n"
520 );
521 ss << " return totallength/nVal;\n";
522 ss << "}";
523}
524
525void OpConfidence::BinInlineFun(std::set<std::string>& decls,
526 std::set<std::string>& funs)
527{
528 decls.insert(gaussinvDecl);
529 funs.insert(gaussinv);
530}
531
533 const std::string &sSymName, SubArguments& vSubArguments)
534{
535 CHECK_PARAMETER_COUNT( 3, 3 );
536 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
537 ss << "{\n";
538 ss << " double tmp = " << GetBottom() <<";\n";
539 ss << " int gid0 = get_global_id(0);\n";
540 GenerateArg( "alpha", 0, vSubArguments, ss );
541 GenerateArg( "sigma", 1, vSubArguments, ss );
542 GenerateArg( "size", 2, vSubArguments, ss );
543 ss << " double rn = floor(size);\n";
544 ss << " if(sigma <= 0.0 || alpha <= 0.0 || alpha >= 1.0";
545 ss << "|| rn < 1.0)\n";
546 ss << " return CreateDoubleError(IllegalArgument);\n";
547 ss << " else\n";
548 ss << " tmp = gaussinv(1.0 - alpha / 2.0) * sigma / sqrt( rn );\n";
549 ss << " return tmp;\n";
550 ss << "}";
551}
552
553void OpCritBinom::BinInlineFun(std::set<std::string>& decls,
554 std::set<std::string>& funs)
555{
556 decls.insert(MinDecl);
557 funs.insert("");
558}
559
561 const std::string &sSymName, SubArguments& vSubArguments)
562{
563 CHECK_PARAMETER_COUNT( 3, 3 );
564 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
565 ss << "{\n";
566 ss << " double tmp = " << GetBottom() <<";\n";
567 ss << " int gid0 = get_global_id(0);\n";
568 GenerateArg( "n", 0, vSubArguments, ss );
569 GenerateArg( "p", 1, vSubArguments, ss );
570 GenerateArg( "alpha", 2, vSubArguments, ss );
571 ss << " double rn = floor(n);\n";
572 ss << " if (rn < 0.0 || alpha < 0.0 || alpha > 1.0 || p < 0.0";
573 ss << " || p > 1.0)\n";
574 ss << " return CreateDoubleError(IllegalArgument);\n";
575 ss << " else if ( alpha == 0 )\n";
576 ss << " return 0;\n";
577 ss << " else if ( alpha == 1 )\n";
578 ss << " return p == 0 ? 0 : rn;\n";
579 ss << " else\n";
580 ss << " {\n";
581 ss << " double rq = (0.5 - p) + 0.5;\n";
582 ss << " double fFactor = pow(rq, rn);\n";
583 ss << " if (fFactor <= Min)\n";
584 ss << " {\n";
585 ss << " fFactor = pow(p, rn);\n";
586 ss << " if (fFactor <= Min)\n";
587 ss << " return CreateDoubleError(NoValue);\n";
588 ss << " else\n";
589 ss << " {\n";
590 ss << " double fSum = 1.0 - fFactor;\n";
591 ss << " uint max =(uint)(rn), i;\n";
592 ss << " for (i = 0; i < max && fSum >= alpha; i++)\n";
593 ss << " {\n";
594 ss << " fFactor *= (rn - i) / (double)(i + 1) * rq / p;\n";
595 ss << " fSum -= fFactor;\n";
596 ss << " }\n";
597 ss << " tmp = (rn - i);\n";
598 ss << " }\n";
599 ss << " }\n";
600 ss << " else\n";
601 ss << " {\n";
602 ss << " double fSum = fFactor;\n";
603 ss << " uint max = (uint)(rn), i;\n";
604 ss << " for (i = 0; i < max && fSum < alpha; i++)\n";
605 ss << " {\n";
606 ss << " fFactor *= (rn - i) / (double)(i + 1) *";
607 ss << " p / rq;\n";
608 ss << " fSum += fFactor;\n";
609 ss << " }\n";
610 ss << " tmp = (i);\n";
611 ss << " }\n";
612 ss << " }\n";
613 ss << " return tmp;\n";
614 ss << "}";
615}
616
617void OpChiInv::BinInlineFun(std::set<std::string>& decls,
618 std::set<std::string>& funs)
619{
620 decls.insert(fMachEpsDecl);
621 funs.insert("");
622 decls.insert(fBigInvDecl);
623 funs.insert("");
624 decls.insert(fHalfMachEpsDecl);
625 funs.insert("");
626 decls.insert(lcl_IterateInverseChiInvDecl);
627 funs.insert(lcl_IterateInverseChiInv);
628 decls.insert(GetChiDistDecl);
629 funs.insert(GetChiDist);
630 decls.insert(lcl_HasChangeOfSignDecl);
631 funs.insert(lcl_HasChangeOfSign);
632 decls.insert(GetUpRegIGammaDecl);
633 funs.insert(GetUpRegIGamma);
634 decls.insert(GetGammaContFractionDecl);
635 funs.insert(GetGammaContFraction);
636 decls.insert(GetGammaSeriesDecl);
637 funs.insert(GetGammaSeries);
638}
640 outputstream &ss,const std::string &sSymName,
641 SubArguments &vSubArguments)
642{
643 CHECK_PARAMETER_COUNT( 2, 2 );
644 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
645 ss << "{\n";
646 ss << " double tmp;\n";
647 ss << " int gid0=get_global_id(0);\n";
648 ss <<"\n ";
649 GenerateArg( "tmp0", 0, vSubArguments, ss );
650 GenerateArg( "tmp1", 1, vSubArguments, ss );
651 ss << " tmp1 = floor(tmp1);";
652 ss << " if (tmp1 < 1.0 || tmp0 <= 0.0 || tmp0 > 1.0 )\n";
653 ss << " {\n";
654 ss << " return CreateDoubleError(IllegalArgument);\n";
655 ss << " }\n";
656 ss << " bool bConvError;\n";
657 ss << " double fVal = lcl_IterateInverseChiInv";
658 ss << "(tmp0, tmp1, tmp1*0.5, tmp1, &bConvError);\n";
659 ss << " if(bConvError)\n";
660 ss << " return CreateDoubleError(NoConvergence);\n";
661 ss << " return fVal;\n";
662 ss << "}\n";
663}
665 outputstream &ss, const std::string &sSymName,
666 SubArguments &vSubArguments)
667{
669 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
670 ss << "{\n";
671 ss << " int gid0=get_global_id(0);\n";
672 GenerateArg( "x", 0, vSubArguments, ss );
673 GenerateArg( "mue", 1, vSubArguments, ss );
674 GenerateArg( "sigma", 2, vSubArguments, ss );
675 GenerateArg( "c", 3, vSubArguments, ss );
676 ss << "if(sigma <= 0)\n";
677 ss << " return CreateDoubleError(IllegalArgument);\n";
678 ss << "double mid,tmp;\n";
679 ss << "mid = (x - mue)/sigma;\n";
680 ss << "if(c)\n";
681 ss << " tmp = 0.5 *erfc(-mid * 0.7071067811865475);\n";
682 ss << "else \n";
683 ss <<" tmp=(0.39894228040143268*exp(-pow(mid,2)/2.0))/sigma;\n";
684 ss << "return tmp;\n";
685 ss << "}\n";
686}
688 outputstream &ss,const std::string &sSymName,
689 SubArguments &vSubArguments)
690{
691 CHECK_PARAMETER_COUNT( 1, 1 );
692 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
693 ss << "{\n";
694 ss << " int gid0=get_global_id(0);\n";
695 GenerateArg( "x", 0, vSubArguments, ss );
696 ss << " double tmp = 0.5 * erfc((-1)*x * 0.7071067811865475);\n";
697 ss << " return tmp;\n";
698 ss << "}\n";
699}
700
702 outputstream &ss,const std::string &sSymName,
703 SubArguments &vSubArguments)
704{
705 CHECK_PARAMETER_COUNT( 2, 2 );
706 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
707 ss << "{\n";
708 ss <<" int gid0=get_global_id(0);\n";
709 ss <<" double tmp = 1 ;\n";
710 GenerateArg( "inA", 0, vSubArguments, ss );
711 GenerateArg( "inB", 1, vSubArguments, ss );
712 ss << " inA = floor( inA );\n";
713 ss << " inB = floor( inB );\n";
714 ss << " if (inA < 0.0 || inB < 0.0 || inB > inA)\n";
715 ss << " return CreateDoubleError(IllegalArgument);\n";
716 ss << " for( int i = 0; i<inB; i++)\n";
717 ss << " {\n";
718 ss << " tmp *= inA ;\n";
719 ss << " inA = inA - 1.0;\n";
720 ss << " }\n";
721 ss << " return tmp;\n";
722 ss << "}\n";
723}
725 outputstream &ss,const std::string &sSymName,
726 SubArguments &vSubArguments)
727{
728 CHECK_PARAMETER_COUNT( 2, 2 );
729 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
730 ss << "{\n";
731 ss <<" int gid0=get_global_id(0);\n";
732 ss <<" double tmp = 1.0;\n";
733 GenerateArg( "inA", 0, vSubArguments, ss );
734 GenerateArg( "inB", 1, vSubArguments, ss );
735 ss << " inA = floor( inA );\n";
736 ss << " inB = floor( inB );\n";
737 ss << " if (inA < 0.0 || inB < 0.0)\n";
738 ss << " return CreateDoubleError(IllegalArgument);\n";
739 ss << " return pow(inA, inB);\n";
740 ss << "}\n";
741}
742
744 outputstream &ss,const std::string &sSymName,
745 SubArguments &vSubArguments)
746{
747 CHECK_PARAMETER_COUNT( 1, 1 );
748 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
749 ss << "{\n";
750 ss << " int gid0=get_global_id(0);\n";
751 GenerateArg( "x", 0, vSubArguments, ss );
752 ss << " double tmp = 0.39894228040143268 * exp((-1)*pow(x,2) / 2.0);\n";
753 ss << " return tmp;\n";
754 ss << "}\n";
755}
756
757void OpNorminv::BinInlineFun(std::set<std::string>& decls,
758 std::set<std::string>& funs)
759{
760 decls.insert(gaussinvDecl);
761 funs.insert(gaussinv);
762}
763
765 outputstream &ss,const std::string &sSymName,
766 SubArguments &vSubArguments)
767{
768 CHECK_PARAMETER_COUNT( 3, 3 );
769 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
770 ss << "{\n";
771 ss << " int gid0=get_global_id(0);\n";
772 GenerateArg( "x", 0, vSubArguments, ss );
773 GenerateArg( "mue", 1, vSubArguments, ss );
774 GenerateArg( "sigma", 2, vSubArguments, ss );
775 ss << " if (sigma <= 0.0 || x < 0.0 || x > 1.0)\n";
776 ss << " return CreateDoubleError(IllegalArgument);\n";
777 ss << " else if (x == 0.0 || x == 1.0)\n";
778 ss << " return CreateDoubleError(NoValue);\n";
779 ss << " return gaussinv(x)*sigma + mue;\n";
780 ss << "}\n";
781}
782
783void OpNormsinv::BinInlineFun(std::set<std::string>& decls,
784 std::set<std::string>& funs)
785{
786 decls.insert(gaussinvDecl);
787 funs.insert(gaussinv);
788}
789
791 (outputstream &ss,const std::string &sSymName,
792 SubArguments &vSubArguments)
793{
794 CHECK_PARAMETER_COUNT( 1, 1 );
795 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
796 ss << "{\n";
797 ss << " int gid0=get_global_id(0);\n";
798 GenerateArg( "x", 0, vSubArguments, ss );
799 ss << " if (x < 0.0 || x > 1.0)\n";
800 ss << " return CreateDoubleError(IllegalArgument);\n";
801 ss << " else if (x == 0.0 || x == 1.0)\n";
802 ss << " return CreateDoubleError(NoValue);\n";
803 ss << " return gaussinv(x);\n";
804 ss << "}\n";
805}
806
807void OpLogInv::BinInlineFun(std::set<std::string>& decls,
808 std::set<std::string>& funs)
809{
810 decls.insert(gaussinvDecl);
811 funs.insert(gaussinv);
812}
813
815 const std::string &sSymName, SubArguments &vSubArguments)
816{
817 CHECK_PARAMETER_COUNT( 3, 3 );
818 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
819 ss << "{\n";
820 ss << " int gid0=get_global_id(0);\n";
821 ss << " double tmp;\n";
822 GenerateArg( "x", 0, vSubArguments, ss );
823 GenerateArgWithDefault( "mue", 1, 0, vSubArguments, ss );
824 GenerateArgWithDefault( "sigma", 2, 1, vSubArguments, ss );
825 ss << " if ( sigma <= 0.0 || x <= 0.0 || x >= 1.0 )\n";
826 ss << " return CreateDoubleError(IllegalArgument);\n";
827 ss << " return exp(mue+sigma*gaussinv(x));\n";
828 ss << "}\n";
829}
830
832 const std::string &sSymName, SubArguments &vSubArguments)
833{
834 CHECK_PARAMETER_COUNT( 2, 4 );
835 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
836 ss << "{\n";
837 ss << " int gid0=get_global_id(0);\n";
838 GenerateArg( "x", 0, vSubArguments, ss );
839 GenerateArgWithDefault( "mue", 1, 0, vSubArguments, ss );
840 GenerateArgWithDefault( "sigma", 2, 1, vSubArguments, ss );
841 GenerateArgWithDefault( "fCumulative", 3, 1, vSubArguments, ss );
842 ss << " if (sigma <= 0.0)\n";
843 ss << " return CreateDoubleError(IllegalArgument);\n";
844 ss << " double tmp;\n";
845 ss << " double temp = (log(x)-mue)/sigma;\n";
846 ss << " if(fCumulative != 0)\n";
847 ss << " {\n";
848 ss << " if(x<=0)\n";
849 ss << " tmp = 0.0;\n";
850 ss << " else\n";
851 ss << " tmp = 0.5 * erfc(-temp * 0.7071067811865475);\n";
852 ss << " }\n";
853 ss << " else\n";
854 ss << " if(x<=0)\n";
855 ss << " return CreateDoubleError(IllegalArgument);\n";
856 ss << " else\n";
857 ss << " tmp = (0.39894228040143268 * exp((-1)*pow(temp, 2)";
858 ss << " / 2.0))/(sigma*x);\n";
859 ss << " return tmp;\n";
860 ss << "}\n";
861}
862
863void OpGammaDist::BinInlineFun(std::set<std::string>& decls,
864 std::set<std::string>& funs)
865{
866 decls.insert(fBigInvDecl);decls.insert(fLogDblMaxDecl);
867 decls.insert(fHalfMachEpsDecl);decls.insert(fMaxGammaArgumentDecl);
868 decls.insert(GetGammaSeriesDecl);decls.insert(GetGammaContFractionDecl);
869 decls.insert(GetLowRegIGammaDecl);decls.insert(GetGammaDistDecl);
870 decls.insert(GetGammaDistPDFDecl);
871 funs.insert(GetGammaSeries);funs.insert(GetGammaContFraction);
872 funs.insert(GetLowRegIGamma);funs.insert(GetGammaDist);
873 funs.insert(GetGammaDistPDF);
874}
875
877 const std::string &sSymName, SubArguments &vSubArguments)
878{
879 CHECK_PARAMETER_COUNT( 3, 4 );
880 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
881 ss << "{\n";
882 ss << " int gid0=get_global_id(0);\n";
883 GenerateArg( 0, vSubArguments, ss );
884 GenerateArg( 1, vSubArguments, ss );
885 GenerateArg( 2, vSubArguments, ss );
886 GenerateArgWithDefault( "arg3", 3, 1, vSubArguments, ss );
887 ss << " if(arg1 <= 0 || arg2 <= 0)\n";
888 ss << " return CreateDoubleError(IllegalArgument);\n";
889 ss << " double tmp;\n";
890 ss << " if (arg3)\n";
891 ss << " tmp=GetGammaDist( arg0, arg1, arg2);\n";
892 ss << " else\n";
893 ss << " tmp=GetGammaDistPDF( arg0, arg1, arg2);\n";
894 ss << " return tmp;\n";
895 ss << "}\n";
896}
897void OpChiDist::BinInlineFun(std::set<std::string>& decls,
898 std::set<std::string>& funs)
899{
900 decls.insert(fBigInvDecl);
901 funs.insert("");
902 decls.insert(fHalfMachEpsDecl);
903 funs.insert("");
904 decls.insert(GetUpRegIGammaDecl);
905 funs.insert(GetUpRegIGamma);
906 decls.insert(GetGammaSeriesDecl);
907 funs.insert(GetGammaSeries);
908 decls.insert(GetGammaContFractionDecl);
909 funs.insert(GetGammaContFraction);
910 decls.insert(GetChiDistDecl);
911 funs.insert(GetChiDist);
912}
914 outputstream &ss,const std::string &sSymName,
915 SubArguments &vSubArguments)
916{
917 CHECK_PARAMETER_COUNT( 2, 2 );
918 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
919 ss << "{\n";
920 ss << " double fx,fDF,tmp=0;\n";
921 ss << " int gid0=get_global_id(0);\n";
922 GenerateArg( "tmp0", 0, vSubArguments, ss );
923 GenerateArg( "tmp1", 1, vSubArguments, ss );
924 ss << " fx = tmp0;\n";
925 ss << " fDF = floor(tmp1);\n";
926 ss << " if(fDF < 1.0)\n";
927 ss << " {\n";
928 ss << " return CreateDoubleError(IllegalArgument);\n";
929 ss << " }\n";
930 ss << " tmp = GetChiDist( fx, fDF);\n";
931 ss << " return tmp;\n";
932 ss << "}\n";
933}
934void OpBinomdist::BinInlineFun(std::set<std::string>& decls,
935 std::set<std::string>& funs)
936{
937 decls.insert(fMachEpsDecl);
938 funs.insert("");
939 decls.insert(MinDecl);
940 funs.insert("");
941 decls.insert(fMaxGammaArgumentDecl);
942 funs.insert("");
943 decls.insert(GetBinomDistPMFDecl);
944 funs.insert(GetBinomDistPMF);
945 decls.insert(GetBetaDistDecl);
946 funs.insert(GetBetaDist);
947 decls.insert(lcl_GetBinomDistRangeDecl);
948 funs.insert(lcl_GetBinomDistRange);
949 decls.insert(lcl_GetBetaHelperContFracDecl);
950 funs.insert(lcl_GetBetaHelperContFrac);
951 decls.insert(GetBetaDistPDFDecl);
952 funs.insert(GetBetaDistPDF);
953 decls.insert(GetLogBetaDecl);
954 funs.insert(GetLogBeta);
955 decls.insert(GetBetaDecl);
956 funs.insert(GetBeta);
957 decls.insert(lcl_getLanczosSumDecl);
958 funs.insert(lcl_getLanczosSum);
959}
961 outputstream &ss,const std::string &sSymName,
962 SubArguments &vSubArguments)
963{
964 CHECK_PARAMETER_COUNT( 4, 4 );
965 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
966 ss << "{\n";
967 ss << " int gid0=get_global_id(0);\n";
968 GenerateArg( "tmp0", 0, vSubArguments, ss );
969 GenerateArg( "tmp1", 1, vSubArguments, ss );
970 GenerateArg( "tmp2", 2, vSubArguments, ss );
971 GenerateArg( "tmp3", 3, vSubArguments, ss );
972 ss << " tmp0 = floor(tmp0);\n";
973 ss << " tmp1 = floor(tmp1);\n";
974 ss << " double rq = (0.5 - tmp2) + 0.5;\n";
975 ss << " if (tmp1 < 0.0 || tmp0 < 0.0 || tmp0 > tmp1 ||";
976 ss << "tmp2 < 0.0 || tmp2 > 1.0)\n";
977 ss << " {\n";
978 ss << " return CreateDoubleError(IllegalArgument);\n";
979 ss << " }\n";
980 ss << " if(tmp2 == 0.0)\n";
981 ss << " return ( (tmp0 == 0.0 || tmp3) ? 1.0 : 0.0 );\n";
982 ss << " if(tmp2 == 1.0)\n";
983 ss << " return ( (tmp0 == tmp1) ? 1.0 : 0.0);\n";
984 ss << " if(!tmp3)\n";
985 ss << " return ( GetBinomDistPMF(tmp0, tmp1, tmp2));\n";
986 ss << " else \n";
987 ss << " {\n";
988 ss << " if(tmp0 == tmp1)\n";
989 ss << " return 1.0;\n";
990 ss << " else\n";
991 ss << " {\n";
992 ss << " double fFactor = pow(rq,tmp1);\n";
993 ss << " if(tmp0 == 0.0)\n";
994 ss << " return (fFactor);\n";
995 ss << " else if(fFactor <= Min)\n";
996 ss << " {\n";
997 ss << " fFactor = pow(tmp2,tmp1);\n";
998 ss << " if(fFactor <= Min)\n";
999 ss << " return GetBetaDist";
1000 ss << "(rq, tmp1 - tmp0, tmp0 + 1.0);\n";
1001 ss << " else\n";
1002 ss << " {\n";
1003 ss << " if(fFactor > fMachEps)\n";
1004 ss << " {\n";
1005 ss << " double fSum = 1.0 - fFactor;\n";
1006 ss << " unsigned int max = ";
1007 ss << "(unsigned int)((tmp1 - tmp0)-1);\n";
1008 ss << " for (uint i = 0; i < max && fFactor > 0.0;";
1009 ss << " i++)\n";
1010 ss << " {\n";
1011 ss << " fFactor *= (tmp1 - i)/(i + 1)*rq/tmp2;\n";
1012 ss << " fSum -= fFactor;\n";
1013 ss << " }\n";
1014 ss << " return ( (fSum < 0.0) ? 0.0 : fSum );\n";
1015 ss << " }\n";
1016 ss << " else \n";
1017 ss << " return (lcl_GetBinomDistRange";
1018 ss << "(tmp1, tmp1 - tmp0, tmp1, fFactor, rq, tmp2));\n";
1019 ss << " }\n";
1020 ss << " }\n";
1021 ss << " else\n";
1022 ss << " {\n";
1023 ss << " double rtmp = ( lcl_GetBinomDistRange";
1024 ss << "(tmp1, 0.0, tmp0, fFactor, tmp2, rq));\n";
1025 ss << " return rtmp;\n";
1026 ss << " }\n";
1027 ss << " }\n";
1028 ss << " }\n";
1029 ss << "}\n";
1030}
1031
1032void OpChiSqDist::BinInlineFun(std::set<std::string>& decls,
1033 std::set<std::string>& funs)
1034{
1035 decls.insert(fMaxGammaArgumentDecl);decls.insert(GetChiSqDistCDFDecl);
1036 decls.insert(GetChiSqDistPDFDecl);decls.insert(GetLowRegIGammaDecl);
1037 decls.insert(GetGammaContFractionDecl);decls.insert(GetGammaSeriesDecl);
1038 decls.insert(fHalfMachEpsDecl);
1039 decls.insert(fBigInvDecl);
1040
1041 funs.insert(GetGammaContFraction);funs.insert(GetChiSqDistCDF);
1042 funs.insert(GetChiSqDistPDF);funs.insert(GetLowRegIGamma);
1043 funs.insert(GetGammaSeries);
1044}
1045
1047 outputstream &ss, const std::string &sSymName, SubArguments &
1048vSubArguments)
1049{
1050 CHECK_PARAMETER_COUNT( 2, 3 );
1051 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
1052 ss << "{\n";
1053 ss << " int gid0 = get_global_id(0);\n";
1054 ss << " double result = 0;\n";
1055 GenerateArg( "tmp0", 0, vSubArguments, ss );
1056 GenerateArg( "tmp1", 1, vSubArguments, ss );
1057 GenerateArgWithDefault( "tmp2", 2, 1, vSubArguments, ss );
1058 ss << " tmp1 = floor(tmp1);\n";
1059 ss << " if(tmp1 < 1.0)\n";
1060 ss << " return CreateDoubleError(IllegalArgument);\n";
1061 ss << " else\n";
1062 ss << " {\n";
1063 ss << " if(tmp2)\n";
1064 ss << " result =GetChiSqDistCDF(tmp0,tmp1);\n";
1065 ss << " else\n";
1066 ss << " result =GetChiSqDistPDF(tmp0,tmp1);\n";
1067 ss << " }\n";
1068 ss << " return result;\n";
1069 ss << "}";
1070}
1071
1072 void OpChiSqInv::BinInlineFun(std::set<std::string>& decls,
1073 std::set<std::string>& funs)
1074{
1075 decls.insert(fMaxGammaArgumentDecl);decls.insert(GetChiSqDistCDFDecl);
1076 decls.insert(GetLowRegIGammaDecl);decls.insert(lcl_IterateInverseChiSQInvDecl);
1077 decls.insert(GetGammaContFractionDecl);decls.insert(GetGammaSeriesDecl);
1078 decls.insert(fHalfMachEpsDecl);
1079 decls.insert(fBigInvDecl);decls.insert(lcl_HasChangeOfSignDecl);
1080 decls.insert(fMachEpsDecl);
1081
1082 funs.insert(GetGammaContFraction);funs.insert(GetChiSqDistCDF);
1083 funs.insert(GetLowRegIGamma);funs.insert(lcl_HasChangeOfSign);
1084 funs.insert(GetGammaSeries);funs.insert(lcl_IterateInverseChiSQInv);
1085}
1086
1088 outputstream &ss, const std::string &sSymName, SubArguments &
1089vSubArguments)
1090{
1091 CHECK_PARAMETER_COUNT( 2, 2 );
1092 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
1093 ss << "{\n";
1094 ss << " int gid0 = get_global_id(0);\n";
1095 ss << " double result = 0;\n";
1096 GenerateArg( "tmp0", 0, vSubArguments, ss );
1097 GenerateArg( "tmp1", 1, vSubArguments, ss );
1098 ss << " tmp1 = floor(tmp1);\n";
1099 ss << " bool bConvError;\n";
1100 ss << " if(tmp1 < 1.0 || tmp0 < 0 || tmp0>=1.0)\n";
1101 ss << " return CreateDoubleError(IllegalArgument);\n";
1102 ss << " else\n";
1103 ss << " {\n";
1104 ss << " result =lcl_IterateInverseChiSQInv( tmp0, tmp1,";
1105 ss << "tmp1*0.5, tmp1, &bConvError );\n";
1106 ss << " }\n";
1107 ss << " if(bConvError)\n";
1108 ss << " return CreateDoubleError(NoConvergence);\n";
1109 ss << " return result;\n";
1110 ss << "}";
1111}
1112
1113void OpGammaInv::BinInlineFun(std::set<std::string>& decls,
1114 std::set<std::string>& funs)
1115{
1116 decls.insert(fBigInvDecl);decls.insert(fHalfMachEpsDecl);
1117 decls.insert(GetGammaSeriesDecl);decls.insert(GetGammaContFractionDecl);
1118 decls.insert(GetGammaInvValueDecl);
1119 funs.insert(GetGammaSeries);funs.insert(GetGammaContFraction);
1120 funs.insert(GetGammaInvValue);
1121}
1122
1124 const std::string &sSymName, SubArguments &vSubArguments)
1125{
1126 CHECK_PARAMETER_COUNT( 3, 3 );
1127 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
1128 ss << "{\n";
1129 ss << " int gid0=get_global_id(0);\n";
1130 ss << " double tmp;\n";
1131 GenerateArg( 0, vSubArguments, ss );
1132 GenerateArg( 1, vSubArguments, ss );
1133 GenerateArg( 2, vSubArguments, ss );
1134 ss << " if( arg0 < 0 || arg0 >= 1 || arg1 <= 0 || arg2 <= 0 )\n";
1135 ss << " return CreateDoubleError(IllegalArgument);\n";
1136 ss << " if (arg0 == 0.0)\n"
1137 " {\n"
1138 " tmp=0.0;\n"
1139 " return tmp;\n"
1140 " }\n"
1141 " else\n"
1142 " {\n"
1143 " bool bConvError;\n"
1144 " double fStart = arg1 * arg2;\n"
1145 " double fAx=fStart*0.5;\n"
1146 " double fBx=fStart;\n"
1147 " bConvError = false;\n"
1148 " double fYEps = 1.0E-307;\n"
1149 " double fXEps = 2.22045e-016;\n"
1150 " double fAy = arg0-GetGammaInvValue(arg1,arg2,fAx);\n"
1151 " double fBy = arg0-GetGammaInvValue(arg1,arg2,fBx);\n"
1152 " double fTemp;\n"
1153 " unsigned short nCount;\n"
1154 " for (nCount = 0; nCount < 1000 && !((fAy < 0.0 && fBy > 0.0)"
1155 " || (fAy > 0.0 && fBy < 0.0)); nCount++)\n"
1156 " {\n"
1157 " if (fabs(fAy) <= fabs(fBy))\n"
1158 " {\n"
1159 " fTemp = fAx;\n"
1160 " fAx += 2.0 * (fAx - fBx);\n"
1161 " if (fAx < 0.0)\n"
1162 " fAx = 0.0;\n"
1163 " fBx = fTemp;\n"
1164 " fBy = fAy;\n"
1165 " fAy = arg0-GetGammaInvValue(arg1,arg2,fAx);\n"
1166 " }\n"
1167 " else\n"
1168 " {\n"
1169 " fTemp = fBx;\n"
1170 " fBx += 2.0 * (fBx - fAx);\n"
1171 " fAx = fTemp;\n"
1172 " fAy = fBy;\n"
1173 " fBy = arg0-GetGammaInvValue(arg1,arg2,fBx);\n"
1174 " }\n"
1175 " }\n"
1176 " if (fAy == 0.0)\n"
1177 " {\n"
1178 " tmp = fAx;\n"
1179 " return tmp;\n"
1180 " }\n"
1181 " if (fBy == 0.0)\n"
1182 " {\n"
1183 " tmp = fBx;\n"
1184 " return tmp;\n"
1185 " }\n"
1186 " if (!((fAy < 0.0 && fBy > 0.0) || (fAy > 0.0 && fBy < 0.0)))\n"
1187 " {\n"
1188 " bConvError = true;\n"
1189 " tmp = 0.0;\n"
1190 " return tmp;\n"
1191 " }\n"
1192 " double fPx = fAx;\n"
1193 " double fPy = fAy;\n"
1194 " double fQx = fBx;\n"
1195 " double fQy = fBy;\n"
1196 " double fRx = fAx;\n"
1197 " double fRy = fAy;\n"
1198 " double fSx = 0.5 * (fAx + fBx);\n"
1199 " bool bHasToInterpolate = true;\n"
1200 " nCount = 0;\n"
1201 " while ( nCount < 500 && fabs(fRy) > fYEps &&"
1202 "(fBx-fAx) > fmax( fabs(fAx), fabs(fBx)) * fXEps )\n"
1203 " {\n"
1204 " if (bHasToInterpolate)\n"
1205 " {\n"
1206 " if (fPy!=fQy && fQy!=fRy && fRy!=fPy)\n"
1207 " {\n"
1208 " fSx = fPx * fRy * fQy / (fRy-fPy) / (fQy-fPy)"
1209 "+ fRx * fQy * fPy / (fQy-fRy) / (fPy-fRy)"
1210 "+ fQx * fPy * fRy / (fPy-fQy) / (fRy-fQy);\n"
1211 " bHasToInterpolate = (fAx < fSx) && (fSx < fBx);\n"
1212 " }\n"
1213 " else\n"
1214 " bHasToInterpolate = false;\n"
1215 " }\n"
1216 " if(!bHasToInterpolate)\n"
1217 " {\n"
1218 " fSx = 0.5 * (fAx + fBx);\n"
1219 " fPx = fAx; fPy = fAy;\n"
1220 " fQx = fBx; fQy = fBy;\n"
1221 " bHasToInterpolate = true;\n"
1222 " }\n"
1223 " fPx = fQx; fQx = fRx; fRx = fSx;\n"
1224 " fPy = fQy; fQy = fRy;\n"
1225 " fRy = arg0-GetGammaInvValue(arg1,arg2,fSx);\n"
1226 " if ((fAy < 0.0 && fRy > 0.0) || (fAy > 0.0 && fRy < 0.0))\n"
1227 " {\n"
1228 " fBx = fRx;\n"
1229 " fBy = fRy;\n"
1230 " }\n"
1231 " else\n"
1232 " {\n"
1233 " fAx = fRx;\n"
1234 " fAy = fRy;\n"
1235 " }\n"
1236 " bHasToInterpolate = bHasToInterpolate && (fabs(fRy)"
1237 " * 2.0 <= fabs(fQy));\n"
1238 " ++nCount;\n"
1239 " }\n"
1240 " tmp = fRx;\n"
1241 " return tmp;\n"
1242 " }\n"
1243 "}\n";
1244}
1245void OpFInv::BinInlineFun(std::set<std::string>& decls,
1246 std::set<std::string>& funs)
1247{
1248 decls.insert(fMachEpsDecl);decls.insert(fMaxGammaArgumentDecl);
1249 decls.insert(lcl_getLanczosSumDecl);decls.insert(GetBetaDecl);
1250 decls.insert(GetLogBetaDecl);decls.insert(GetBetaDistPDFDecl);
1251 decls.insert(lcl_GetBetaHelperContFracDecl);decls.insert(GetFInvValueDecl);
1252 funs.insert(lcl_getLanczosSum);funs.insert(GetBeta);
1253 funs.insert(GetLogBeta);funs.insert(GetBetaDistPDF);
1254 funs.insert(lcl_GetBetaHelperContFrac);funs.insert(GetFInvValue);
1255}
1256
1258 const std::string &sSymName, SubArguments &vSubArguments)
1259{
1260 CHECK_PARAMETER_COUNT( 3, 3 );
1261 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
1262 ss << "{\n";
1263 ss << " int gid0=get_global_id(0);\n";
1264 ss << " double tmp;\n";
1265 GenerateArg( 0, vSubArguments, ss );
1266 GenerateArg( 1, vSubArguments, ss );
1267 GenerateArg( 2, vSubArguments, ss );
1268 ss << " double fF2=floor(arg2);\n"
1269 " double fF1=floor(arg1);\n"
1270 " if( arg0 <= 0 || arg1 < 1 || arg2 < 1 || arg1 >= 1.0e10 || arg2 >= 1.0e10 || arg0 > 1 )\n"
1271 " return CreateDoubleError(IllegalArgument);\n"
1272 " double fAx=fF1*0.5;\n"
1273 " double fBx=fF1;\n"
1274 " const double fYEps = 1.0E-307;\n"
1275 " const double fXEps = 2.22045e-016;\n"
1276 " double fAy = arg0-GetFInvValue(fF1,fF2,fAx);\n"
1277 " double fBy = arg0-GetFInvValue(fF1,fF2,fBx);\n"
1278 " double fTemp;\n"
1279 " unsigned short nCount;\n"
1280 " for (nCount = 0; nCount < 1000 && !((fAy < 0.0 && fBy > 0.0)"
1281 " || (fAy > 0.0 && fBy < 0.0)); nCount++)\n"
1282 " {\n"
1283 " if (fabs(fAy) <= fabs(fBy))\n"
1284 " {\n"
1285 " fTemp = fAx;\n"
1286 " fAx += 2.0 * (fAx - fBx);\n"
1287 " if (fAx < 0.0)\n"
1288 " fAx = 0.0;\n"
1289 " fBx = fTemp;\n"
1290 " fBy = fAy;\n"
1291 " fAy = arg0-GetFInvValue(fF1,fF2,fAx);\n"
1292 " }\n"
1293 " else\n"
1294 " {\n"
1295 " fTemp = fBx;\n"
1296 " fBx += 2.0 * (fBx - fAx);\n"
1297 " fAx = fTemp;\n"
1298 " fAy = fBy;\n"
1299 " fBy = arg0-GetFInvValue(fF1,fF2,fBx);\n"
1300 " }\n"
1301 " }\n"
1302 " if (fAy == 0.0)\n"
1303 " {\n"
1304 " tmp = fAx;\n"
1305 " return tmp;\n"
1306 " }\n"
1307 " if (fBy == 0.0)\n"
1308 " {\n"
1309 " tmp = fBx;\n"
1310 " return tmp;\n"
1311 " }\n"
1312 " if (!((fAy < 0.0 && fBy > 0.0) || (fAy > 0.0 && fBy < 0.0)))\n"
1313 " return CreateDoubleError(NoConvergence);\n"
1314 " double fPx = fAx;\n"
1315 " double fPy = fAy;\n"
1316 " double fQx = fBx;\n"
1317 " double fQy = fBy;\n"
1318 " double fRx = fAx;\n"
1319 " double fRy = fAy;\n"
1320 " double fSx = 0.5 * (fAx + fBx);\n"
1321 " bool bHasToInterpolate = true;\n"
1322 " nCount = 0;\n"
1323 " while ( nCount < 500 && fabs(fRy) > fYEps &&"
1324 "(fBx-fAx) > fmax( fabs(fAx), fabs(fBx)) * fXEps )\n"
1325 " {\n"
1326 " if (bHasToInterpolate)\n"
1327 " {\n"
1328 " if (fPy!=fQy && fQy!=fRy && fRy!=fPy)\n"
1329 " {\n"
1330 " fSx = fPx * fRy * fQy / (fRy-fPy)"
1331 " / (fQy-fPy)+fRx * fQy * fPy / (fQy-fRy)"
1332 " / (fPy-fRy)+ fQx * fPy * fRy / (fPy-fQy)"
1333 " / (fRy-fQy);\n"
1334 " bHasToInterpolate = (fAx < fSx) && (fSx < fBx);\n"
1335 " }\n"
1336 " else\n"
1337 " bHasToInterpolate = false;\n"
1338 " }\n"
1339 " if(!bHasToInterpolate)\n"
1340 " {\n"
1341 " fSx = 0.5 * (fAx + fBx);\n"
1342 " fPx = fAx; fPy = fAy;\n"
1343 " fQx = fBx; fQy = fBy;\n"
1344 " bHasToInterpolate = true;\n"
1345 " }\n"
1346 " fPx = fQx; fQx = fRx; fRx = fSx;\n"
1347 " fPy = fQy; fQy = fRy;\n"
1348 " fRy = arg0-GetFInvValue(fF1,fF2,fSx);\n"
1349 " if ((fAy < 0.0 && fRy > 0.0) || (fAy > 0.0 && fRy < 0.0))\n"
1350 " {\n"
1351 " fBx = fRx; fBy = fRy;\n"
1352 " }\n"
1353 " else\n"
1354 " {\n"
1355 " fAx = fRx; fAy = fRy;\n"
1356 " }\n"
1357 " bHasToInterpolate = bHasToInterpolate && (fabs(fRy)"
1358 " * 2.0 <= fabs(fQy));\n"
1359 " ++nCount;\n"
1360 " }\n"
1361 " tmp = fRx;\n"
1362 " return tmp;"
1363 "}";
1364}
1365void OpFTest::BinInlineFun(std::set<std::string>& decls,
1366 std::set<std::string>& funs)
1367{
1368 decls.insert(fMachEpsDecl);decls.insert(fMaxGammaArgumentDecl);
1369 decls.insert(lcl_getLanczosSumDecl);decls.insert(GetBetaDecl);
1370 decls.insert(GetLogBetaDecl);decls.insert(GetBetaDistPDFDecl);
1371 decls.insert(lcl_GetBetaHelperContFracDecl);decls.insert(GetBetaDistDecl);
1372 decls.insert(GetFDistDecl);
1373 funs.insert(lcl_getLanczosSum);funs.insert(GetBeta);
1374 funs.insert(GetLogBeta);funs.insert(GetBetaDistPDF);
1375 funs.insert(lcl_GetBetaHelperContFrac);funs.insert(GetBetaDist);
1376 funs.insert(GetFDist);
1377}
1379 const std::string &sSymName, SubArguments &vSubArguments)
1380{
1381 CHECK_PARAMETER_COUNT( 2, 2 );
1382 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
1383 ss << "{\n";
1384 ss << " int gid0 = get_global_id(0);\n";
1385 ss << " double fSum1 = 0.0;\n";
1386 ss << " double fSumSqr1 = 0.0;\n";
1387 ss << " double fSum2 = 0.0;\n";
1388 ss << " double fSumSqr2 = 0.0;\n";
1389 ss << " double fLength1 = 0.0;\n";
1390 ss << " double fLength2 = 0.0;\n";
1391 ss << " double tmp = 0;\n";
1392 GenerateRangeArg( 0, vSubArguments, ss, SkipEmpty,
1393 " fSum1 += arg;\n"
1394 " fSumSqr1 += arg * arg;\n"
1395 " fLength1 += 1;\n"
1396 );
1397 GenerateRangeArg( 1, vSubArguments, ss, SkipEmpty,
1398 " fSum2 += arg;\n"
1399 " fSumSqr2 += arg * arg;\n"
1400 " fLength2 += 1;\n"
1401 );
1402 ss << " if(fLength1 < 2 || fLength2 < 2)\n"
1403 " return CreateDoubleError(NoValue);\n"
1404 " double fS1 = (fSumSqr1-fSum1*fSum1/fLength1)/(fLength1-1.0);\n"
1405 " double fS2 = (fSumSqr2-fSum2*fSum2/fLength2)/(fLength2-1.0);\n"
1406 " if(fS1 == 0 || fS2 == 0)\n"
1407 " return CreateDoubleError(NoValue);\n"
1408 " double fF, fF1, fF2;\n"
1409 " if (fS1 > fS2)\n"
1410 " {\n"
1411 " fF = fS1/fS2;\n"
1412 " fF1 = fLength1-1.0;\n"
1413 " fF2 = fLength2-1.0;\n"
1414 " }\n"
1415 " else\n"
1416 " {\n"
1417 " fF = fS2/fS1;\n"
1418 " fF1 = fLength2-1.0;\n"
1419 " fF2 = fLength1-1.0;\n"
1420 " }\n"
1421 " double fFcdf = GetFDist(fF, fF1, fF2);\n"
1422 " return 2.0*min(fFcdf, 1 - fFcdf);\n";
1423 ss << "}";
1424}
1425void OpB::BinInlineFun(std::set<std::string>& decls,
1426 std::set<std::string>& funs)
1427{
1428 //decls.insert(fBigInvDecl);decls.insert(fLogDblMaxDecl);
1429 decls.insert(GetBinomDistPMFDecl);decls.insert(MinDecl);
1430 decls.insert(fMachEpsDecl);decls.insert(fMaxGammaArgumentDecl);
1431 decls.insert(GetBetaDistDecl);decls.insert(GetBetaDistPDFDecl);
1432 decls.insert(lcl_GetBetaHelperContFracDecl);decls.insert(GetLogBetaDecl);
1433 decls.insert(lcl_getLanczosSumDecl); decls.insert(GetBetaDecl);
1434 funs.insert(GetBinomDistPMF);funs.insert(lcl_GetBinomDistRange);
1435 funs.insert(GetBetaDist);funs.insert(GetBetaDistPDF);
1436 funs.insert(lcl_GetBetaHelperContFrac);funs.insert(GetLogBeta);
1437 funs.insert(lcl_getLanczosSum);funs.insert(GetBeta);
1438}
1439
1441 const std::string &sSymName, SubArguments &vSubArguments)
1442{
1443 CHECK_PARAMETER_COUNT( 4, 4 );
1444 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
1445 ss << "{\n";
1446 ss << " int gid0=get_global_id(0);\n";
1447 ss << " double min = 2.22507e-308;\n";
1448 ss << " double tmp;\n";
1449 GenerateArg( 0, vSubArguments, ss );
1450 GenerateArg( 1, vSubArguments, ss );
1451 GenerateArg( 2, vSubArguments, ss );
1452 GenerateArg( 3, vSubArguments, ss );
1453 ss << " double rxs = floor(arg2);\n"
1454 " double rxe = floor(arg3);\n"
1455 " double rn = floor(arg0);\n"
1456 " double rq = (0.5 - arg1) + 0.5;\n"
1457 " bool bIsValidX = (0.0 <= rxs && rxs <= rxe && rxe <= rn);\n"
1458 " if (bIsValidX && 0.0 < arg1 && arg1 < 1.0)\n"
1459 " {\n"
1460 " if (rxs == rxe)\n"
1461 " tmp = GetBinomDistPMF(rxs, rn, arg1);\n"
1462 " else\n"
1463 " {\n"
1464 " double fFactor = pow(rq, rn);\n"
1465 " if (fFactor > min)\n"
1466 " tmp ="
1467 " lcl_GetBinomDistRange(rn, rxs, rxe, fFactor, arg1, rq);\n"
1468 " else\n"
1469 " {\n"
1470 " fFactor = pow(arg1, rn);\n"
1471 " if (fFactor > min)\n"
1472 " {\n"
1473 " tmp ="
1474 "lcl_GetBinomDistRange(rn, rn - rxe, rn - rxs, fFactor, rq, arg1);\n"
1475 " }\n"
1476 " else\n"
1477 " tmp ="
1478 "GetBetaDist(rq, rn - rxe, rxe + 1.0)"
1479 "- GetBetaDist(rq, rn - rxs + 1, rxs);\n"
1480 " }\n"
1481 " }\n"
1482 " }\n"
1483 " else\n"
1484 " {\n"
1485 " if (bIsValidX)\n"
1486 " {\n"
1487 " if (arg1 == 0.0)\n"
1488 " tmp = (rxs == 0.0 ? 1.0 : 0.0);\n"
1489 " else if (arg1 == 1.0)\n"
1490 " tmp = (rxe == rn ? 1.0 : 0.0);\n"
1491 " else\n"
1492 " {\n"
1493 " tmp = DBL_MIN;\n"
1494 " }\n"
1495 " }\n"
1496 " else\n"
1497 " {\n"
1498 " tmp = DBL_MIN;\n"
1499 " }\n"
1500 " }\n"
1501 " return tmp;"
1502 "}\n";
1503}
1504void OpBetaDist::BinInlineFun(std::set<std::string>& decls,
1505 std::set<std::string>& funs)
1506{
1507 decls.insert(fMachEpsDecl);decls.insert(fMaxGammaArgumentDecl);
1508 decls.insert(GetBetaDistDecl);decls.insert(GetBetaDistPDFDecl);
1509 decls.insert(lcl_GetBetaHelperContFracDecl);decls.insert(GetLogBetaDecl);
1510 decls.insert(GetBetaDecl);decls.insert(lcl_getLanczosSumDecl);
1511 funs.insert(GetBetaDist);funs.insert(GetBetaDistPDF);
1512 funs.insert(lcl_GetBetaHelperContFrac);funs.insert(GetLogBeta);
1513 funs.insert(GetBeta);funs.insert(lcl_getLanczosSum);
1514}
1515void OpPoisson::BinInlineFun(std::set<std::string>& decls,
1516 std::set<std::string>& funs)
1517{
1518 decls.insert(fHalfMachEpsDecl);
1519 funs.insert("");
1520 decls.insert(fMaxGammaArgumentDecl);
1521 funs.insert("");
1522 decls.insert(fBigInvDecl);
1523 funs.insert("");
1524 decls.insert(GetLogGammaDecl);
1525 funs.insert(GetLogGamma);
1526 decls.insert(lcl_GetLogGammaHelperDecl);
1527 funs.insert(lcl_GetLogGammaHelper);
1528 decls.insert(lcl_GetGammaHelperDecl);
1529 funs.insert(lcl_GetGammaHelper);
1530 decls.insert(lcl_getLanczosSumDecl);
1531 funs.insert(lcl_getLanczosSum);
1532 decls.insert(GetUpRegIGammaDecl);
1533 funs.insert(GetUpRegIGamma);
1534 decls.insert(GetGammaContFractionDecl);
1535 funs.insert(GetGammaContFraction);
1536 decls.insert(GetGammaSeriesDecl);
1537 funs.insert(GetGammaSeries);
1538}
1540 outputstream &ss,const std::string &sSymName,
1541 SubArguments &vSubArguments)
1542{
1543 CHECK_PARAMETER_COUNT( 2, 3 );
1544 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
1545 ss << "{\n";
1546 ss << " double tmp;\n";
1547 ss << " int gid0=get_global_id(0);\n";
1548 GenerateArg( "x", 0, vSubArguments, ss );
1549 GenerateArg( "lambda", 1, vSubArguments, ss );
1550 GenerateArgWithDefault( "bCumulative", 2, 1, vSubArguments, ss );
1551 ss << " x = floor(x);\n";
1552 ss << " if (lambda <= 0.0 || x < 0.0)\n";
1553 ss << " return CreateDoubleError(IllegalArgument);\n";
1554 ss << " if (!bCumulative)\n";
1555 ss << " {\n";
1556 ss << " if(lambda == 0.0)\n";
1557 ss << " {\n";
1558 ss << " return 0;\n";
1559 ss << " }\n";
1560 ss << " else\n";
1561 ss << " {\n";
1562 ss << " if (lambda >712)\n";
1563 ss << " {\n";
1564 ss << " tmp = (exp(x*log(lambda)-lambda-GetLogGamma(x+1.0)));\n";
1565 ss << " return tmp;\n";
1566 ss << " }\n";
1567 ss << " else\n";
1568 ss << " {\n";
1569 ss << " double fPoissonVar = 1.0;\n";
1570 ss << " for ( int f = 0; f < x; ++f )\n";
1571 ss << " fPoissonVar *= lambda / ( (double)f + 1.0 );\n";
1572 ss << " tmp = ( fPoissonVar * exp( -lambda ) );\n";
1573 ss << " return tmp;\n";
1574 ss << " }\n";
1575 ss << " }\n";
1576 ss << " } \n";
1577 ss << " else\n";
1578 ss << " {\n";
1579 ss << " if (lambda == 0.0)\n";
1580 ss << " {\n";
1581 ss << " return 1;\n";
1582 ss << " }\n";
1583 ss << " else\n";
1584 ss << " {\n";
1585 ss << " if (lambda > 712 )\n";
1586 ss << " {\n";
1587 ss << " tmp = (GetUpRegIGamma(x+1.0,lambda));\n";
1588 ss << " return tmp;\n";
1589 ss << " }\n";
1590 ss << " else\n";
1591 ss << " {\n";
1592 ss << " if (x >= 936.0)\n";
1593 ss << " {\n";
1594 ss << " return 1;\n";
1595 ss << " }\n";
1596 ss << " else\n";
1597 ss << " {\n";
1598 ss << " double fSummand = exp(-lambda);\n";
1599 ss << " double fSum = fSummand;\n";
1600 ss << " int nEnd = (int) (x + 0.5);\n";
1601 ss << " for (int i = 1; i <= nEnd; i++)\n";
1602 ss << " {\n";
1603 ss << " fSummand = (fSummand*lambda)/((double)i);\n";
1604 ss << " fSum += fSummand;\n";
1605 ss << " }\n";
1606 ss << " tmp = fSum;\n";
1607 ss << " return tmp;\n";
1608 ss << " }\n";
1609 ss << " }\n";
1610 ss << " }\n";
1611 ss << " }\n";
1612 ss << "}\n";
1613}
1614
1616 const std::string &sSymName, SubArguments &vSubArguments)
1617{
1618 CHECK_PARAMETER_COUNT( 3, 6 );
1619 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
1620 ss << "{\n";
1621 ss << " int gid0=get_global_id(0);\n";
1622 ss << " double tmp;\n";
1623 GenerateArg( 0, vSubArguments, ss );
1624 GenerateArg( 1, vSubArguments, ss );
1625 GenerateArg( 2, vSubArguments, ss );
1626 GenerateArgWithDefault( "arg3", 3, 0, vSubArguments, ss );
1627 GenerateArgWithDefault( "arg4", 4, 1, vSubArguments, ss );
1628 GenerateArgWithDefault( "arg5", 5, 1, vSubArguments, ss );
1629 ss << " double fScale = arg4 - arg3;\n"
1630 " if (fScale <= 0.0 || arg1 <= 0.0 || arg2 <= 0.0)\n"
1631 " return CreateDoubleError(IllegalArgument);\n"
1632 " if (arg5)\n"
1633 " {\n"
1634 " if (arg0< arg3)\n"
1635 " {\n"
1636 " tmp = 0.0;\n"
1637 " return tmp;\n"
1638 " }\n"
1639 " if (arg0 > arg4)\n"
1640 " {\n"
1641 " tmp = 1.0;\n"
1642 " return tmp;\n"
1643 " }\n"
1644 " arg0 = (arg0-arg3)/fScale;\n"
1645 " tmp = GetBetaDist(arg0, arg1, arg2);\n"
1646 " }\n"
1647 " else\n"
1648 " {\n"
1649 " if (arg0 < arg3 || arg0 > arg4 )\n"
1650 " {\n"
1651 " tmp = 0.0;\n"
1652 " return tmp;\n"
1653 " }\n"
1654 " arg0 = (arg0 - arg3)/fScale;\n"
1655 " tmp = GetBetaDistPDF(arg0, arg1, arg2)/fScale;\n"
1656 " }\n";
1657 ss << " return tmp;\n";
1658 ss << "}\n";
1659}
1660void OpBetainv::BinInlineFun(std::set<std::string>& decls,
1661 std::set<std::string>& funs)
1662{
1663 decls.insert(fMachEpsDecl);
1664 funs.insert("");
1665 decls.insert(fMaxGammaArgumentDecl);
1666 funs.insert("");
1667 decls.insert(lcl_IterateInverseBetaInvDecl);
1668 funs.insert(lcl_IterateInverseBetaInv);
1669 decls.insert(GetBetaDistDecl);
1670 funs.insert(GetBetaDist);
1671 decls.insert(lcl_HasChangeOfSignDecl);
1672 funs.insert(lcl_HasChangeOfSign);
1673 decls.insert(lcl_HasChangeOfSignDecl);
1674 funs.insert(lcl_HasChangeOfSign);
1675 decls.insert(lcl_HasChangeOfSignDecl);
1676 funs.insert(lcl_HasChangeOfSign);
1677 decls.insert(lcl_GetBetaHelperContFracDecl);
1678 funs.insert(lcl_GetBetaHelperContFrac);
1679 decls.insert(GetBetaDistPDFDecl);
1680 funs.insert(GetBetaDistPDF);
1681 decls.insert(GetLogBetaDecl);
1682 funs.insert(GetLogBeta);
1683 decls.insert(GetBetaDecl);
1684 funs.insert(GetBeta);
1685 decls.insert(lcl_getLanczosSumDecl);
1686 funs.insert(lcl_getLanczosSum);
1687}
1689 outputstream &ss,const std::string &sSymName,
1690 SubArguments &vSubArguments)
1691{
1692 CHECK_PARAMETER_COUNT( 3, 5 );
1693 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
1694 ss << "{\n";
1695 ss << " int gid0=get_global_id(0);\n";
1696 GenerateArg( "tmp0", 0, vSubArguments, ss );
1697 GenerateArg( "tmp1", 1, vSubArguments, ss );
1698 GenerateArg( "tmp2", 2, vSubArguments, ss );
1699 GenerateArgWithDefault( "tmp3", 3, 0, vSubArguments, ss );
1700 GenerateArgWithDefault( "tmp4", 4, 1, vSubArguments, ss );
1701 ss << " if (tmp0 < 0.0 || tmp0 > 1.0 ||";
1702 ss << "tmp3 >= tmp4 || tmp1 <= 0.0 || tmp2 <= 0.0)\n";
1703 ss << " {\n";
1704 ss << " return CreateDoubleError(IllegalArgument);\n";
1705 ss << " }\n";
1706 ss << " bool bConvError;\n";
1707 ss << " double fVal = lcl_IterateInverseBetaInv";
1708 ss << "(tmp0, tmp1, tmp2, 0.0, 1.0, &bConvError);\n";
1709 ss << " if(bConvError)\n";
1710 ss << " return CreateDoubleError(NoConvergence);\n";
1711 ss << " return (tmp3 + fVal*(tmp4 - tmp3));\n";
1712 ss << "}\n";
1713}
1715 const std::string &sSymName, SubArguments& vSubArguments)
1716{
1717 CHECK_PARAMETER_COUNT( 1, 30 );
1718 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
1719 ss << "{\n";
1720 ss << " int gid0 = get_global_id(0);\n";
1721 ss << " double vSum = 0.0;\n";
1722 ss << " double vMean = 0.0;\n";
1723 ss << " int cnt = 0;\n";
1724 GenerateRangeArgs( vSubArguments, ss, SkipEmpty,
1725 " vSum += arg;\n"
1726 " ++cnt;\n"
1727 );
1728 ss << " vMean = vSum / cnt;\n";
1729 ss << " vSum = 0.0;\n";
1730 GenerateRangeArgs( vSubArguments, ss, SkipEmpty,
1731 " vSum += ( arg - vMean ) * ( arg - vMean );\n"
1732 );
1733 ss << " return vSum;\n";
1734 ss << "}\n";
1735}
1736
1738 const std::string &sSymName, SubArguments &vSubArguments)
1739{
1740 CHECK_PARAMETER_COUNT( 4, 5 );
1741 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
1742 ss << "{\n";
1743 ss << " int gid0=get_global_id(0);\n";
1744 GenerateArg( "x", 0, vSubArguments, ss );
1745 GenerateArg( "n", 1, vSubArguments, ss );
1746 GenerateArg( "M", 2, vSubArguments, ss );
1747 GenerateArg( "N", 3, vSubArguments, ss );
1748 GenerateArgWithDefault( "fCumulative", 4, 0, vSubArguments, ss );
1749 ss <<
1750 " x = floor(x);\n"
1751 " n = floor(n);\n"
1752 " M = floor(M);\n"
1753 " N = floor(N);\n"
1754 " double num[9];\n"
1755 " double tmp = 0;\n"
1756 " if( (x < 0.0) || (n < x) || (N < n) ||"
1757 "(N < M) || (M < 0.0) )\n"
1758 " {\n"
1759 " return CreateDoubleError(IllegalArgument);\n"
1760 " }\n"
1761 " for(int i = (fCumulative ? 0 : x); i <= x; ++i )\n"
1762 " {\n"
1763 " if( (M < i) || (i < n - N + M) )\n"
1764 " continue;\n"
1765 " num[0]=M;\n"
1766 " num[1]=i;\n"
1767 " num[2]=M-i;\n"
1768 " num[3]=N-M;\n"
1769 " num[4]=n-i;\n"
1770 " num[5]=N-M-n+i;\n"
1771 " num[6]=N;\n"
1772 " num[7]=n;\n"
1773 " num[8]=N-n;\n"
1774 " for(int i=0;i<9;i++)\n"
1775 " {\n"
1776 " if(num[i]<171)\n"
1777 " {\n"
1778 " if(num[i]==0)\n"
1779 " num[i]=0;\n"
1780 " else\n"
1781 " num[i]=log(tgamma(num[i])*num[i]);\n"
1782 " }\n"
1783 " else\n"
1784 " num[i]=0.5*log(2.0*M_PI)+(num[i]+0.5)*log(num[i])-num[i]+\n"
1785 " (1.0/(12.0*num[i])-1.0/(360*pow(num[i],3)));\n"
1786 " }\n"
1787 " tmp+=pow(M_E,(num[0]+num[3]+num[7]+num[8]-num[1]-num[2]-num[4]-num[5]-num[6]));\n"
1788 " }\n"
1789 " return tmp;\n";
1790 ss << "}\n";
1791}
1792
1794 const std::string &sSymName, SubArguments &vSubArguments)
1795{
1796 CHECK_PARAMETER_COUNT( 1, 30 );
1797 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
1798 ss << "{\n";
1799 ss << " int gid0 = get_global_id(0);\n";
1800 ss << " double sum=0.0;\n";
1801 ss << " double totallength=0;\n";
1802 GenerateRangeArgs( vSubArguments, ss, SkipEmpty,
1803 " sum += arg;\n"
1804 " ++totallength;\n"
1805 );
1806 ss << " double mean = sum / totallength;\n";
1807 ss << " sum = 0.0;\n";
1808 GenerateRangeArgs( vSubArguments, ss, SkipEmpty,
1809 " sum += fabs(arg-mean);\n"
1810 );
1811 ss << " return sum/totallength;\n";
1812 ss << "}";
1813}
1814
1815// Block of functions OpCovar-OpRsq that are rather similar.
1816
1818 const std::string &sSymName, SubArguments& vSubArguments)
1819{
1820 CHECK_PARAMETER_COUNT( 2, 2 );
1823 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
1824 ss << "{\n";
1825 ss << " int gid0 = get_global_id(0);\n";
1826 ss << " double fSumX = 0.0;\n";
1827 ss << " double fSumY = 0.0;\n";
1828 ss << " double fMeanX = 0.0;\n";
1829 ss << " double fMeanY = 0.0;\n";
1830 ss << " double fSumDeltaXDeltaY = 0.0;\n";
1831 ss << " double fCount = 0.0;\n";
1832 GenerateRangeArgPair( 0, 1, vSubArguments, ss, SkipEmpty,
1833 " fSumX += arg1;\n"
1834 " fSumY += arg2;\n"
1835 " fCount += 1.0;\n"
1836 );
1837 ss << " if( fCount < 1 )\n";
1838 ss << " return CreateDoubleError(NoValue);\n";
1839 ss << " fMeanX = fSumX / fCount;\n";
1840 ss << " fMeanY = fSumY / fCount;\n";
1841 GenerateRangeArgPair( 0, 1, vSubArguments, ss, SkipEmpty,
1842 " fSumDeltaXDeltaY +=(arg1-fMeanX)*(arg2-fMeanY);\n"
1843 );
1844 ss << " return fSumDeltaXDeltaY / fCount;\n";
1845 ss << "}\n";
1846}
1847
1849 const std::string &sSymName, SubArguments &vSubArguments)
1850{
1851 CHECK_PARAMETER_COUNT( 3, 3 );
1854 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
1855 ss << "{\n";
1856 ss << " int gid0 = get_global_id(0);\n";
1857 ss << " double fSumX = 0.0;\n";
1858 ss << " double fSumY = 0.0;\n";
1859 ss << " double fMeanX = 0.0;\n";
1860 ss << " double fMeanY = 0.0;\n";
1861 ss << " double fSumDeltaXDeltaY = 0.0;\n";
1862 ss << " double fSumSqrDeltaX = 0.0;\n";
1863 ss << " double fCount = 0.0;\n";
1864 GenerateArg( "arg0", 0, vSubArguments, ss );
1865 GenerateRangeArgPair( 1, 2, vSubArguments, ss, SkipEmpty,
1866 // note that arg1 -> Y, arg2 -> X
1867 " fSumX += arg2;\n"
1868 " fSumY += arg1;\n"
1869 " fCount += 1.0;\n"
1870 );
1871 ss << " if( fCount < 1 )\n";
1872 ss << " return CreateDoubleError(NoValue);\n";
1873 ss << " fMeanX = fSumX / fCount;\n";
1874 ss << " fMeanY = fSumY / fCount;\n";
1875 GenerateRangeArgPair( 1, 2, vSubArguments, ss, SkipEmpty,
1876 " fSumDeltaXDeltaY +=(arg2-fMeanX)*(arg1-fMeanY);\n"
1877 " fSumSqrDeltaX += (arg2-fMeanX)*(arg2-fMeanX);\n"
1878 );
1879 ss << " if(fSumSqrDeltaX == 0.0)\n";
1880 ss << " return CreateDoubleError(DivisionByZero);\n";
1881 ss << " return fMeanY + fSumDeltaXDeltaY / fSumSqrDeltaX * (arg0 - fMeanX);\n";
1882 ss << "}\n";
1883}
1884
1885void OpInterceptSlopeBase::GenerateCode( outputstream &ss, const std::string &sSymName,
1886 SubArguments &vSubArguments, const char* finalComputeCode )
1887{
1888 CHECK_PARAMETER_COUNT( 2, 2 );
1891 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
1892 ss << "{\n";
1893 ss << " int gid0 = get_global_id(0);\n";
1894 ss << " double fSumX = 0.0;\n";
1895 ss << " double fSumY = 0.0;\n";
1896 ss << " double fMeanX = 0.0;\n";
1897 ss << " double fMeanY = 0.0;\n";
1898 ss << " double fSumDeltaXDeltaY = 0.0;\n";
1899 ss << " double fSumSqrDeltaX = 0.0;\n";
1900 ss << " double fCount = 0.0;\n";
1901 GenerateRangeArgPair( 0, 1, vSubArguments, ss, SkipEmpty,
1902 // note that arg1 -> Y, arg2 -> X
1903 " fSumX += arg2;\n"
1904 " fSumY += arg1;\n"
1905 " fCount += 1.0;\n"
1906 );
1907 ss << " if( fCount < 1 )\n";
1908 ss << " return CreateDoubleError(NoValue);\n";
1909 ss << " fMeanX = fSumX / fCount;\n";
1910 ss << " fMeanY = fSumY / fCount;\n";
1911 GenerateRangeArgPair( 0, 1, vSubArguments, ss, SkipEmpty,
1912 " fSumDeltaXDeltaY +=(arg2-fMeanX)*(arg1-fMeanY);\n"
1913 " fSumSqrDeltaX += (arg2-fMeanX)*(arg2-fMeanX);\n"
1914 );
1915 ss << finalComputeCode;
1916 ss << "}\n";
1917}
1918
1920 const std::string &sSymName, SubArguments &vSubArguments)
1921{
1922 GenerateCode( ss, sSymName, vSubArguments,
1923 " if(fSumSqrDeltaX == 0.0)\n"
1924 " return CreateDoubleError(DivisionByZero);\n"
1925 " return fMeanY - (fSumDeltaXDeltaY/fSumSqrDeltaX)*fMeanX;\n"
1926 );
1927}
1928
1930 const std::string &sSymName, SubArguments &vSubArguments)
1931{
1932 GenerateCode( ss, sSymName, vSubArguments,
1933 " if(fSumSqrDeltaX == 0.0)\n"
1934 " return CreateDoubleError(DivisionByZero);\n"
1935 " return fSumDeltaXDeltaY / fSumSqrDeltaX;\n"
1936 );
1937}
1938
1939void OpPearsonCovarBase::GenerateCode( outputstream &ss, const std::string &sSymName,
1940 SubArguments &vSubArguments, double minimalCountValue, const char* finalComputeCode )
1941{
1942 CHECK_PARAMETER_COUNT( 2, 2 );
1945 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
1946 ss << "{\n";
1947 ss << " int gid0 = get_global_id(0);\n";
1948 ss << " double fSumX = 0.0;\n";
1949 ss << " double fSumY = 0.0;\n";
1950 ss << " double fMeanX = 0.0;\n";
1951 ss << " double fMeanY = 0.0;\n";
1952 ss << " double fSumDeltaXDeltaY = 0.0;\n";
1953 ss << " double fSumSqrDeltaX = 0.0;\n";
1954 ss << " double fSumSqrDeltaY = 0.0;\n";
1955 ss << " double fCount = 0.0;\n";
1956 GenerateRangeArgPair( 0, 1, vSubArguments, ss, SkipEmpty,
1957 // note that arg1 -> Y, arg2 -> X
1958 " fSumX += arg2;\n"
1959 " fSumY += arg1;\n"
1960 " fCount += 1.0;\n"
1961 );
1962 ss << " if( fCount < " << minimalCountValue <<" )\n";
1963 ss << " return CreateDoubleError(NoValue);\n";
1964 ss << " fMeanX = fSumX / fCount;\n";
1965 ss << " fMeanY = fSumY / fCount;\n";
1966 GenerateRangeArgPair( 0, 1, vSubArguments, ss, SkipEmpty,
1967 " fSumDeltaXDeltaY +=(arg2-fMeanX)*(arg1-fMeanY);\n"
1968 " fSumSqrDeltaX += (arg2-fMeanX)*(arg2-fMeanX);\n"
1969 " fSumSqrDeltaY += (arg1-fMeanY)*(arg1-fMeanY);\n"
1970 );
1971 ss << finalComputeCode;
1972 ss << "}\n";
1973}
1974
1976 outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments)
1977{
1978 GenerateCode( ss, sSymName, vSubArguments, 1,
1979 " if (fSumSqrDeltaX == 0 || fSumSqrDeltaY == 0)\n"
1980 " return CreateDoubleError(DivisionByZero);\n"
1981 " return ( fSumDeltaXDeltaY / sqrt( fSumSqrDeltaX * fSumSqrDeltaY));\n"
1982 );
1983}
1984
1986 const std::string &sSymName, SubArguments &vSubArguments)
1987{
1988 GenerateCode( ss, sSymName, vSubArguments, 3,
1989 " if(fSumSqrDeltaX == 0.0)\n"
1990 " return CreateDoubleError(DivisionByZero);\n"
1991 " return sqrt((fSumSqrDeltaY - fSumDeltaXDeltaY * \n"
1992 " fSumDeltaXDeltaY / fSumSqrDeltaX)\n"
1993 " /(fCount - 2.0));\n"
1994 );
1995}
1996
1998 outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments)
1999{
2000 // This is basically pow(OpPearson,2)
2001 GenerateCode( ss, sSymName, vSubArguments, 1,
2002 " if (fSumSqrDeltaX == 0 || fSumSqrDeltaY == 0)\n"
2003 " return CreateDoubleError(DivisionByZero);\n"
2004 " return ( fSumDeltaXDeltaY * fSumDeltaXDeltaY / (fSumSqrDeltaX * fSumSqrDeltaY));\n"
2005 );
2006}
2007
2008void OpVarStDevBase::BinInlineFun(std::set<std::string>& decls,std::set<std::string>& funs)
2009{
2010 decls.insert(is_representable_integerDecl);
2011 funs.insert(is_representable_integer);
2012 decls.insert(approx_equalDecl);
2013 funs.insert(approx_equal);
2014 decls.insert(fsub_approxDecl);
2015 funs.insert(fsub_approx);
2016}
2017
2018void OpVarStDevBase::GenerateCode(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments)
2019{
2020 CHECK_PARAMETER_COUNT( 1, 30 );
2021 GenerateFunctionDeclaration( sSymName, vSubArguments, ss );
2022 ss << "{\n"; // this must be matched by whoever calls this function
2023 ss << " int gid0 = get_global_id(0);\n";
2024 ss << " double fSum = 0.0;\n";
2025 ss << " double fCount = 0.0;\n";
2026 GenerateRangeArgs( vSubArguments, ss, SkipEmpty,
2027 " fSum += arg;\n"
2028 " fCount += 1.0;\n"
2029 );
2030 ss << " if (fCount == 0)\n";
2031 ss << " return CreateDoubleError(DivisionByZero);\n";
2032 ss << " double fMean = fSum / fCount;\n";
2033 ss << " double vSum = 0.0;\n";
2034 GenerateRangeArgs( vSubArguments, ss, SkipEmpty,
2035 " vSum += pown( fsub_approx(arg, fMean), 2 );\n"
2036 );
2037}
2038
2040 const std::string &sSymName, SubArguments &vSubArguments)
2041{
2042 GenerateCode( ss, sSymName, vSubArguments );
2043 ss << " if (fCount <= 1.0)\n";
2044 ss << " return CreateDoubleError(DivisionByZero);\n";
2045 ss << " else\n";
2046 ss << " return vSum / (fCount - 1.0);\n";
2047 ss << "}\n";
2048}
2049
2051 const std::string &sSymName, SubArguments &vSubArguments)
2052{
2053 GenerateCode( ss, sSymName, vSubArguments );
2054 ss << " if (fCount == 0.0)\n";
2055 ss << " return CreateDoubleError(DivisionByZero);\n";
2056 ss << " else\n";
2057 ss << " return vSum / fCount;\n";
2058 ss << "}\n";
2059}
2060
2062 const std::string &sSymName, SubArguments &vSubArguments)
2063{
2064 GenerateCode( ss, sSymName, vSubArguments );
2065 ss << " if (fCount <= 1.0)\n";
2066 ss << " return CreateDoubleError(DivisionByZero);\n";
2067 ss << " else\n";
2068 ss << " return sqrt(vSum / (fCount - 1.0));\n";
2069 ss << "}\n";
2070}
2071
2073 const std::string &sSymName, SubArguments &vSubArguments)
2074{
2075 GenerateCode( ss, sSymName, vSubArguments );
2076 ss << " if (fCount <= 0.0)\n";
2077 ss << " return CreateDoubleError(DivisionByZero);\n";
2078 ss << " else\n";
2079 ss << " return sqrt(vSum / fCount);\n";
2080 ss << "}\n";
2081}
2082
2084 const std::string &sSymName, SubArguments &vSubArguments)
2085{
2086 GenerateCode( ss, sSymName, vSubArguments );
2087 ss << " if(fCount <= 2.0)\n";
2088 ss << " return CreateDoubleError(DivisionByZero);\n";
2089 ss << " double fStdDev = sqrt(vSum / (fCount - 1.0));\n";
2090 ss << " double dx = 0.0;\n";
2091 ss << " double xcube = 0.0;\n";
2092 ss << " if(fStdDev == 0.0)\n";
2093 ss << " return CreateDoubleError(IllegalArgument);\n";
2094 GenerateRangeArgs( vSubArguments, ss, SkipEmpty,
2095 " dx = fsub_approx(arg, fMean) / fStdDev;\n"
2096 " xcube = xcube + dx * dx * dx;\n"
2097 );
2098 ss << " return ((xcube * fCount) / (fCount - 1.0)) / (fCount - 2.0);\n";
2099 ss << "}\n";
2100}
2101
2103 const std::string &sSymName, SubArguments &vSubArguments)
2104{
2105 GenerateCode( ss, sSymName, vSubArguments );
2106 ss << " if(fCount <= 2.0)\n";
2107 ss << " return CreateDoubleError(DivisionByZero);\n";
2108 ss << " double fStdDev = sqrt(vSum / fCount);\n";
2109 ss << " double dx = 0.0;\n";
2110 ss << " double xcube = 0.0;\n";
2111 ss << " if(fStdDev == 0.0)\n";
2112 ss << " return CreateDoubleError(IllegalArgument);\n";
2113 GenerateRangeArgs( vSubArguments, ss, SkipEmpty,
2114 " dx = fsub_approx(arg, fMean) / fStdDev;\n"
2115 " xcube = xcube + dx * dx * dx;\n"
2116 );
2117 ss << " return xcube / fCount;\n";
2118 ss << "}\n";
2119}
2120
2122 const std::string &sSymName, SubArguments &vSubArguments)
2123{
2124 GenerateCode( ss, sSymName, vSubArguments );
2125 ss << " if( fCount < 4 )\n";
2126 ss << " return CreateDoubleError(DivisionByZero);\n";
2127 ss << " double fStdDev = sqrt(vSum / (fCount - 1.0));\n";
2128 ss << " double dx = 0.0;\n";
2129 ss << " double xpower4 = 0.0;\n";
2130 GenerateRangeArgs( vSubArguments, ss, SkipEmpty,
2131 " dx = (arg -fMean) / fStdDev;\n"
2132 " xpower4 = xpower4 + (dx * dx * dx * dx);\n"
2133 );
2134 ss<< " double k_d = (fCount - 2.0) * (fCount - 3.0);\n";
2135 ss<< " double k_l = fCount * (fCount + 1.0) / ((fCount - 1.0) * k_d);\n";
2136 ss<< " double k_t = 3.0 * (fCount - 1.0) * (fCount - 1.0) / k_d;\n";
2137 ss<< " return xpower4 * k_l - k_t;\n";
2138 ss << "}";
2139}
2140
2141void OpMin::BinInlineFun(std::set<std::string>& decls,std::set<std::string>& funs)
2142{
2143 decls.insert(fmin_countDecl);
2144 funs.insert(fmin_count);
2145}
2146
2147void OpMax::BinInlineFun(std::set<std::string>& decls,std::set<std::string>& funs)
2148{
2149 decls.insert(fmax_countDecl);
2150 funs.insert(fmax_count);
2151}
2152
2153void OpAverage::BinInlineFun(std::set<std::string>& decls,std::set<std::string>& funs)
2154{
2155 decls.insert(fsum_countDecl);
2156 funs.insert(fsum_count);
2157}
2158
2159}
2160
2161/* vim:set shiftwidth=4 softtabstop=4 expandtab: */
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void BinInlineFun(std::set< std::string > &decls, std::set< std::string > &funs) override
void BinInlineFun(std::set< std::string > &decls, std::set< std::string > &funs) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
void BinInlineFun(std::set< std::string > &decls, std::set< std::string > &funs) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void BinInlineFun(std::set< std::string > &, std::set< std::string > &) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void BinInlineFun(std::set< std::string > &, std::set< std::string > &) override
virtual void BinInlineFun(std::set< std::string > &, std::set< std::string > &) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void BinInlineFun(std::set< std::string > &, std::set< std::string > &) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void BinInlineFun(std::set< std::string > &, std::set< std::string > &) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void BinInlineFun(std::set< std::string > &, std::set< std::string > &) override
virtual std::string GetBottom() override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void BinInlineFun(std::set< std::string > &, std::set< std::string > &) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void BinInlineFun(std::set< std::string > &, std::set< std::string > &) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual std::string GetBottom() override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
void BinInlineFun(std::set< std::string > &decls, std::set< std::string > &funs) override
void BinInlineFun(std::set< std::string > &decls, std::set< std::string > &funs) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void BinInlineFun(std::set< std::string > &, std::set< std::string > &) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
void BinInlineFun(std::set< std::string > &decls, std::set< std::string > &funs) override
void BinInlineFun(std::set< std::string > &decls, std::set< std::string > &funs) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void BinInlineFun(std::set< std::string > &, std::set< std::string > &) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
void GenerateCode(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments, const char *finalComputeCode)
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void BinInlineFun(std::set< std::string > &, std::set< std::string > &) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void BinInlineFun(std::set< std::string > &decls, std::set< std::string > &funs) override
virtual void BinInlineFun(std::set< std::string > &decls, std::set< std::string > &funs) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void BinInlineFun(std::set< std::string > &, std::set< std::string > &) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void BinInlineFun(std::set< std::string > &, std::set< std::string > &) override
void GenerateCode(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments, double minimalCountValue, const char *finalComputeCode)
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void BinInlineFun(std::set< std::string > &, std::set< std::string > &) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void BinInlineFun(std::set< std::string > &, std::set< std::string > &) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void BinInlineFun(std::set< std::string > &, std::set< std::string > &) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void BinInlineFun(std::set< std::string > &, std::set< std::string > &) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
void GenerateCode(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments)
virtual void BinInlineFun(std::set< std::string > &decls, std::set< std::string > &funs) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void GenSlidingWindowFunction(outputstream &ss, const std::string &sSymName, SubArguments &vSubArguments) override
virtual void BinInlineFun(std::set< std::string > &, std::set< std::string > &) override
void GenerateArg(const char *name, int arg, SubArguments &vSubArguments, outputstream &ss, EmptyArgType empty=EmptyIsZero, GenerateArgTypeType generateType=DoNotGenerateArgType)
Definition: opbase.cxx:229
static void GenerateRangeArg(int arg, SubArguments &vSubArguments, outputstream &ss, EmptyArgType empty, const char *code)
Definition: opbase.cxx:417
void GenerateArgWithDefault(const char *name, int arg, double def, SubArguments &vSubArguments, outputstream &ss, EmptyArgType empty=EmptyIsZero)
Definition: opbase.cxx:307
void GenerateFunctionDeclaration(const std::string &sSymName, SubArguments &vSubArguments, outputstream &ss)
Definition: opbase.cxx:566
static void GenerateRangeArgPair(int arg1, int arg2, SubArguments &vSubArguments, outputstream &ss, EmptyArgType empty, const char *code, const char *firstElementDiff=nullptr)
Definition: opbase.cxx:423
static void GenerateRangeArgs(int firstArg, int lastArg, SubArguments &vSubArguments, outputstream &ss, EmptyArgType empty, const char *code)
Definition: opbase.cxx:316
std::vector< DynamicKernelArgumentRef > SubArguments
Definition: opbase.hxx:347
const char is_representable_integerDecl[]
const char fsub_approx[]
const char fsub_approxDecl[]
const char approx_equalDecl[]
const char approx_equal[]
const char is_representable_integer[]
const char lcl_GetGammaHelperDecl[]
const char GetGammaContFraction[]
const char GetFInvValue[]
const char lcl_GetLogGammaHelperDecl[]
const char GetGammaInvValueDecl[]
const char GetLogGamma[]
const char fMachEpsDecl[]
const char GetLogBetaDecl[]
const char lcl_IterateInverseBetaInv[]
const char GetBetaDecl[]
const char GetChiSqDistPDF[]
const char lcl_IterateInverseDecl[]
const char taylor[]
const char lcl_IterateInverseBetaInvDecl[]
const char lcl_GetGammaHelper[]
const char taylorDecl[]
const char GetLogBeta[]
const char GetBetaDistPDFDecl[]
const char GetUpRegIGamma[]
const char fMaxGammaArgumentDecl[]
const char gaussinvDecl[]
const char fsum_countDecl[]
const char GetFInvValueDecl[]
const char lcl_GetBetaHelperContFracDecl[]
const char GetBetaDistDecl[]
const char GetGammaInvValue[]
const char gaussinv[]
const char GetLogGammaDecl[]
const char GetBinomDistPMF[]
const char GetLowRegIGamma[]
const char gaussDecl[]
const char lcl_IterateInverseChiInv[]
const char MinDecl[]
const char GetTDist[]
const char GetChiDistDecl[]
const char lcl_IterateInverseChiSQInvDecl[]
const char fmax_count[]
const char GetValueDecl[]
const char GetBeta[]
const char GetChiSqDistPDFDecl[]
const char lcl_IterateInverse[]
const char lcl_GetBinomDistRange[]
const char phi[]
const char lcl_IterateInverseChiInvDecl[]
const char gauss[]
const char GetLowRegIGammaDecl[]
const char lcl_getLanczosSumDecl[]
const char fsum_count[]
const char GetGammaSeriesDecl[]
const char fmin_countDecl[]
const char GetGammaDistPDFDecl[]
const char lcl_HasChangeOfSignDecl[]
const char lcl_GetLogGammaHelper[]
const char GetFDistDecl[]
const char lcl_GetBinomDistRangeDecl[]
const char GetValue[]
const char fHalfMachEpsDecl[]
const char lcl_HasChangeOfSign[]
const char GetChiSqDistCDFDecl[]
const char fLogDblMaxDecl[]
const char GetTDistDecl[]
const char GetGammaContFractionDecl[]
const char lcl_IterateInverseChiSQInv[]
const char GetBetaDist[]
const char fmin_count[]
const char lcl_getLanczosSum[]
const char GetGammaDistPDF[]
const char GetChiSqDistCDF[]
const char GetGammaSeries[]
const char fBigInvDecl[]
const char GetGammaDistDecl[]
const char GetChiDist[]
const char GetUpRegIGammaDecl[]
const char GetBinomDistPMFDecl[]
const char fmax_countDecl[]
const char phiDecl[]
const char GetGammaDist[]
const char GetBetaDistPDF[]
const char GetFDist[]
const char lcl_GetBetaHelperContFrac[]
#define CHECK_PARAMETER_DOUBLEVECTORREF(arg)
Definition: opbase.hxx:97
#define CHECK_PARAMETER_COUNT(min, max)
Definition: opbase.hxx:85