randlib plugged.
[scilab.git] / scilab / modules / randlib / sci_gateway / cpp / sci_grand.cpp
1 /*
2 * Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
3 * Copyright (C) 2011 - DIGITEO - Cedric DELAMARRE
4 *
5 * This file must be used under the terms of the CeCILL.
6 * This source file is licensed as described in the file COPYING, which
7 * you should have received as part of this distribution.  The terms
8 * are also available at
9 * http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt
10 *
11 */
12 /*--------------------------------------------------------------------------*/
13 #include "randlib_gw.hxx"
14 #include "function.hxx"
15 #include "double.hxx"
16 #include "string.hxx"
17 #include "configvariable.hxx"
18
19 extern "C"
20 {
21 #include "localization.h"
22 #include "Scierror.h"
23 #include "sciprint.h"
24 #include "grand.h"
25 #include "clcg4.h"
26 #include "others_generators.h"
27 }
28 /*--------------------------------------------------------------------------*/
29
30 types::Function::ReturnValue sci_grand(types::typed_list &in, int _iRetCount, types::typed_list &out)
31 {
32     enum {MT, KISS, CLCG4, CLCG2, URAND, FSULTRA};
33
34     //  names at the scilab level
35     wchar_t* names_gen[6] = {L"mt", L"kiss", L"clcg4", L"clcg2", L"urand", L"fsultra"};
36
37     types::String* pStrMethod   = NULL;
38     types::String* pStrGenOrPhr = NULL;
39
40     std::vector<types::Double*> vectpDblInput;
41
42     int iStrPos     = 0;
43     int iPos        = 0;
44     int meth        = -1;
45     int iRows       = -1;
46     int iCols       = -1;
47     int iNumIter    = -1;
48
49     int current_base_gen = ConfigVariable::getCurrentBaseGen();
50
51 // *** check the maximal number of input args. ***
52     if(in.size() > 6)
53     {
54         ScierrorW(77, _W("%ls: Wrong number of input argument(s): %d to %d expected.\n"), L"grand", 1, 6);
55         return types::Function::Error;
56     }
57
58 // *** check number of output args. ***
59     if(_iRetCount > 1)
60     {
61         ScierrorW(78, _W("%ls: Wrong number of output argument(s): %d expected.\n"), L"grand", 1);
62         return types::Function::Error;
63     }
64
65 // *** find the mothod string. ***
66     for(int i = 0; i < in.size(); i++)
67     {
68         if(in[i]->isString())
69         {
70             pStrMethod = in[i]->getAs<types::String>();
71             iStrPos = i;
72             break;
73         }
74     }
75
76     if(pStrMethod == NULL)
77     {
78         ScierrorW(78, _W("%ls: Wrong number of output argument(s): At least %d string expected.\n"), L"grand", 1);
79         return types::Function::Error;
80     }
81
82     wchar_t* wcsMeth = pStrMethod->get(0);
83     int iNumInputArg = 5;
84     if(wcscmp(wcsMeth, L"bet") == 0) // beta
85     {
86         meth = 1;
87     }
88     else if(wcscmp(wcsMeth,L"bin") == 0) // binomial
89     {
90         meth = 2;
91     }
92     else if(wcscmp(wcsMeth, L"nbn") == 0) // negative binomial
93     {
94         meth = 3;
95     }
96     else if(wcscmp(wcsMeth, L"chi") == 0) // chisquare
97     {
98         iNumInputArg = 4;
99         meth = 4;
100     }
101     else if(wcscmp(wcsMeth, L"nch") == 0) // non central chisquare
102     {
103         meth = 5;
104     }
105     else if(wcscmp(wcsMeth, L"exp") == 0) // exponential
106     {
107         iNumInputArg = 4;
108         meth = 6;
109     }
110     else if(wcscmp(wcsMeth, L"f") == 0) // F variance ratio
111     {
112         meth = 7;
113     }
114     else if(wcscmp(wcsMeth, L"nf") == 0) // non central F variance ratio
115     {
116         iNumInputArg = 6;
117         meth = 8;
118     }
119     else if(wcscmp(wcsMeth, L"gam") == 0) // gamma
120     {
121         meth = 9;
122     }
123     else if(wcscmp(wcsMeth, L"nor") == 0) // Gauss Laplace (normal)
124     {
125         meth = 10;
126     }
127     else if(wcscmp(wcsMeth, L"mn") == 0) // multivariate gaussian (multivariate normal)
128     {
129         iNumInputArg = 4;
130         meth = 11;
131     }
132     else if(wcscmp(wcsMeth, L"geom") == 0) // geometric
133     {
134         iNumInputArg = 4;
135         meth = 12;
136     }
137     else if(wcscmp(wcsMeth, L"markov") == 0) // markov
138     {
139         iNumInputArg = 4;
140         meth = 13;
141     }
142     else if(wcscmp(wcsMeth, L"mul") == 0) // multinomial
143     {
144         iNumInputArg = 4;
145         meth = 14;
146     }
147     else if(wcscmp(wcsMeth, L"poi") == 0) // Poisson
148     {
149         iNumInputArg = 4;
150         meth = 15;
151     }
152     else if(wcscmp(wcsMeth, L"prm") == 0) // random permutations
153     {
154         iNumInputArg = 3;
155         meth = 16;
156     }
157     else if(wcscmp(wcsMeth, L"def") == 0) // uniform (def)
158     {
159         iNumInputArg = 3;
160         meth = 17;
161     }
162     else if(wcscmp(wcsMeth, L"unf") == 0) // uniform (unf)
163     {
164         meth = 18;
165     }
166     else if(wcscmp(wcsMeth, L"uin") == 0) // uniform (uin)
167     {
168         meth = 19;
169     }
170     else if(wcscmp(wcsMeth, L"lgi") == 0) // uniform (lgi)
171     {
172         iNumInputArg = 3;
173         meth = 20;
174     }
175     else if(wcscmp(wcsMeth, L"getgen") == 0) // getgen
176     {
177         iNumInputArg = 1;
178         meth = 21;
179     }
180     else if(wcscmp(wcsMeth, L"setgen") == 0) // setgen
181     {
182         iNumInputArg = 2;
183         meth = 22;
184     }
185     else if(wcscmp(wcsMeth, L"getsd") == 0) // getsd
186     {
187         iNumInputArg = 1;
188         meth = 23;
189     }
190     else if(wcscmp(wcsMeth, L"setsd") == 0) // setsd
191     {
192         switch(current_base_gen)
193         {
194             case MT:        iNumInputArg = 2; break;
195             case KISS:      iNumInputArg = 5; break;
196             case CLCG2:     iNumInputArg = 3; break;
197             case CLCG4:     iNumInputArg = 5; break;
198             case URAND:     iNumInputArg = 2; break;
199             case FSULTRA:
200             {
201                 if(in.size() != 2 && in.size() != 3)
202                 {
203                     ScierrorW(77, _W("%ls: Wrong number of input argument(s) for method %ls: %d or %d expected.\n"), L"grand", wcsMeth, 2, 3);
204                     return types::Function::Error;
205                 }
206
207                 iNumInputArg = in.size();
208                 break;
209             }
210         }
211         meth = 24;
212     }
213     else if(wcscmp(wcsMeth, L"phr2sd") == 0) // phr2sd
214     {
215         iNumInputArg = 2;
216         meth = 25;
217     }
218     else if(wcscmp(wcsMeth, L"setcgn") == 0) // setcgn
219     {
220         iNumInputArg = 2;
221         meth = 26;
222     }
223     else if(wcscmp(wcsMeth, L"getcgn") == 0) // getcgn
224     {
225         iNumInputArg = 1;
226         meth = 27;
227     }
228     else if(wcscmp(wcsMeth, L"initgn") == 0) // initgn
229     {
230         iNumInputArg = 2;
231         meth = 28;
232     }
233     else if(wcscmp(wcsMeth, L"setall") == 0) // setall
234     {
235         iNumInputArg = 5;
236         meth = 29;
237     }
238     else if(wcscmp(wcsMeth, L"advnst") == 0) // advnst
239     {
240         iNumInputArg = 2;
241         meth = 30;
242     }
243     else
244     {
245         ScierrorW(78, _W("%ls: Wrong method for input argument #%d: %ls unknown.\n"), L"grand", iStrPos+1, wcsMeth);
246         return types::Function::Error;
247     }
248
249 // *** get arguments before methode string. ***
250     if(meth == 11 || meth == 13 || meth == 14 || meth == 16) // grand(n, "...", ...
251     {
252         if(iStrPos != 1)
253         {
254             ScierrorW(999, _W("%ls: Wrong position for input argument #%d : Must be in position %d.\n"), L"grand", iStrPos+1, 2);
255             return types::Function::Error;
256         }
257
258         if(in[iPos]->isDouble() == false)
259         {
260             ScierrorW(999, _W("%ls: Wrong type for input argument #%d : A matrix expected.\n"), L"grand", iPos+1);
261             return types::Function::Error;
262         }
263
264         types::Double* pDblTemp = in[iPos]->getAs<types::Double>();
265
266         if(pDblTemp->isScalar() == false)
267         {
268             ScierrorW(999, _W("%ls: Wrong type for input argument #%d : A scalar expected.\n"), L"grand", iPos+1);
269             return types::Function::Error;
270         }
271
272         iNumIter = pDblTemp->get(0);
273         iPos++;
274     }
275     else if(meth < 21) // grand(m, n, "...", ... || grand(matrix, "...", ...
276     {
277         if(iStrPos != 1 && iStrPos != 2)
278         {
279             ScierrorW(999, _W("%ls: Wrong position for input argument #%d : Must be in position %d or %d.\n"), L"grand", iStrPos+1, 2, 3);
280             return types::Function::Error;
281         }
282
283         std::vector<types::Double*> vectpDblTemp;
284         for(int i = 0; i < iStrPos; i++)
285         {
286             if(in[iPos]->isDouble() == false)
287             {
288                 ScierrorW(999, _W("%ls: Wrong type for input argument #%d : A scalar expected.\n"), L"grand", iPos+1);
289                 return types::Function::Error;
290             }
291
292             vectpDblTemp.push_back(in[iPos]->getAs<types::Double>());
293
294             if(iStrPos == 3 && vectpDblTemp[i]->isScalar() == false)
295             {
296                 ScierrorW(999, _W("%ls: Wrong type for input argument #%d : A scalar expected.\n"), L"grand", iPos+1);
297                 return types::Function::Error;
298             }
299             iPos++;
300         }
301
302         if(iStrPos == 2)
303         {
304             iRows = static_cast<int>(vectpDblTemp[0]->get(0));
305             iCols = static_cast<int>(vectpDblTemp[1]->get(0));;
306         }
307         else // iStrPos == 2
308         {
309             iRows = vectpDblTemp[0]->getRows();
310             iCols = vectpDblTemp[0]->getCols();
311             iNumInputArg--;
312         }
313     }
314
315     if(in[iPos]->isString() == false)
316     {
317         ScierrorW(999, _W("%ls: Wrong type for input argument #%d : A string expected.\n"), L"grand", iPos+1);
318         return types::Function::Error;
319     }
320     iPos++; // method string has been already got.
321
322 // *** check the number of input args according the methode. ***
323     if(in.size() != iNumInputArg)
324     {
325         ScierrorW(77, _W("%ls: Wrong number of input argument(s) for method %ls: %d expected.\n"), L"grand", wcsMeth, iNumInputArg);
326         return types::Function::Error;
327     }
328
329 // *** get arguments after methode string. ***
330     if(meth == 22 || meth == 25) // setgen || phr2sd
331     {
332         if(in[iPos]->isString() == false)
333         {
334             ScierrorW(999, _W("%ls: Wrong type for input argument #%d : A string expected.\n"), L"grand", iPos+1);
335             return types::Function::Error;
336         }
337
338         pStrGenOrPhr = in[iPos]->getAs<types::String>();
339
340         if(pStrGenOrPhr->isScalar() == false)
341         {
342             ScierrorW(999, _W("%ls: Wrong size for input argument #%d : Only one string expected.\n"), L"grand", iPos+1);
343             return types::Function::Error;
344         }
345     }
346     else
347     {
348         for(int i = iPos; i < in.size(); i++)
349         {
350             if(in[i]->isDouble() == false)
351             {
352                 ScierrorW(999, _W("%ls: Wrong type for input argument #%d : A scalar expected.\n"), L"grand", iPos+1);
353                 return types::Function::Error;
354             }
355
356             vectpDblInput.push_back(in[i]->getAs<types::Double>());
357         }
358     }
359
360 // *** perform operation in according method string and return result. ***
361
362     switch(meth)
363     {
364         case 1: // beta
365         {
366             types::Double* pDblOut = new types::Double(iRows, iCols);
367             double minlog   = 1.e-37;
368
369             for(int i = 0; i < 2; i++)
370             {
371                 if(vectpDblInput[i]->isScalar() == false)
372                 {
373                     ScierrorW(999, _W("%ls: Wrong type for input argument #%d : A scalar expected.\n"), L"grand", i+4);
374                     return types::Function::Error;
375                 }
376
377                 if(vectpDblInput[i]->get(0) < minlog)
378                 {
379                     ScierrorW(999, _W("%ls: Wrong value for input argument #%d : At least %lf expected.\n"), L"grand", iPos+1, minlog);
380                     return types::Function::Error;
381                 }
382             }
383
384             for(int i = 0; i < pDblOut->getSize(); i++)
385             {
386                 pDblOut->set(i, C2F(genbet)(vectpDblInput[0]->get(),vectpDblInput[1]->get()));
387             }
388
389             out.push_back(pDblOut);
390             break;
391         }
392         case 2: // binomial
393         case 3: // negative binomial
394         {
395             types::Double* pDblOut = new types::Double(iRows, iCols);
396
397             for(int i = 0; i < 2; i++)
398             {
399                 if(vectpDblInput[i]->isScalar() == false)
400                 {
401                     ScierrorW(999, _W("%ls: Wrong type for input argument #%d : A scalar expected.\n"), L"grand", i+4);
402                     return types::Function::Error;
403                 }
404             }
405
406             if(vectpDblInput[0]->get(0) < 0.0) // N
407             {
408                 ScierrorW(999, _W("%ls: Wrong value for input argument #%d : Positive integer expected.\n"), L"grand", 4);
409                 return types::Function::Error;
410             }
411
412             if(vectpDblInput[1]->get(0) < 0.0 || vectpDblInput[1]->get(0) > 1.0) // p
413             {
414                 ScierrorW(999, _W("%ls: Wrong value for input argument #%d : A value  expected.\n"), L"grand", 5);
415                 return types::Function::Error;
416             }
417
418             int N = static_cast<int>(vectpDblInput[0]->get(0));
419             if(meth == 2)
420             {
421                 for(int i = 0; i < pDblOut->getSize(); i++)
422                 {
423                     pDblOut->set(i, static_cast<double>(C2F(ignbin)(&N, vectpDblInput[1]->get())));
424                 }
425             }
426             else // meth == 3
427             {
428                 for(int i = 0; i < pDblOut->getSize(); i++)
429                 {
430                     pDblOut->set(i, static_cast<double>(C2F(ignnbn)(&N, vectpDblInput[1]->get())));
431                 }
432             }
433
434             out.push_back(pDblOut);
435             break;
436         }
437         case 4: // chisquare
438         {
439             types::Double* pDblOut = new types::Double(iRows, iCols);
440
441             if(vectpDblInput[0]->isScalar() == false)
442             {
443                 ScierrorW(999, _W("%ls: Wrong type for input argument #%d : A scalar expected.\n"), L"grand", 4);
444                 return types::Function::Error;
445             }
446
447             if(vectpDblInput[0]->get(0) <= 0.0) // Df
448             {
449                 ScierrorW(999, _W("%ls: Wrong value for input argument #%d : Positive no null value expected.\n"), L"grand", 4);
450                 return types::Function::Error;
451             }
452
453             for(int i = 0; i < pDblOut->getSize(); i++)
454             {
455                 pDblOut->set(i, C2F(genchi)(vectpDblInput[0]->get()));
456             }
457
458             out.push_back(pDblOut);
459             break;
460         }
461         case 5: // non central chisquare
462         {
463             types::Double* pDblOut = new types::Double(iRows, iCols);
464
465             for(int i = 0; i < 2; i++)
466             {
467                 if(vectpDblInput[i]->isScalar() == false)
468                 {
469                     ScierrorW(999, _W("%ls: Wrong type for input argument #%d : A scalar expected.\n"), L"grand", i+4);
470                     return types::Function::Error;
471                 }
472             }
473
474             if(vectpDblInput[0]->get(0) < 1.0) // Df
475             {
476                 ScierrorW(999, _W("%ls: Wrong value for input argument #%d : value >1.0 expected.\n"), L"grand", 4);
477                 return types::Function::Error;
478             }
479
480             if(vectpDblInput[1]->get(0) < 0.0) // Xnon
481             {
482                 ScierrorW(999, _W("%ls: Wrong value for input argument #%d : Positive value expected.\n"), L"grand", 5);
483                 return types::Function::Error;
484             }
485
486             for(int i = 0; i < pDblOut->getSize(); i++)
487             {
488                 pDblOut->set(i, C2F(gennch)(vectpDblInput[0]->get(), vectpDblInput[1]->get()));
489             }
490
491             out.push_back(pDblOut);
492             break;
493         }
494         case 6: // exponential
495         {
496             types::Double* pDblOut = new types::Double(iRows, iCols);
497
498             if(vectpDblInput[0]->isScalar() == false)
499             {
500                 ScierrorW(999, _W("%ls: Wrong type for input argument #%d : A scalar expected.\n"), L"grand", 4);
501                 return types::Function::Error;
502             }
503
504             if(vectpDblInput[0]->get(0) < 0.0) // Av
505             {
506                 ScierrorW(999, _W("%ls: Wrong value for input argument #%d : Positive value expected.\n"), L"grand", 4);
507                 return types::Function::Error;
508             }
509
510             for(int i = 0; i < pDblOut->getSize(); i++)
511             {
512                 pDblOut->set(i, C2F(genexp)(vectpDblInput[0]->get()));
513             }
514
515             out.push_back(pDblOut);
516             break;
517         }
518         case 7: // F variance ratio
519         {
520             types::Double* pDblOut = new types::Double(iRows, iCols);
521
522             for(int i = 0; i < 2; i++)
523             {
524                 if(vectpDblInput[i]->isScalar() == false)
525                 {
526                     ScierrorW(999, _W("%ls: Wrong type for input argument #%d : A scalar expected.\n"), L"grand", i+4);
527                     return types::Function::Error;
528                 }
529
530                 if(vectpDblInput[i]->get(0) <= 0.0) // Dfn Dfd
531                 {
532                     ScierrorW(999, _W("%ls: Wrong value for input argument #%d : Positive no null value expected.\n"), L"grand", i+4);
533                     return types::Function::Error;
534                 }
535             }
536
537             for(int i = 0; i < pDblOut->getSize(); i++)
538             {
539                 pDblOut->set(i, C2F(genf)(vectpDblInput[0]->get(), vectpDblInput[1]->get()));
540             }
541
542             out.push_back(pDblOut);
543             break;
544         }
545         case 8: // non F variance ratio
546         {
547             types::Double* pDblOut = new types::Double(iRows, iCols);
548
549             for(int i = 0; i < 3; i++)
550             {
551                 if(vectpDblInput[i]->isScalar() == false)
552                 {
553                     ScierrorW(999, _W("%ls: Wrong type for input argument #%d : A scalar expected.\n"), L"grand", i+4);
554                     return types::Function::Error;
555                 }
556             }
557
558             if(vectpDblInput[0]->get(0) < 1.0) // Dfn
559             {
560                 ScierrorW(999, _W("%ls: Wrong value for input argument #%d : value > 1.0 expected.\n"), L"grand", 4);
561                 return types::Function::Error;
562             }
563
564             if(vectpDblInput[1]->get(0) <= 0.0) // Dfd
565             {
566                 ScierrorW(999, _W("%ls: Wrong value for input argument #%d : Positive non null value expected.\n"), L"grand", 5);
567                 return types::Function::Error;
568             }
569
570             if(vectpDblInput[2]->get(0) < 0.0) // Xnon
571             {
572                 ScierrorW(999, _W("%ls: Wrong value for input argument #%d : Positive value expected.\n"), L"grand", 6);
573                 return types::Function::Error;
574             }
575
576             for(int i = 0; i < pDblOut->getSize(); i++)
577             {
578                 pDblOut->set(i, C2F(gennf)(vectpDblInput[0]->get(), vectpDblInput[1]->get(), vectpDblInput[2]->get()));
579             }
580
581             out.push_back(pDblOut);
582             break;
583         }
584         case 9: // gamma
585         {
586             types::Double* pDblOut = new types::Double(iRows, iCols);
587
588             for(int i = 0; i < 2; i++)
589             {
590                 if(vectpDblInput[i]->isScalar() == false)
591                 {
592                     ScierrorW(999, _W("%ls: Wrong type for input argument #%d : A scalar expected.\n"), L"grand", i+4);
593                     return types::Function::Error;
594                 }
595
596                 if(vectpDblInput[i]->get(0) <= 0.0)
597                 {
598                     ScierrorW(999, _W("%ls: Wrong value for input argument #%d : Positive non null value expected.\n"), L"grand", i+4);
599                     return types::Function::Error;
600                 }
601             }
602
603             for(int i = 0; i < pDblOut->getSize(); i++)
604             {
605                 //** WARNING ** : order is changed in parameters for
606                 //compatibility between Rand(...'gam',..) and cdfgam
607
608                 pDblOut->set(i, C2F(gengam)(vectpDblInput[1]->get(), vectpDblInput[0]->get()));
609             }
610
611             out.push_back(pDblOut);
612             break;
613         }
614         case 10: // Gauss Laplace (normal)
615         {
616             types::Double* pDblOut = new types::Double(iRows, iCols);
617
618             for(int i = 0; i < 2; i++)
619             {
620                 if(vectpDblInput[i]->isScalar() == false)
621                 {
622                     ScierrorW(999, _W("%ls: Wrong type for input argument #%d : A scalar expected.\n"), L"grand", i+4);
623                     return types::Function::Error;
624                 }
625             }
626
627             if(vectpDblInput[1]->get(0) < 0.0)
628             {
629                 ScierrorW(999, _W("%ls: Wrong value for input argument #%d : Positive value expected.\n"), L"grand", 5);
630                 return types::Function::Error;
631             }
632
633             for(int i = 0; i < pDblOut->getSize(); i++)
634             {
635                 pDblOut->set(i, C2F(gennor)(vectpDblInput[0]->get(), vectpDblInput[1]->get()));
636             }
637
638             out.push_back(pDblOut);
639             break;
640         }
641         case 11: // multivariate gaussian (multivariate normal)
642         {
643             if(vectpDblInput[0]->getCols() != 1 || vectpDblInput[0]->getSize() == 0)
644             {
645                 ScierrorW(999, _W("%ls: Wrong size for input argument #%d : A matrix of size m x 1 expected.(m > 0)\n"), L"grand", 4);
646                 return types::Function::Error;
647             }
648
649             if(vectpDblInput[0]->getRows() != vectpDblInput[1]->getRows())
650             {
651                 ScierrorW(999, _W("%ls: Wrong size for input argument #%d and #%d: Mean and Cov have incompatible dimensions.\n"), L"grand", 4, 5);
652                 return types::Function::Error;
653             }
654
655             if(vectpDblInput[1]->getRows() != vectpDblInput[1]->getCols())
656             {
657                 ScierrorW(999, _W("%ls: Wrong size for input argument #%d : A square symmetric positive definite matrix expected.\n"), L"grand", 5);
658                 return types::Function::Error;
659             }
660
661             int ierr = 0;
662             int size = vectpDblInput[0]->getRows();
663             int mp   = (size * (size + 3)) / 2 + 1;
664
665             types::Double* pDblOut = new types::Double(size, iNumIter);
666             double* work  = (double*)malloc(size * sizeof(double));
667             double* param = (double*)malloc(mp   * sizeof(double));
668
669             types::Double* pDblMean = vectpDblInput[0]->clone()->getAs<types::Double>();
670             types::Double* pDblCov  = vectpDblInput[1]->clone()->getAs<types::Double>();
671
672             C2F(setgmn)(pDblMean->get(), pDblCov->get(), &size, &size, param, &ierr);
673
674             delete pDblMean;
675             delete pDblCov;
676
677             if(ierr == 1)
678             {
679                 free(work);
680                 free(param);
681                 ScierrorW(999, _W("%ls: setgmn return with state %d.\n"), L"grand", ierr);
682                 return types::Function::Error;
683             }
684
685             for(int i = 0; i < iNumIter; i++)
686             {
687                 C2F(genmn)(param, pDblOut->get() + (size * i), work);
688             }
689
690             free(work);
691             free(param);
692
693             out.push_back(pDblOut);
694             break;
695         }
696         case 12: // geometric
697         {
698             types::Double* pDblOut = new types::Double(iRows, iCols);
699             double pmin = 1.3e-307;
700
701             if(vectpDblInput[0]->isScalar() == false)
702             {
703                 ScierrorW(999, _W("%ls: Wrong type for input argument #%d : A scalar expected.\n"), L"grand", 4);
704                 return types::Function::Error;
705             }
706
707             if(vectpDblInput[0]->get(0) < pmin || vectpDblInput[0]->get(0) > 1.0)
708             {
709                 ScierrorW(999, _W("%ls: Wrong value for input argument #%d : Must be between %lf and %d.\n"), L"grand", 4, pmin, 1);
710                 return types::Function::Error;
711             }
712
713             for(int i = 0; i < pDblOut->getSize(); i++)
714             {
715                 pDblOut->set(i, igngeom(vectpDblInput[0]->get(0)));
716             }
717
718             out.push_back(pDblOut);
719             break;
720         }
721         case 13: // markov
722         {
723             if( vectpDblInput[0]->getRows() != vectpDblInput[0]->getCols() &&
724                 vectpDblInput[0]->getRows() != 1)
725             {
726                 ScierrorW(999, _W("%ls: Wrong size for input argument #%d : A square matrix or a row vector expected.\n"), L"grand", 4);
727                 return types::Function::Error;
728             }
729
730             if(vectpDblInput[1]->getSize() == 0)
731             {
732                 ScierrorW(999, _W("%ls: Wrong size for input argument #%d: No empty matrix expected.\n"), L"grand", 5);
733                 return types::Function::Error;
734             }
735
736             int sizeOfX0 = vectpDblInput[1]->getSize();
737
738             for( int i = 0; i < sizeOfX0; i++)
739             {
740                 if(vectpDblInput[1]->get(i) < 1 || (vectpDblInput[1]->get(i) - 1) >= vectpDblInput[0]->getCols())
741                 {
742                     ScierrorW(999, _W("%ls: X0(%d) must be in the range [1,%d[.\n"), L"grand", i+1, vectpDblInput[0]->getCols() + 1);
743                     return types::Function::Error;
744                 }
745             }
746
747             types::Double* pDblOut = new types::Double(sizeOfX0, iNumIter);
748             int size = vectpDblInput[0]->getSize() + vectpDblInput[0]->getCols();
749             double* work = (double*)malloc(size * sizeof(double));
750
751             for(int i = 0; i < vectpDblInput[0]->getRows(); i++)
752             {
753                 double ptot = 0.0;
754                 for(int j = 0; j < vectpDblInput[0]->getCols(); j++)
755                 {
756                     int position = vectpDblInput[0]->getRows() * j + i;
757
758                     if(vectpDblInput[0]->get(position) < 0 || vectpDblInput[0]->get(position) > 1)
759                     {
760                         ScierrorW(999, _W("%ls: Wrong value for input argument #%d: P(%d,%d) must be in the range [0 1].\n"), L"grand", i+1, j+1);
761                         return types::Function::Error;
762                     }
763
764                     ptot += vectpDblInput[0]->get(position);
765                 }
766
767                 if(fabs(ptot - 1.0) > 1e-8)
768                 {
769                     ScierrorW(999, _W("%ls: Sum of P(%d,1:%d)=%lf ~= 1.\n"), L"grand", i+1, vectpDblInput[0]->getCols(), ptot);
770                     return types::Function::Error;
771                 }
772             }
773             // ** Computing the cumulative sum of the P matrix **
774
775             for(int i = 0; i < vectpDblInput[0]->getRows(); i++)
776             {
777                 double cumsum = 0.0;
778                 work[i] = 0.0;
779                 for(int j = 1; j < vectpDblInput[0]->getCols() + 1; j++)
780                 {
781                     cumsum += vectpDblInput[0]->get(vectpDblInput[0]->getRows() * (j - 1) + i);
782                     work[vectpDblInput[0]->getRows() * j + i] = cumsum;
783                 }
784             }
785
786             for(int i = 0; i < vectpDblInput[1]->getSize(); i++)
787             {
788                 int iCur = static_cast<int>(vectpDblInput[1]->get(i)) - 1;
789                 for(int j = 0; j < iNumIter; j++)
790                 {
791                     int niv = 0;
792                     double rr = C2F(ranf)();
793
794                     if(vectpDblInput[0]->getRows() == 1)
795                     {
796                         iCur = 0;
797                     }
798
799                     while(  rr >= work[iCur + vectpDblInput[0]->getRows() * niv] &&
800                             niv < (vectpDblInput[0]->getRows() + 1))
801                     {
802                         niv++;
803                     }
804
805                     // ** projection to avoid boundaries **
806                     niv = Max(Min(niv, vectpDblInput[0]->getCols()), 1);
807                     pDblOut->set(vectpDblInput[1]->getSize() * j + i, static_cast<double>(niv));
808                     iCur = niv - 1;
809                 }
810             }
811
812             out.push_back(pDblOut);
813             break;
814         }
815         case 14: // multinomial
816         {
817             if(vectpDblInput[0]->isScalar() == false)
818             {
819                 ScierrorW(999, _W("%ls: Wrong size for input argument #%d: A scalar expected.\n"), L"grand", 4);
820                 return types::Function::Error;
821             }
822
823             if(vectpDblInput[0]->get(0) < 0)
824             {
825                 ScierrorW(999, _W("%ls: Wrong value for input argument #%d: A positive scalar expected.\n"), L"grand", 4);
826                 return types::Function::Error;
827             }
828
829             if(vectpDblInput[1]->getCols() != 1 || vectpDblInput[1]->getRows() <= 0)
830             {
831                 ScierrorW(999, _W("%ls: Wrong size for input argument #%d: A colomn vector expected.\n"), L"grand", 5);
832                 return types::Function::Error;
833             }
834
835             double ptot = 0.0;
836             int ncat = vectpDblInput[1]->getRows() + 1;
837             types::Double* pDblOut = new types::Double(ncat, iNumIter);
838
839             for(int i = 0; i < vectpDblInput[1]->getCols(); i++)
840             {
841                 if(vectpDblInput[1]->get(i) < 0.0 || vectpDblInput[1]->get(i) > 1.0)
842                 {
843                     ScierrorW(999, _W("%ls: Wrong value for input argument #%d: P(%d) must be in the range [0 1].\n"), L"grand", 5, i+1);
844                     return types::Function::Error;
845                 }
846                 ptot += vectpDblInput[1]->get(i);
847             }
848
849             if(ptot > 1.0)
850             {
851                 ScierrorW(999, _W("%ls: Wrong value for input argument #%d: Sum of P(i) > 1.\n"), L"grand", 5);
852                 return types::Function::Error;
853             }
854
855             int* piP    = (int*)malloc(vectpDblInput[0]->getSize() * sizeof(int));
856             int* piOut  = (int*)malloc(pDblOut->getSize() * sizeof(int));
857
858             for(int i = 0; i < vectpDblInput[0]->getSize(); i++)
859             {
860                 piP[i] = static_cast<int>(vectpDblInput[0]->get(i));
861             }
862
863             for(int i = 0; i < iNumIter; i++)
864             {
865                 C2F(genmul)(piP, vectpDblInput[1]->get(), &ncat,  + &piOut[ncat * i]);
866             }
867
868             for(int i = 0; i < pDblOut->getSize(); i++)
869             {
870                 pDblOut->set(i, static_cast<double>(piOut[i]));
871             }
872
873             free(piP);
874             free(piOut);
875
876             out.push_back(pDblOut);
877             break;
878         }
879         case 15: // Poisson
880         {
881             types::Double* pDblOut = new types::Double(iRows, iCols);
882
883             if(vectpDblInput[0]->isScalar() == false)
884             {
885                 ScierrorW(999, _W("%ls: Wrong type for input argument #%d : A scalar expected.\n"), L"grand", 4);
886                 return types::Function::Error;
887             }
888
889             if(vectpDblInput[0]->get(0) < 0.0)
890             {
891                 ScierrorW(999, _W("%ls: Wrong value for input argument #%d : Positive value expected.\n"), L"grand", 4);
892                 return types::Function::Error;
893             }
894
895             for(int i = 0; i < pDblOut->getSize(); i++)
896             {
897                 pDblOut->set(i, static_cast<double>(C2F(ignpoi)(vectpDblInput[0]->get())));
898             }
899
900             out.push_back(pDblOut);
901             break;
902         }
903         case 16: // random permutations
904         {
905             if(vectpDblInput[0]->getCols() != 1)
906             {
907                 ScierrorW(999, _W("%ls: Wrong type for input argument #%d : A colomn vector expected.\n"), L"grand", 4);
908                 return types::Function::Error;
909             }
910
911             int iVectRows = vectpDblInput[0]->getRows();
912             types::Double* pDblOut = new types::Double(iVectRows, iNumIter);
913
914             for(int i = 0; i < iNumIter; i++)
915             {
916                 for(int j = 0; j < iVectRows; j++)
917                 {
918                     pDblOut->set(iVectRows * i + j, vectpDblInput[0]->get(j));
919                 }
920                 C2F(genprm)(pDblOut->get() + (iVectRows * i), &iVectRows);
921             }
922
923             out.push_back(pDblOut);
924             break;
925         }
926         case 17: // uniform (def)
927         {
928             types::Double* pDblOut = new types::Double(iRows, iCols);
929
930             for(int i = 0; i < pDblOut->getSize(); i++)
931             {
932                 pDblOut->set(i, C2F(ranf)());
933             }
934
935             out.push_back(pDblOut);
936             break;
937         }
938         case 18: // uniform (unf)
939         {
940             types::Double* pDblOut = new types::Double(iRows, iCols);
941
942             for(int i = 0; i < 2; i++)
943             {
944                 if(vectpDblInput[i]->isScalar() == false)
945                 {
946                     ScierrorW(999, _W("%ls: Wrong type for input argument #%d : A scalar expected.\n"), L"grand", i+4);
947                     return types::Function::Error;
948                 }
949             }
950
951             double low  = vectpDblInput[0]->get(0);
952             double high = vectpDblInput[1]->get(0);
953
954             if(low > high)
955             {
956                 ScierrorW(999, _W("%ls: Wrong value for input argument #%d and #%d: Low < High expected.\n"), L"grand", 4, 5);
957                 return types::Function::Error;
958             }
959
960             for(int i = 0; i < pDblOut->getSize(); i++)
961             {
962                 pDblOut->set(i, low + (high - low) * C2F(ranf)());
963             }
964
965             out.push_back(pDblOut);
966             break;
967         }
968         case 19: // uniform (uin)
969         {
970             types::Double* pDblOut = new types::Double(iRows, iCols);
971
972             for(int i = 0; i < 2; i++)
973             {
974                 if(vectpDblInput[i]->isScalar() == false)
975                 {
976                     ScierrorW(999, _W("%ls: Wrong type for input argument #%d : A scalar expected.\n"), L"grand", i+4);
977                     return types::Function::Error;
978                 }
979             }
980
981             int low  = static_cast<int>(vectpDblInput[0]->get(0));
982             int high = static_cast<int>(vectpDblInput[1]->get(0));
983
984             if(low > high)
985             {
986                 ScierrorW(999, _W("%ls: Wrong value for input argument #%d and #%d: Low < High expected.\n"), L"grand", 4, 5);
987                 return types::Function::Error;
988             }
989
990             if( low  != vectpDblInput[0]->get(0) ||
991                 high != vectpDblInput[1]->get(0) ||
992                 (high - low + 1) > 2147483561)
993             {
994                 ScierrorW(999, _W("%ls: Wrong value for input argument #%d and #%d: Low and High must be integers and (high - low + 1) <=  2147483561.\n"), L"grand", 4, 5);
995                 return types::Function::Error;
996             }
997
998             for(int i = 0; i < pDblOut->getSize(); i++)
999             {
1000                 pDblOut->set(i, C2F(ignuin)(vectpDblInput[0]->get(), vectpDblInput[1]->get()));
1001             }
1002
1003             out.push_back(pDblOut);
1004             break;
1005         }
1006         case 20: // uniform (lgi)
1007         {
1008             types::Double* pDblOut = new types::Double(iRows, iCols);
1009
1010             for(int i = 0; i < pDblOut->getSize(); i++)
1011             {
1012                 pDblOut->set(i, static_cast<double>(ignlgi()));
1013             }
1014
1015             out.push_back(pDblOut);
1016             break;
1017         }
1018         case 21: // getgen
1019         {
1020             types::String* pStrOut = new types::String(names_gen[current_base_gen]);
1021             out.push_back(pStrOut);
1022             break;
1023         }
1024         case 22: // setgen
1025         {
1026             wchar_t* wcsGen = pStrGenOrPhr->get(0);
1027
1028             if(wcscmp(wcsGen, L"mt") == 0)
1029             {
1030                 ConfigVariable::setCurrentBaseGen(MT);
1031             }
1032             else if(wcscmp(wcsGen, L"kiss") == 0)
1033             {
1034                 ConfigVariable::setCurrentBaseGen(KISS);
1035             }
1036             else if(wcscmp(wcsGen, L"clcg4") == 0)
1037             {
1038                 ConfigVariable::setCurrentBaseGen(CLCG4);
1039             }
1040             else if(wcscmp(wcsGen, L"clcg2") == 0)
1041             {
1042                 ConfigVariable::setCurrentBaseGen(CLCG2);
1043             }
1044             else if(wcscmp(wcsGen, L"urand") == 0)
1045             {
1046                 ConfigVariable::setCurrentBaseGen(URAND);
1047             }
1048             else if(wcscmp(wcsGen, L"fsultra") == 0)
1049             {
1050                 ConfigVariable::setCurrentBaseGen(FSULTRA);
1051             }
1052             else
1053             {
1054                 ScierrorW(999, _W("%ls: Wrong value for input argument #%d: '%ls', '%ls', '%ls', '%ls', '%ls' or '%ls' expected.\n"), L"grand", 2, L"mt", L"kiss", L"clcg4", L"clcg2", L"urand", L"fsultra");
1055                 return types::Function::Error;
1056             }
1057
1058             int current_gen = ConfigVariable::getCurrentBaseGen();
1059             types::String* pStrOut = new types::String(names_gen[current_gen]);
1060             out.push_back(pStrOut);
1061             break;
1062         }
1063         case 23: // getsd
1064         {
1065             types::Double* pDblOut = NULL;
1066
1067             switch(current_base_gen)
1068             {
1069                 case MT:
1070                 {
1071                     pDblOut = new types::Double(625, 1);
1072                     get_state_mt(pDblOut->get());
1073                     break;
1074                 }
1075                 case KISS:
1076                 {
1077                     pDblOut = new types::Double(4, 1);
1078                     get_state_kiss(pDblOut->get());
1079                     break;
1080                 }
1081                 case CLCG4:
1082                 {
1083                     pDblOut = new types::Double(4, 1);
1084                     int current_clcg4 = ConfigVariable::getCurrentClcg4();
1085                     get_state_clcg4(current_clcg4, pDblOut->get());
1086                     break;
1087                 }
1088                 case CLCG2:
1089                 {
1090                     pDblOut = new types::Double(2, 1);
1091                     get_state_clcg2(pDblOut->get());
1092                     break;
1093                 }
1094                 case URAND:
1095                 {
1096                     pDblOut = new types::Double(1, 1);
1097                     get_state_urand(pDblOut->get());
1098                     break;
1099                 }
1100                 case FSULTRA:
1101                 {
1102                     pDblOut = new types::Double(40, 1);
1103                     get_state_fsultra(pDblOut->get());
1104                     break;
1105                 }
1106             }
1107
1108             out.push_back(pDblOut);
1109             break;
1110         }
1111         case 24: // setsd
1112         {
1113             int ierr = 0;
1114             switch(current_base_gen)
1115             {
1116                 case MT:
1117                 {
1118                     if(vectpDblInput[0]->isScalar())
1119                     {
1120                         ierr = set_state_mt_simple(vectpDblInput[0]->get(0));
1121                     }
1122                     else if(vectpDblInput[0]->getSize() != 625)
1123                     {
1124                         ierr = set_state_mt(vectpDblInput[0]->get());
1125                     }
1126                     else
1127                     {
1128                         ScierrorW(999, _W("%ls: Wrong size for input argument #%d : A scalar or a vector of size %d expected.\n"), L"grand", 4, 625);
1129                         return types::Function::Error;
1130                     }
1131
1132                     break;
1133                 }
1134                 case KISS:
1135                 case CLCG4:
1136                 {
1137                     for(int i = 0; i < 4; i++)
1138                     {
1139                         if(vectpDblInput[i]->isScalar() == false)
1140                         {
1141                             ScierrorW(999, _W("%ls: Wrong type for input argument #%d : A scalar expected.\n"), L"grand", i+2);
1142                             return types::Function::Error;
1143                         }
1144                     }
1145
1146                     if(current_base_gen == KISS)
1147                     {
1148                         ierr = set_state_kiss(  vectpDblInput[0]->get(0),
1149                                                 vectpDblInput[1]->get(0),
1150                                                 vectpDblInput[2]->get(0),
1151                                                 vectpDblInput[3]->get(0));
1152                     }
1153                     else // CLCG4
1154                     {
1155                         ierr = set_seed_clcg4(  ConfigVariable::getCurrentClcg4(),
1156                                                 vectpDblInput[0]->get(0),
1157                                                 vectpDblInput[1]->get(0),
1158                                                 vectpDblInput[2]->get(0),
1159                                                 vectpDblInput[3]->get(0));
1160                     }
1161
1162                     break;
1163                 }
1164                 case CLCG2:
1165                 {
1166                     for(int i = 0; i < 2; i++)
1167                     {
1168                         if(vectpDblInput[i]->isScalar() == false)
1169                         {
1170                             ScierrorW(999, _W("%ls: Wrong type for input argument #%d : A scalar expected.\n"), L"grand", i+2);
1171                             return types::Function::Error;
1172                         }
1173                     }
1174
1175                     ierr = set_state_clcg2(vectpDblInput[0]->get(0), vectpDblInput[1]->get(0));
1176                     break;
1177                 }
1178                 case URAND:
1179                 {
1180                     if(vectpDblInput[0]->isScalar() == false)
1181                     {
1182                         ScierrorW(999, _W("%ls: Wrong type for input argument #%d : A scalar expected.\n"), L"grand", 2);
1183                         return types::Function::Error;
1184                     }
1185
1186                     ierr = set_state_urand(vectpDblInput[0]->get(0));
1187                     break;
1188                 }
1189                 case FSULTRA:
1190                 {
1191                     if(in.size() == 2)
1192                     {
1193                         if(vectpDblInput[0]->getRows() != 40 || vectpDblInput[0]->getCols() != 1)
1194                         {
1195                             ScierrorW(999, _W("%ls: Wrong size for input argument #%d : A vector of size %d x %d expected.\n"), L"grand", 2, 40, 1);
1196                             return types::Function::Error;
1197                         }
1198
1199                         ierr = set_state_fsultra(vectpDblInput[0]->get());
1200                     }
1201                     else // in.size() == 3
1202                     {
1203                         for(int i = 0; i < 2; i++)
1204                         {
1205                             if(vectpDblInput[i]->isScalar() == false)
1206                             {
1207                                 ScierrorW(999, _W("%ls: Wrong type for input argument #%d : A scalar expected.\n"), L"grand", i+2);
1208                                 return types::Function::Error;
1209                             }
1210                         }
1211
1212                         ierr = set_state_fsultra_simple(vectpDblInput[0]->get(0), vectpDblInput[1]->get(0));
1213                     }
1214
1215                     break;
1216                 }
1217             }
1218
1219             if(ierr == 0)
1220             {
1221                 ScierrorW(999, _W("%ls: Wrong value for the last %d input argument(s).\n"), L"grand", in.size() - 1);
1222                 return types::Function::Error;
1223             }
1224
1225             break;
1226         }
1227         case 25: // phr2sd
1228         {
1229             if(pStrGenOrPhr->isScalar() == false)
1230             {
1231                 ScierrorW(999, _W("%ls: Wrong type for input argument #%d : One string expected.\n"), L"grand", 2);
1232                 return types::Function::Error;
1233             }
1234
1235             types::Double* pDblOut = new types::Double(1, 2);
1236             int size = wcslen(pStrGenOrPhr->get(0));
1237             int piOut[2];
1238             char* strPhr = wide_string_to_UTF8(pStrGenOrPhr->get(0));
1239
1240             C2F(phrtsd)(strPhr, &size, piOut, piOut + 1, size);
1241
1242             pDblOut->set(0, static_cast<double>(piOut[0]));
1243             pDblOut->set(1, static_cast<double>(piOut[1]));
1244
1245             FREE(strPhr);
1246             out.push_back(pDblOut);
1247             break;
1248         }
1249         case 26: // setcgn
1250         {
1251             if(current_base_gen != CLCG4)
1252             {
1253                 sciprintW(_W("The %ls option affects only the %ls generator\n"), L"setcgn", L"clcg4");
1254             }
1255
1256             if(vectpDblInput[0]->isScalar() == false)
1257             {
1258                 ScierrorW(999, _W("%ls: Wrong type for input argument #%d : A scalar expected.\n"), L"grand", 2);
1259                 return types::Function::Error;
1260             }
1261
1262             if(vectpDblInput[0]->get(0) < 0 || vectpDblInput[0]->get(0) > Maxgen)
1263             {
1264                ScierrorW(999, _W("%ls: Wrong value for input argument #%d : Must be between %d and %d.\n"), L"grand", 0, Maxgen);
1265                 return types::Function::Error;
1266             }
1267
1268             ConfigVariable::setCurrentClcg4(static_cast<int>(vectpDblInput[0]->get(0)));
1269             double dOut = static_cast<double>(ConfigVariable::getCurrentClcg4());
1270             out.push_back(new types::Double(dOut));
1271             break;
1272         }
1273         case 27: // getcgn
1274         {
1275             double dOut = static_cast<double>(ConfigVariable::getCurrentClcg4());
1276             out.push_back(new types::Double(dOut));
1277             break;
1278         }
1279         case 28: // initgn
1280         {
1281             SeedType where;
1282             if(current_base_gen != CLCG4)
1283             {
1284                 sciprintW(_W("The %ls option affects only the %ls generator\n"), L"initgn", L"clcg4");
1285             }
1286
1287             if(vectpDblInput[0]->isScalar() == false)
1288             {
1289                 ScierrorW(999, _W("%ls: Wrong type for input argument #%d : A scalar expected.\n"), L"grand", 2);
1290                 return types::Function::Error;
1291             }
1292
1293             if( vectpDblInput[0]->get(0) != 0 &&
1294                 vectpDblInput[0]->get(0) != -1 &&
1295                 vectpDblInput[0]->get(0) != 1)
1296             {
1297                ScierrorW(999, _W("%ls: Wrong value for input argument #%d : Must be between %d, %d or %d.\n"), L"grand", 2, -1, 0, 1);
1298                 return types::Function::Error;
1299             }
1300
1301             where = (SeedType)(int)(vectpDblInput[0]->get(0) + 1);
1302             init_generator_clcg4(ConfigVariable::getCurrentClcg4(), where);
1303             out.push_back(vectpDblInput[0]);
1304             break;
1305         }
1306         case 29: // setall
1307         {
1308             if(current_base_gen != CLCG4)
1309             {
1310                 sciprintW(_W("The %ls option affects only the %ls generator\n"), L"setall", L"clcg4");
1311             }
1312
1313             for(int i = 0; i < 4; i++)
1314             {
1315                 if(vectpDblInput[i]->isScalar() == false)
1316                 {
1317                     ScierrorW(999, _W("%ls: Wrong type for input argument #%d : A scalar expected.\n"), L"grand", i+2);
1318                     return types::Function::Error;
1319                 }
1320             }
1321
1322             int ierr = set_initial_seed_clcg4(  vectpDblInput[0]->get(0),
1323                                                 vectpDblInput[1]->get(0),
1324                                                 vectpDblInput[2]->get(0),
1325                                                 vectpDblInput[3]->get(0));
1326             if(ierr == 0)
1327             {
1328                 ScierrorW(999, _W("%ls: Wrong value for the last %d input argument(s).\n"), L"grand", 4);
1329                 return types::Function::Error;
1330             }
1331
1332             out.push_back(pStrMethod);
1333             break;
1334         }
1335         case 30: // advnst
1336         {
1337             if(current_base_gen != CLCG4)
1338             {
1339                 sciprintW(_W("The %ls option affects only the %ls generator\n"), L"advnst", L"clcg4");
1340             }
1341
1342             if(vectpDblInput[0]->isScalar() == false)
1343             {
1344                 ScierrorW(999, _W("%ls: Wrong type for input argument #%d : A scalar expected.\n"), L"grand", 2);
1345                 return types::Function::Error;
1346             }
1347
1348             int k = static_cast<int>(vectpDblInput[0]->get(0));
1349
1350             if(k < 1)
1351             {
1352                 ScierrorW(999, _W("%ls: Wrong value for input argument #%d : Must be > %d.\n"), L"grand", 2, 0);
1353                 return types::Function::Error;
1354             }
1355
1356             advance_state_clcg4(ConfigVariable::getCurrentClcg4(), k);
1357             out.push_back(new types::Double(static_cast<double>(k)));
1358             break;
1359         }
1360     }
1361
1362     return types::Function::OK;
1363 }
1364 /*--------------------------------------------------------------------------*/
1365