* Bug #15243 fixed - grand() did not check its inputs
[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 * Copyright (C) 2014 - Scilab Enterprises - Sylvain GENIN
5 *
6  * Copyright (C) 2012 - 2016 - Scilab Enterprises
7  *
8  * This file is hereby licensed under the terms of the GNU GPL v2.0,
9  * pursuant to article 5.3.4 of the CeCILL v.2.1.
10  * This file was originally licensed under the terms of the CeCILL v2.1,
11  * and continues to be available under such terms.
12  * For more information, see the COPYING file which you should have received
13  * along with this program.
14 *
15 */
16 /*--------------------------------------------------------------------------*/
17 #include "randlib_gw.hxx"
18 #include "function.hxx"
19 #include "double.hxx"
20 #include "string.hxx"
21 #include "configvariable.hxx"
22 #include "int.hxx"
23 #include "polynom.hxx"
24 #include "sparse.hxx"
25 #include "overload.hxx"
26
27 extern "C"
28 {
29 #include "sci_malloc.h"
30 #include "localization.h"
31 #include "Scierror.h"
32 #include "sciprint.h"
33 #include "grand.h"
34 #include "clcg4.h"
35 #include "others_generators.h"
36 }
37 /*--------------------------------------------------------------------------*/
38
39 template<class U>
40 void sci_grand_prm(int iNumIter, U *pIn, types::InternalType** pOut);
41
42 types::Function::ReturnValue sci_grand(types::typed_list &in, int _iRetCount, types::typed_list &out)
43 {
44     enum
45     {
46         MT, KISS, CLCG4, CLCG2, URAND
47     };
48
49     //  names at the scilab level
50     const wchar_t* names_gen[6] = { L"mt", L"kiss", L"clcg4", L"clcg2", L"urand" };
51
52     types::String* pStrMethod = NULL;
53     types::String* pStrGenOrPhr = NULL;
54
55     std::vector<types::Double*> vectpDblInput;
56
57     int iStrPos = 0;
58     int iPos = 0;
59     int meth = -1;
60     int iRows = -1;
61     int iCols = -1;
62     int iNumIter = -1;
63
64
65     int current_base_gen = ConfigVariable::getCurrentBaseGen();
66
67     // *** check the maximal number of input args. ***
68     /* if (in.size() > 6)
69     {
70     Scierror(77, _("%s: Wrong number of input argument(s): %d to %d expected.\n"), "grand", 1, 6);
71     return types::Function::Error;
72     }*/
73
74     // *** check number of output args. ***
75     if (_iRetCount > 1)
76     {
77         Scierror(78, _("%s: Wrong number of output argument(s): %d expected.\n"), "grand", 1);
78         return types::Function::Error;
79     }
80
81     // *** find the method string. ***
82     for (int i = 0; i < in.size(); i++)
83     {
84         if (in[i]->isString())
85         {
86             pStrMethod = in[i]->getAs<types::String>();
87             iStrPos = i;
88             break;
89         }
90     }
91
92
93     if (pStrMethod == NULL)
94     {
95         for (int i = 0; i < in.size(); i++)
96         {
97             if (in[i]->isDouble() == false)
98             {
99                 std::wstring wstFuncName = L"%" + in[0]->getShortTypeStr() + L"_grand";
100                 return Overload::call(wstFuncName, in, _iRetCount, out);
101             }
102         }
103
104         Scierror(999, _("%s: Wrong type for input argument #%d: A string expected.\n"), "grand", in.size() + 1);
105         return types::Function::Error;
106     }
107
108     int iDims = iStrPos > 1 ? iStrPos : 2;
109     int *itab = new int[iDims];
110
111     //all variables are matrix
112     itab[0] = 1;
113     itab[1] = 1;
114
115     wchar_t* wcsMeth = pStrMethod->get(0);
116     int iNumInputArg = 2;
117     if (wcscmp(wcsMeth, L"bet") == 0) // beta
118     {
119         meth = 1;
120     }
121     else if (wcscmp(wcsMeth, L"bin") == 0) // binomial
122     {
123         meth = 2;
124     }
125     else if (wcscmp(wcsMeth, L"nbn") == 0) // negative binomial
126     {
127         meth = 3;
128     }
129     else if (wcscmp(wcsMeth, L"chi") == 0) // chisquare
130     {
131         iNumInputArg = 1;
132         meth = 4;
133     }
134     else if (wcscmp(wcsMeth, L"nch") == 0) // non central chisquare
135     {
136         meth = 5;
137     }
138     else if (wcscmp(wcsMeth, L"exp") == 0) // exponential
139     {
140         iNumInputArg = 1;
141         meth = 6;
142     }
143     else if (wcscmp(wcsMeth, L"f") == 0) // F variance ratio
144     {
145         meth = 7;
146     }
147     else if (wcscmp(wcsMeth, L"nf") == 0) // non central F variance ratio
148     {
149         iNumInputArg = 3;
150         meth = 8;
151     }
152     else if (wcscmp(wcsMeth, L"gam") == 0) // gamma
153     {
154         meth = 9;
155     }
156     else if (wcscmp(wcsMeth, L"nor") == 0) // Gauss Laplace (normal)
157     {
158         meth = 10;
159     }
160     else if (wcscmp(wcsMeth, L"mn") == 0) // multivariate gaussian (multivariate normal)
161     {
162         iNumInputArg = 2;
163         meth = 11;
164     }
165     else if (wcscmp(wcsMeth, L"geom") == 0) // geometric
166     {
167         iNumInputArg = 1;
168         meth = 12;
169     }
170     else if (wcscmp(wcsMeth, L"markov") == 0) // markov
171     {
172         iNumInputArg = 2;
173         meth = 13;
174     }
175     else if (wcscmp(wcsMeth, L"mul") == 0) // multinomial
176     {
177         iNumInputArg = 2;
178         meth = 14;
179     }
180     else if (wcscmp(wcsMeth, L"poi") == 0) // Poisson
181     {
182         iNumInputArg = 1;
183         meth = 15;
184     }
185     else if (wcscmp(wcsMeth, L"prm") == 0) // random permutations
186     {
187         iNumInputArg = 1;
188         meth = 16;
189     }
190     else if (wcscmp(wcsMeth, L"def") == 0) // uniform (def)
191     {
192         iNumInputArg = 0;
193         meth = 17;
194     }
195     else if (wcscmp(wcsMeth, L"unf") == 0) // uniform (unf)
196     {
197         meth = 18;
198     }
199     else if (wcscmp(wcsMeth, L"uin") == 0) // uniform (uin)
200     {
201         meth = 19;
202     }
203     else if (wcscmp(wcsMeth, L"lgi") == 0) // uniform (lgi)
204     {
205         iNumInputArg = 0;
206         meth = 20;
207     }
208     else if (wcscmp(wcsMeth, L"getgen") == 0) // getgen
209     {
210         iNumInputArg = 0;
211         meth = 21;
212     }
213     else if (wcscmp(wcsMeth, L"setgen") == 0) // setgen
214     {
215         iNumInputArg = 1;
216         meth = 22;
217     }
218     else if (wcscmp(wcsMeth, L"getsd") == 0) // getsd
219     {
220         iNumInputArg = 0;
221         meth = 23;
222     }
223     else if (wcscmp(wcsMeth, L"setsd") == 0) // setsd
224     {
225         switch (current_base_gen)
226         {
227             case MT:
228                 iNumInputArg = 1;
229                 break;
230             case KISS:
231                 iNumInputArg = 4;
232                 break;
233             case CLCG2:
234                 iNumInputArg = 2;
235                 break;
236             case CLCG4:
237                 iNumInputArg = 4;
238                 break;
239             case URAND:
240                 iNumInputArg = 1;
241                 break;
242         }
243         meth = 24;
244     }
245     else if (wcscmp(wcsMeth, L"phr2sd") == 0) // phr2sd
246     {
247         iNumInputArg = 1;
248         meth = 25;
249     }
250     else if (wcscmp(wcsMeth, L"setcgn") == 0) // setcgn
251     {
252         iNumInputArg = 1;
253         meth = 26;
254     }
255     else if (wcscmp(wcsMeth, L"getcgn") == 0) // getcgn
256     {
257         iNumInputArg = 0;
258         meth = 27;
259     }
260     else if (wcscmp(wcsMeth, L"initgn") == 0) // initgn
261     {
262         iNumInputArg = 1;
263         meth = 28;
264     }
265     else if (wcscmp(wcsMeth, L"setall") == 0) // setall
266     {
267         iNumInputArg = 4;
268         meth = 29;
269     }
270     else if (wcscmp(wcsMeth, L"advnst") == 0) // advnst
271     {
272         iNumInputArg = 1;
273         meth = 30;
274     }
275     else
276     {
277         char* pstMeth = wide_string_to_UTF8(wcsMeth);
278         Scierror(78, _("%s: Wrong method for input argument #%d: %s unknown.\n"), "grand", iStrPos + 1, pstMeth);
279         FREE(pstMeth);
280         delete[] itab;
281         return types::Function::Error;
282     }
283
284     // *** get arguments before methode string. ***
285     if (meth == 11 || meth == 13 || meth == 14 || meth == 16) // grand(n, "...", ...
286     {
287         if (iStrPos != 1)
288         {
289             Scierror(999, _("%s: Wrong position for input argument #%d : Must be in position %d.\n"), "grand", iStrPos + 1, 2);
290             delete[] itab;
291             return types::Function::Error;
292         }
293
294         if (in[iPos]->isDouble() == false)
295         {
296             Scierror(999, _("%s: Wrong type for input argument #%d : A matrix expected.\n"), "grand", iPos + 1);
297             delete[] itab;
298             return types::Function::Error;
299         }
300
301         types::Double* pDblTemp = in[iPos]->getAs<types::Double>();
302
303         if (pDblTemp->isScalar() == false)
304         {
305             Scierror(999, _("%s: Wrong type for input argument #%d : A scalar expected.\n"), "grand", iPos + 1);
306             delete[] itab;
307             return types::Function::Error;
308         }
309
310         iNumIter = (int)pDblTemp->get(0);
311         iPos++;
312     }
313     else if (meth < 21) // grand(m, n, "...", ... || grand(matrix, "...", ...
314     {
315
316         /*if (iStrPos != 1 && iStrPos != 2)
317         {
318         Scierror(999, _("%s: Wrong position for input argument #%d : Must be in position %d or %d.\n"), "grand", iStrPos + 1, 2, 3);
319         return types::Function::Error;
320         }*/
321
322         std::vector<types::Double*> vectpDblTemp;
323         for (int i = 0; i < iStrPos; i++)
324         {
325             if (in[iPos]->isDouble() == false)
326             {
327                 Scierror(999, _("%s: Wrong type for input argument #%d : A scalar expected.\n"), "grand", iPos + 1);
328                 delete[] itab;
329                 return types::Function::Error;
330             }
331
332             vectpDblTemp.push_back(in[iPos]->getAs<types::Double>());
333
334             if (iStrPos == 3 && vectpDblTemp[i]->isScalar() == false)
335             {
336                 Scierror(999, _("%s: Wrong type for input argument #%d : A scalar expected.\n"), "grand", iPos + 1);
337                 delete[] itab;
338                 return types::Function::Error;
339             }
340             iPos++;
341         }
342
343         //get number of dimensions to output
344         for (int i = 0; i < iStrPos; i++)
345         {
346             itab[i] = static_cast<int>(vectpDblTemp[i]->get(0));
347         }
348     }
349
350     if (in[iPos]->isString() == false)
351     {
352         Scierror(999, _("%s: Wrong type for input argument #%d : string expected.\n"), "grand", iPos + 1);
353         delete[] itab;
354         return types::Function::Error;
355     }
356     iPos++; // method string has been already got.
357
358     // *** get arguments after methode string. ***
359     if (meth == 22 || meth == 25) // setgen || phr2sd
360     {
361         if (in[iPos]->isString() == false)
362         {
363             Scierror(999, _("%s: Wrong type for input argument #%d : string expected.\n"), "grand", iPos + 1);
364             delete[] itab;
365             return types::Function::Error;
366         }
367
368         pStrGenOrPhr = in[iPos]->getAs<types::String>();
369
370         if (pStrGenOrPhr->isScalar() == false)
371         {
372             Scierror(999, _("%s: Wrong size for input argument #%d : Only one string expected.\n"), "grand", iPos + 1);
373             delete[] itab;
374             return types::Function::Error;
375         }
376     }
377     else if (meth == 16)
378     {
379         if (in.size() > 3)
380         {
381             Scierror(999, _("Missing vect for random permutation\n"));
382             delete[] itab;
383             return types::Function::Error;
384         }
385     }
386     else
387     {
388         if (in.size() - iPos != iNumInputArg)
389         {
390             Scierror(77, _("%s: Wrong number of input argument(s): %d expected.\n"), "grand", iPos + iNumInputArg);
391             delete[] itab;
392             return types::Function::Error;
393         }
394
395         for (int i = iPos; i < in.size(); i++)
396         {
397             if (in[i]->isDouble() == false)
398             {
399                 Scierror(999, _("%s: Wrong type for input argument #%d : A scalar expected.\n"), "grand", iPos + 1);
400                 delete[] itab;
401                 return types::Function::Error;
402             }
403
404             vectpDblInput.push_back(in[i]->getAs<types::Double>());
405         }
406     }
407
408     // *** perform operation in according method string and return result. ***
409
410     types::Double* pDblOut = NULL;
411
412     if (iStrPos >= 2)
413     {
414         int iProdiTab = 1;
415         for (int i = 0; i < iStrPos; i++)
416         {
417             iProdiTab = iProdiTab * itab[i];
418         }
419
420         if (iProdiTab == 0)
421         {
422             pDblOut = types::Double::Empty();
423         }
424         else if ((itab[1] != itab[0]) && (itab[0] <= -1))
425         {
426             Scierror(999, _("%s: Wrong value for input argument #%d: Positive scalar expected.\n"), "grand", 1);
427             delete[] itab;
428             return types::Function::Error;
429         }
430         else if ((itab[1] != itab[0]) && (itab[1] <= -1))
431         {
432             Scierror(999, _("%s: Wrong value for input argument #%d: Positive scalar expected.\n"), "grand", 2);
433             delete[] itab;
434             return types::Function::Error;
435         }
436         else
437         {
438             pDblOut = new types::Double(iDims, itab);
439         }
440     }
441     else
442     {
443         types::Double* pDblIn = in[0]->getAs<types::Double>();
444         if (meth == 14 || meth == 13)//'mul' or 'markov'
445         {
446             int iDimsArraytempo[2];
447             iDimsArraytempo[0] = in[3]->getAs<types::Double>()->getSize();
448             if (meth == 14)
449             {
450                 iDimsArraytempo[0] = iDimsArraytempo[0] + 1;
451             }
452
453             iDimsArraytempo[1] = iNumIter;
454             pDblOut = new types::Double(pDblIn->getDims(), iDimsArraytempo);
455         }
456         else
457         {
458             pDblOut = new types::Double(pDblIn->getDims(), pDblIn->getDimsArray());
459         }
460
461     }
462
463     delete[] itab;
464
465     switch (meth)
466     {
467         case 1: // beta
468         {
469             double minlog   = 1.e-37;
470
471             for (int i = 0; i < 2; i++)
472             {
473                 if (vectpDblInput[i]->isScalar() == false)
474                 {
475                     delete pDblOut;
476                     Scierror(999, _("%s: Wrong type for input argument #%d : A scalar expected.\n"), "grand", i + 4);
477                     return types::Function::Error;
478                 }
479
480                 if (vectpDblInput[i]->get(0) < minlog)
481                 {
482                     delete pDblOut;
483                     Scierror(999, _("%s: Wrong value for input argument #%d : At least %lf expected.\n"), "grand", iPos + 1, minlog);
484                     return types::Function::Error;
485                 }
486             }
487
488             for (int i = 0; i < pDblOut->getSize(); i++)
489             {
490                 pDblOut->set(i, C2F(genbet)(vectpDblInput[0]->get(), vectpDblInput[1]->get()));
491             }
492
493             out.push_back(pDblOut);
494             break;
495         }
496         case 2: // binomial
497         case 3: // negative binomial
498         {
499             for (int i = 0; i < 2; i++)
500             {
501                 if (vectpDblInput[i]->isScalar() == false)
502                 {
503                     delete pDblOut;
504                     Scierror(999, _("%s: Wrong type for input argument #%d : A scalar expected.\n"), "grand", i + 4);
505                     return types::Function::Error;
506                 }
507             }
508
509             if (vectpDblInput[0]->get(0) < 0.0) // N
510             {
511                 delete pDblOut;
512                 Scierror(999, _("%s: Wrong value for input argument #%d : Positive integer expected.\n"), "grand", 4);
513                 return types::Function::Error;
514             }
515
516             if (vectpDblInput[1]->get(0) < 0.0 || vectpDblInput[1]->get(0) > 1.0) // p
517             {
518                 delete pDblOut;
519                 Scierror(999, _("%s: Wrong value for input argument #%d : A value expected.\n"), "grand", 5);
520                 return types::Function::Error;
521             }
522
523             int N = static_cast<int>(vectpDblInput[0]->get(0));
524             if (meth == 2)
525             {
526                 for (int i = 0; i < pDblOut->getSize(); i++)
527                 {
528                     pDblOut->set(i, static_cast<double>(C2F(ignbin)(&N, vectpDblInput[1]->get())));
529                 }
530             }
531             else // meth == 3
532             {
533                 for (int i = 0; i < pDblOut->getSize(); i++)
534                 {
535                     pDblOut->set(i, static_cast<double>(C2F(ignnbn)(&N, vectpDblInput[1]->get())));
536                 }
537             }
538
539             out.push_back(pDblOut);
540             break;
541         }
542         case 4: // chisquare
543         {
544             if (vectpDblInput[0]->isScalar() == false)
545             {
546                 delete pDblOut;
547                 Scierror(999, _("%s: Wrong type for input argument #%d : A scalar expected.\n"), "grand", 4);
548                 return types::Function::Error;
549             }
550
551             if (vectpDblInput[0]->get(0) <= 0.0) // Df
552             {
553                 delete pDblOut;
554                 Scierror(999, _("%s: Wrong value for input argument #%d : Positive no null value expected.\n"), "grand", 4);
555                 return types::Function::Error;
556             }
557
558             for (int i = 0; i < pDblOut->getSize(); i++)
559             {
560                 pDblOut->set(i, C2F(genchi)(vectpDblInput[0]->get()));
561             }
562
563             out.push_back(pDblOut);
564             break;
565         }
566         case 5: // non central chisquare
567         {
568             for (int i = 0; i < 2; i++)
569             {
570                 if (vectpDblInput[i]->isScalar() == false)
571                 {
572                     delete pDblOut;
573                     Scierror(999, _("%s: Wrong type for input argument #%d : A scalar expected.\n"), "grand", i + 4);
574                     return types::Function::Error;
575                 }
576             }
577
578             if (vectpDblInput[0]->get(0) < 1.0) // Df
579             {
580                 delete pDblOut;
581                 Scierror(999, _("%s: Wrong value for input argument #%d : value > 1.0 expected.\n"), "grand", 4);
582                 return types::Function::Error;
583             }
584
585             if (vectpDblInput[1]->get(0) < 0.0) // Xnon
586             {
587                 delete pDblOut;
588                 Scierror(999, _("%s: Wrong value for input argument #%d : Positive value expected.\n"), "grand", 5);
589                 return types::Function::Error;
590             }
591
592             for (int i = 0; i < pDblOut->getSize(); i++)
593             {
594                 pDblOut->set(i, C2F(gennch)(vectpDblInput[0]->get(), vectpDblInput[1]->get()));
595             }
596
597             out.push_back(pDblOut);
598             break;
599         }
600         case 6: // exponential
601         {
602             if (vectpDblInput[0]->isScalar() == false)
603             {
604                 delete pDblOut;
605                 Scierror(999, _("%s: Wrong type for input argument #%d : A scalar expected.\n"), "grand", 4);
606                 return types::Function::Error;
607             }
608
609             if (vectpDblInput[0]->get(0) < 0.0) // Av
610             {
611                 delete pDblOut;
612                 Scierror(999, _("%s: Wrong value for input argument #%d : Positive value expected.\n"), "grand", 4);
613                 return types::Function::Error;
614             }
615
616             for (int i = 0; i < pDblOut->getSize(); i++)
617             {
618                 pDblOut->set(i, C2F(genexp)(vectpDblInput[0]->get()));
619             }
620
621             out.push_back(pDblOut);
622             break;
623         }
624         case 7: // F variance ratio
625         {
626             for (int i = 0; i < 2; i++)
627             {
628                 if (vectpDblInput[i]->isScalar() == false)
629                 {
630                     delete pDblOut;
631                     Scierror(999, _("%s: Wrong type for input argument #%d : A scalar expected.\n"), "grand", i + 4);
632                     return types::Function::Error;
633                 }
634
635                 if (vectpDblInput[i]->get(0) <= 0.0) // Dfn Dfd
636                 {
637                     delete pDblOut;
638                     Scierror(999, _("%s: Wrong value for input argument #%d : Positive no null value expected.\n"), "grand", i + 4);
639                     return types::Function::Error;
640                 }
641             }
642
643             for (int i = 0; i < pDblOut->getSize(); i++)
644             {
645                 pDblOut->set(i, C2F(genf)(vectpDblInput[0]->get(), vectpDblInput[1]->get()));
646             }
647
648             out.push_back(pDblOut);
649             break;
650         }
651         case 8: // non F variance ratio
652         {
653             for (int i = 0; i < 3; i++)
654             {
655                 if (vectpDblInput[i]->isScalar() == false)
656                 {
657                     delete pDblOut;
658                     Scierror(999, _("%s: Wrong type for input argument #%d : A scalar expected.\n"), "grand", i + 4);
659                     return types::Function::Error;
660                 }
661             }
662
663             if (vectpDblInput[0]->get(0) < 1.0) // Dfn
664             {
665                 delete pDblOut;
666                 Scierror(999, _("%s: Wrong value for input argument #%d : value > 1.0 expected.\n"), "grand", 4);
667                 return types::Function::Error;
668             }
669
670             if (vectpDblInput[1]->get(0) <= 0.0) // Dfd
671             {
672                 delete pDblOut;
673                 Scierror(999, _("%s: Wrong value for input argument #%d : Positive non null value expected.\n"), "grand", 5);
674                 return types::Function::Error;
675             }
676
677             if (vectpDblInput[2]->get(0) < 0.0) // Xnon
678             {
679                 delete pDblOut;
680                 Scierror(999, _("%s: Wrong value for input argument #%d : Positive value expected.\n"), "grand", 6);
681                 return types::Function::Error;
682             }
683
684             for (int i = 0; i < pDblOut->getSize(); i++)
685             {
686                 pDblOut->set(i, C2F(gennf)(vectpDblInput[0]->get(), vectpDblInput[1]->get(), vectpDblInput[2]->get()));
687             }
688
689             out.push_back(pDblOut);
690             break;
691         }
692         case 9: // gamma
693         {
694             for (int i = 0; i < 2; i++)
695             {
696                 if (vectpDblInput[i]->isScalar() == false)
697                 {
698                     delete pDblOut;
699                     Scierror(999, _("%s: Wrong type for input argument #%d : A scalar expected.\n"), "grand", i + 4);
700                     return types::Function::Error;
701                 }
702
703                 if (vectpDblInput[i]->get(0) <= 0.0)
704                 {
705                     delete pDblOut;
706                     Scierror(999, _("%s: Wrong value for input argument #%d : Positive non null value expected.\n"), "grand", i + 4);
707                     return types::Function::Error;
708                 }
709             }
710
711             for (int i = 0; i < pDblOut->getSize(); i++)
712             {
713                 //** WARNING ** : order is changed in parameters for
714                 //compatibility between Rand(...'gam',..) and cdfgam
715
716                 pDblOut->set(i, C2F(gengam)(vectpDblInput[1]->get(), vectpDblInput[0]->get()));
717             }
718
719             out.push_back(pDblOut);
720             break;
721         }
722         case 10: // Gauss Laplace (normal)
723         {
724             for (int i = 0; i < 2; i++)
725             {
726                 if (vectpDblInput[i]->isScalar() == false)
727                 {
728                     delete pDblOut;
729                     Scierror(999, _("%s: Wrong type for input argument #%d : A scalar expected.\n"), "grand", i + 4);
730                     return types::Function::Error;
731                 }
732             }
733
734             if (vectpDblInput[1]->get(0) < 0.0)
735             {
736                 delete pDblOut;
737                 Scierror(999, _("%s: Wrong value for input argument #%d : Positive value expected.\n"), "grand", 5);
738                 return types::Function::Error;
739             }
740
741             for (int i = 0; i < pDblOut->getSize(); i++)
742             {
743                 pDblOut->set(i, C2F(gennor)(vectpDblInput[0]->get(), vectpDblInput[1]->get()));
744             }
745
746             out.push_back(pDblOut);
747             break;
748         }
749         case 11: // multivariate gaussian (multivariate normal)
750         {
751             if (vectpDblInput[0]->getCols() != 1 || vectpDblInput[0]->getSize() == 0)
752             {
753                 delete pDblOut;
754                 Scierror(999, _("%s: Wrong size for input argument #%d : A matrix of size m x 1 expected.(m > 0)\n"), "grand", 4);
755                 return types::Function::Error;
756             }
757
758             if (vectpDblInput[0]->getRows() != vectpDblInput[1]->getRows())
759             {
760                 delete pDblOut;
761                 Scierror(999, _("%s: Wrong size for input arguments #%d and #%d: Mean and Cov have incompatible dimensions.\n"), "grand", 4, 5);
762                 return types::Function::Error;
763             }
764
765             if (vectpDblInput[1]->getRows() != vectpDblInput[1]->getCols())
766             {
767                 delete pDblOut;
768                 Scierror(999, _("%s: Wrong size for input argument #%d : A square symmetric positive definite matrix expected.\n"), "grand", 5);
769                 return types::Function::Error;
770             }
771
772             int ierr = 0;
773             int size = vectpDblInput[0]->getRows();
774             int mp   = (size * (size + 3)) / 2 + 1;
775
776             delete pDblOut;
777             types::Double* pDblOut = new types::Double(size, iNumIter);
778             double* work = new double[size];
779             double* param = new double[mp];
780
781             types::Double* pDblMean = vectpDblInput[0]->clone()->getAs<types::Double>();
782             types::Double* pDblCov  = vectpDblInput[1]->clone()->getAs<types::Double>();
783
784             C2F(setgmn)(pDblMean->get(), pDblCov->get(), &size, &size, param, &ierr);
785
786             delete pDblMean;
787             delete pDblCov;
788
789             if (ierr == 1)
790             {
791                 delete[] work;
792                 delete[] param;
793                 delete pDblOut;
794                 Scierror(999, _("%s: setgmn return with state %d.\n"), "grand", ierr);
795                 return types::Function::Error;
796             }
797
798             for (int i = 0; i < iNumIter; i++)
799             {
800                 C2F(genmn)(param, pDblOut->get() + (size * i), work);
801             }
802
803             delete[] work;
804             delete[] param;
805
806             out.push_back(pDblOut);
807             break;
808         }
809         case 12: // geometric
810         {
811             double pmin = 1.3e-307;
812
813             if (vectpDblInput[0]->isScalar() == false)
814             {
815                 delete pDblOut;
816                 Scierror(999, _("%s: Wrong type for input argument #%d : A scalar expected.\n"), "grand", 4);
817                 return types::Function::Error;
818             }
819
820             if (vectpDblInput[0]->get(0) < pmin || vectpDblInput[0]->get(0) > 1.0)
821             {
822                 delete pDblOut;
823                 Scierror(999, _("%s: Wrong value for input argument #%d : Must be between %lf and %d.\n"), "grand", 4, pmin, 1);
824                 return types::Function::Error;
825             }
826
827             for (int i = 0; i < pDblOut->getSize(); i++)
828             {
829                 pDblOut->set(i, igngeom(vectpDblInput[0]->get(0)));
830             }
831
832             out.push_back(pDblOut);
833             break;
834         }
835         case 13: // markov
836         {
837             if ( vectpDblInput[0]->getRows() != vectpDblInput[0]->getCols() &&
838                     vectpDblInput[0]->getRows() != 1)
839             {
840                 delete pDblOut;
841                 Scierror(999, _("%s: Wrong size for input argument #%d : A square matrix or a row vector expected.\n"), "grand", 4);
842                 return types::Function::Error;
843             }
844
845             if (vectpDblInput[1]->getSize() == 0)
846             {
847                 delete pDblOut;
848                 Scierror(999, _("%s: Wrong size for input argument #%d: No empty matrix expected.\n"), "grand", 5);
849                 return types::Function::Error;
850             }
851
852             int sizeOfX0 = vectpDblInput[1]->getSize();
853
854             for ( int i = 0; i < sizeOfX0; i++)
855             {
856                 if (vectpDblInput[1]->get(i) < 1 || (vectpDblInput[1]->get(i) - 1) >= vectpDblInput[0]->getCols())
857                 {
858                     Scierror(999, _("%s: X0(%d) must be in the range [1,%d[.\n"), "grand", i + 1, vectpDblInput[0]->getCols() + 1);
859                     return types::Function::Error;
860                 }
861             }
862
863             int size = vectpDblInput[0]->getSize() + vectpDblInput[0]->getCols();
864             double* work = new double[(size * sizeof(double))];
865
866             for (int i = 0; i < vectpDblInput[0]->getRows(); i++)
867             {
868                 double ptot = 0.0;
869                 for (int j = 0; j < vectpDblInput[0]->getCols(); j++)
870                 {
871                     int position = vectpDblInput[0]->getRows() * j + i;
872
873                     if (vectpDblInput[0]->get(position) < 0 || vectpDblInput[0]->get(position) > 1)
874                     {
875                         delete pDblOut;
876                         Scierror(999, _("%s: Wrong value for input argument #%d: P(%d,%d) must be in the range [0 1].\n"), "grand", i + 1, j + 1);
877                         delete[] work;
878                         return types::Function::Error;
879                     }
880
881                     ptot += vectpDblInput[0]->get(position);
882                 }
883
884                 if (fabs(ptot - 1.0) > 1e-8)
885                 {
886                     delete pDblOut;
887                     Scierror(999, _("%s: Sum of P(%d,1:%d)=%lf ~= 1.\n"), "grand", i + 1, vectpDblInput[0]->getCols(), ptot);
888                     delete[] work;
889                     return types::Function::Error;
890                 }
891             }
892             // ** Computing the cumulative sum of the P matrix **
893
894             for (int i = 0; i < vectpDblInput[0]->getRows(); i++)
895             {
896                 double cumsum = 0.0;
897                 work[i] = 0.0;
898                 for (int j = 1; j < vectpDblInput[0]->getCols() + 1; j++)
899                 {
900                     cumsum += vectpDblInput[0]->get(vectpDblInput[0]->getRows() * (j - 1) + i);
901                     work[vectpDblInput[0]->getRows() * j + i] = cumsum;
902                 }
903             }
904
905             for (int i = 0; i < vectpDblInput[1]->getSize(); i++)
906             {
907                 int iCur = static_cast<int>(vectpDblInput[1]->get(i)) - 1;
908                 for (int j = 0; j < iNumIter; j++)
909                 {
910                     int niv = 0;
911                     double rr = C2F(ranf)();
912
913                     if (vectpDblInput[0]->getRows() == 1)
914                     {
915                         iCur = 0;
916                     }
917
918                     while (  rr >= work[iCur + vectpDblInput[0]->getRows() * niv] &&
919                              niv < (vectpDblInput[0]->getRows() + 1))
920                     {
921                         niv++;
922                     }
923
924                     // ** projection to avoid boundaries **
925                     niv = std::max(std::min(niv, vectpDblInput[0]->getCols()), 1);
926                     pDblOut->set(vectpDblInput[1]->getSize() * j + i, static_cast<double>(niv));
927                     iCur = niv - 1;
928                 }
929             }
930
931             out.push_back(pDblOut);
932             delete[] work;
933             break;
934         }
935         case 14: // multinomial
936         {
937             if (vectpDblInput[0]->isScalar() == false)
938             {
939                 delete pDblOut;
940                 Scierror(999, _("%s: Wrong size for input argument #%d: A scalar expected.\n"), "grand", 4);
941                 return types::Function::Error;
942             }
943
944             if (vectpDblInput[0]->get(0) < 0)
945             {
946                 delete pDblOut;
947                 Scierror(999, _("%s: Wrong value for input argument #%d: A positive scalar expected.\n"), "grand", 4);
948                 return types::Function::Error;
949             }
950
951             if (vectpDblInput[1]->getCols() != 1 || vectpDblInput[1]->getRows() <= 0)
952             {
953                 delete pDblOut;
954                 Scierror(999, _("%s: Wrong size for input argument #%d: A column vector expected.\n"), "grand", 5);
955                 return types::Function::Error;
956             }
957
958             double ptot = 0.0;
959             int ncat = vectpDblInput[1]->getRows() + 1;
960
961             for (int i = 0; i < vectpDblInput[1]->getCols(); i++)
962             {
963                 if (vectpDblInput[1]->get(i) < 0.0 || vectpDblInput[1]->get(i) > 1.0)
964                 {
965                     delete pDblOut;
966                     Scierror(999, _("%s: Wrong value for input argument #%d: P(%d) must be in the range [0 1].\n"), "grand", 5, i + 1);
967                     return types::Function::Error;
968                 }
969                 ptot += vectpDblInput[1]->get(i);
970             }
971
972             if (ptot > 1.0)
973             {
974                 delete pDblOut;
975                 Scierror(999, _("%s: Wrong value for input argument #%d: Sum of P(i) > 1.\n"), "grand", 5);
976                 return types::Function::Error;
977             }
978
979             int* piP = new int[vectpDblInput[0]->getSize()];
980             int* piOut = new int[pDblOut->getSize()];
981
982             for (int i = 0; i < vectpDblInput[0]->getSize(); i++)
983             {
984                 piP[i] = static_cast<int>(vectpDblInput[0]->get(i));
985             }
986
987             for (int i = 0; i < iNumIter; i++)
988             {
989                 C2F(genmul)(piP, vectpDblInput[1]->get(), &ncat,  + &piOut[ncat * i]);
990             }
991
992             for (int i = 0; i < pDblOut->getSize(); i++)
993             {
994                 pDblOut->set(i, static_cast<double>(piOut[i]));
995             }
996
997             delete piP;
998             delete piOut;
999
1000             out.push_back(pDblOut);
1001             break;
1002         }
1003         case 15: // Poisson
1004         {
1005             if (vectpDblInput[0]->isScalar() == false)
1006             {
1007                 delete pDblOut;
1008                 Scierror(999, _("%s: Wrong type for input argument #%d : A scalar expected.\n"), "grand", 4);
1009                 return types::Function::Error;
1010             }
1011
1012             if (vectpDblInput[0]->get(0) < 0.0)
1013             {
1014                 delete pDblOut;
1015                 Scierror(999, _("%s: Wrong value for input argument #%d : Positive value expected.\n"), "grand", 4);
1016                 return types::Function::Error;
1017             }
1018
1019             for (int i = 0; i < pDblOut->getSize(); i++)
1020             {
1021                 pDblOut->set(i, static_cast<double>(C2F(ignpoi)(vectpDblInput[0]->get())));
1022             }
1023
1024             out.push_back(pDblOut);
1025             break;
1026         }
1027         case 16: // random permutations
1028         {
1029             delete pDblOut;
1030             types::InternalType* pITOut = NULL;
1031
1032             // get dims and isComplex
1033             int iDims = 0;
1034             int* piDimsArray = NULL;
1035             bool bIsComplex = false;
1036             if (in[2]->isGenericType())
1037             {
1038                 types::GenericType* pGT = in[2]->getAs<types::GenericType>();
1039                 iDims = pGT->getDims();
1040                 piDimsArray = pGT->getDimsArray();
1041                 bIsComplex = pGT->isComplex();
1042             }
1043
1044             switch (in[2]->getType())
1045             {
1046                 case types::InternalType::ScilabInt8:
1047                     pITOut = new types::Int8(iDims, piDimsArray);
1048                     sci_grand_prm(iNumIter, in[2]->getAs<types::Int8>(), &pITOut);
1049                     break;
1050                 case types::InternalType::ScilabUInt8:
1051                     pITOut = new types::UInt8(iDims, piDimsArray);
1052                     sci_grand_prm(iNumIter, in[2]->getAs<types::UInt8>(), &pITOut);
1053                     break;
1054                 case types::InternalType::ScilabInt16:
1055                     pITOut = new types::Int16(iDims, piDimsArray);
1056                     sci_grand_prm(iNumIter, in[2]->getAs<types::Int16>(), &pITOut);
1057                     break;
1058                 case types::InternalType::ScilabUInt16:
1059                     pITOut = new types::UInt16(iDims, piDimsArray);
1060                     sci_grand_prm(iNumIter, in[2]->getAs<types::UInt16>(), &pITOut);
1061                     break;
1062                 case types::InternalType::ScilabInt32:
1063                     pITOut = new types::Int32(iDims, piDimsArray);
1064                     sci_grand_prm(iNumIter, in[2]->getAs<types::Int32>(), &pITOut);
1065                     break;
1066                 case types::InternalType::ScilabUInt32:
1067                     pITOut = new types::UInt32(iDims, piDimsArray);
1068                     sci_grand_prm(iNumIter, in[2]->getAs<types::UInt32>(), &pITOut);
1069                     break;
1070                 case types::InternalType::ScilabInt64:
1071                     pITOut = new types::Int64(iDims, piDimsArray);
1072                     sci_grand_prm(iNumIter, in[2]->getAs<types::Int64>(), &pITOut);
1073                     break;
1074                 case types::InternalType::ScilabUInt64:
1075                     pITOut = new types::UInt64(iDims, piDimsArray);
1076                     sci_grand_prm(iNumIter, in[2]->getAs<types::UInt64>(), &pITOut);
1077                     break;
1078                 case types::InternalType::ScilabDouble:
1079                     pITOut = new types::Double(iDims, piDimsArray, bIsComplex);
1080                     sci_grand_prm(iNumIter, in[2]->getAs<types::Double>(), &pITOut);
1081                     break;
1082                 case types::InternalType::ScilabBool:
1083                     pITOut = new types::Bool(iDims, piDimsArray);
1084                     sci_grand_prm(iNumIter, in[2]->getAs<types::Bool>(), &pITOut);
1085                     break;
1086                 case types::InternalType::ScilabString:
1087                     pITOut = new types::String(iDims, piDimsArray);
1088                     sci_grand_prm(iNumIter, in[2]->getAs<types::String>(), &pITOut);
1089                     break;
1090                 case types::InternalType::ScilabPolynom:
1091                     pITOut = new types::Polynom(in[2]->getAs<types::Polynom>()->getVariableName(), iDims, piDimsArray);
1092                     sci_grand_prm(iNumIter, in[2]->getAs<types::Polynom>(), &pITOut);
1093                     break;
1094                 case types::InternalType::ScilabSparse:
1095                 {
1096                     std::complex<double> cplxDbl;
1097                     types::InternalType* pITOutTempo = NULL;
1098                     types::Double* pDblTempo = NULL;
1099                     types::Sparse* pSP = in[2]->getAs<types::Sparse>();
1100                     int isize = pSP->getSize();
1101                     pITOut = new types::Sparse(pSP->getRows(), pSP->getCols(), pSP->isComplex());
1102                     pDblTempo = new types::Double(isize, 1, pSP->isComplex());
1103                     pITOutTempo = new types::Double(isize, iNumIter, pSP->isComplex());
1104
1105                     if (pDblTempo->isComplex())
1106                     {
1107                         for (int i = 0; i < isize; i++)
1108                         {
1109                             cplxDbl = in[2]->getAs<types::Sparse>()->getImg(i);
1110                             pDblTempo->set(i, cplxDbl.real());
1111                             pDblTempo->set(i, cplxDbl.imag());
1112                         }
1113                     }
1114                     else
1115                     {
1116                         for (int i = 0; i < isize; i++)
1117                         {
1118                             pDblTempo->set(i, in[2]->getAs<types::Sparse>()->get(i));
1119                         }
1120                     }
1121                     sci_grand_prm(iNumIter, pDblTempo, &pITOutTempo);
1122
1123                     if (pDblTempo->isComplex())
1124                     {
1125                         for (int i = 0; i < isize; i++)
1126                         {
1127                             cplxDbl.real(pITOutTempo->getAs<types::Double>()->get(i));
1128                             cplxDbl.imag(pITOutTempo->getAs<types::Double>()->getImg(i));
1129                             pITOutTempo->getAs<types::Sparse>()->set(i, cplxDbl);
1130                         }
1131                     }
1132                     else
1133                     {
1134                         for (int i = 0; i < isize; i++)
1135                         {
1136                             pITOut->getAs<types::Sparse>()->set(i, pITOutTempo->getAs<types::Double>()->get(i));
1137                         }
1138                     }
1139
1140                     delete pITOutTempo;
1141                     delete pDblTempo;
1142                     break;
1143                 }
1144                 default:
1145                 {
1146                     Scierror(999, _("%s: Wrong type for input argument: Matrix (full or sparse) or Hypermatrix of Reals, Complexes, Integers, Booleans, Strings or Polynomials expected.\n"), "grand");
1147                     return types::Function::Error;
1148                 }
1149             }
1150
1151             out.push_back(pITOut);
1152             break;
1153         }
1154         case 17: // uniform (def)
1155         {
1156             for (int i = 0; i < pDblOut->getSize(); i++)
1157             {
1158                 pDblOut->set(i, C2F(ranf)());
1159             }
1160
1161             out.push_back(pDblOut);
1162             break;
1163         }
1164         case 18: // uniform (unf)
1165         {
1166             for (int i = 0; i < 2; i++)
1167             {
1168                 if (vectpDblInput[i]->isScalar() == false)
1169                 {
1170                     delete pDblOut;
1171                     Scierror(999, _("%s: Wrong type for input argument #%d : A scalar expected.\n"), "grand", i + 4);
1172                     return types::Function::Error;
1173                 }
1174             }
1175
1176             double low  = vectpDblInput[0]->get(0);
1177             double high = vectpDblInput[1]->get(0);
1178
1179             if (low > high)
1180             {
1181                 delete pDblOut;
1182                 Scierror(999, _("%s: Wrong value for input arguments #%d and #%d: Low < High expected.\n"), "grand", 4, 5);
1183                 return types::Function::Error;
1184             }
1185
1186             for (int i = 0; i < pDblOut->getSize(); i++)
1187             {
1188                 pDblOut->set(i, low + (high - low) * C2F(ranf)());
1189             }
1190
1191             out.push_back(pDblOut);
1192             break;
1193         }
1194         case 19: // uniform (uin)
1195         {
1196             for (int i = 0; i < 2; i++)
1197             {
1198                 if (vectpDblInput[i]->isScalar() == false)
1199                 {
1200                     delete pDblOut;
1201                     Scierror(999, _("%s: Wrong type for input argument #%d : A scalar expected.\n"), "grand", i + 4);
1202                     return types::Function::Error;
1203                 }
1204             }
1205
1206             int low  = static_cast<int>(vectpDblInput[0]->get(0));
1207             int high = static_cast<int>(vectpDblInput[1]->get(0));
1208
1209             if (low > high)
1210             {
1211                 delete pDblOut;
1212                 Scierror(999, _("%s: Wrong value for input arguments #%d and #%d: Low < High expected.\n"), "grand", 4, 5);
1213                 return types::Function::Error;
1214             }
1215
1216             if ( low  != vectpDblInput[0]->get(0) ||
1217                     high != vectpDblInput[1]->get(0) ||
1218                     (high - low + 1) > 2147483561)
1219             {
1220                 delete pDblOut;
1221                 Scierror(999, _("%s: Wrong value for input arguments #%d and #%d: Low and High must be integers and (high - low + 1) <=  2147483561.\n"), "grand", 4, 5);
1222                 return types::Function::Error;
1223             }
1224
1225             for (int i = 0; i < pDblOut->getSize(); i++)
1226             {
1227                 pDblOut->set(i, C2F(ignuin)(vectpDblInput[0]->get(), vectpDblInput[1]->get()));
1228             }
1229
1230             out.push_back(pDblOut);
1231             break;
1232         }
1233         case 20: // uniform (lgi)
1234         {
1235             for (int i = 0; i < pDblOut->getSize(); i++)
1236             {
1237                 pDblOut->set(i, static_cast<double>(ignlgi()));
1238             }
1239
1240             out.push_back(pDblOut);
1241             break;
1242         }
1243         case 21: // getgen
1244         {
1245             delete pDblOut;
1246             types::String* pStrOut = new types::String(names_gen[current_base_gen]);
1247             out.push_back(pStrOut);
1248             break;
1249         }
1250         case 22: // setgen
1251         {
1252             delete pDblOut;
1253             wchar_t* wcsGen = pStrGenOrPhr->get(0);
1254
1255             if (wcscmp(wcsGen, L"mt") == 0)
1256             {
1257                 ConfigVariable::setCurrentBaseGen(MT);
1258             }
1259             else if (wcscmp(wcsGen, L"kiss") == 0)
1260             {
1261                 ConfigVariable::setCurrentBaseGen(KISS);
1262             }
1263             else if (wcscmp(wcsGen, L"clcg4") == 0)
1264             {
1265                 ConfigVariable::setCurrentBaseGen(CLCG4);
1266             }
1267             else if (wcscmp(wcsGen, L"clcg2") == 0)
1268             {
1269                 ConfigVariable::setCurrentBaseGen(CLCG2);
1270             }
1271             else if (wcscmp(wcsGen, L"urand") == 0)
1272             {
1273                 ConfigVariable::setCurrentBaseGen(URAND);
1274             }
1275             else
1276             {
1277                 Scierror(999, _("%s: Wrong value for input argument #%d: '%s', '%s', '%s', '%s', '%s' or '%s' expected.\n"), "grand", 2, "mt", "kiss", "clcg4", "clcg2", "urand");
1278                 return types::Function::Error;
1279             }
1280
1281             int current_gen = ConfigVariable::getCurrentBaseGen();
1282             types::String* pStrOut = new types::String(names_gen[current_gen]);
1283             out.push_back(pStrOut);
1284             break;
1285         }
1286         case 23: // getsd
1287         {
1288             switch (current_base_gen)
1289             {
1290                 case MT:
1291                 {
1292                     pDblOut = new types::Double(625, 1);
1293                     get_state_mt(pDblOut->get());
1294                     break;
1295                 }
1296                 case KISS:
1297                 {
1298                     pDblOut = new types::Double(4, 1);
1299                     get_state_kiss(pDblOut->get());
1300                     break;
1301                 }
1302                 case CLCG4:
1303                 {
1304                     pDblOut = new types::Double(4, 1);
1305                     int current_clcg4 = ConfigVariable::getCurrentClcg4();
1306                     get_state_clcg4(current_clcg4, pDblOut->get());
1307                     break;
1308                 }
1309                 case CLCG2:
1310                 {
1311                     pDblOut = new types::Double(2, 1);
1312                     get_state_clcg2(pDblOut->get());
1313                     break;
1314                 }
1315                 case URAND:
1316                 {
1317                     pDblOut = new types::Double(1, 1);
1318                     get_state_urand(pDblOut->get());
1319                     break;
1320                 }
1321             }
1322
1323             out.push_back(pDblOut);
1324             break;
1325         }
1326         case 24: // setsd
1327         {
1328             delete pDblOut;
1329             int ierr = 0;
1330             switch (current_base_gen)
1331             {
1332                 case MT:
1333                 {
1334                     if (vectpDblInput[0]->isScalar())
1335                     {
1336                         ierr = set_state_mt_simple(vectpDblInput[0]->get(0));
1337                     }
1338                     else if (vectpDblInput[0]->getSize() == 625)
1339                     {
1340                         ierr = set_state_mt(vectpDblInput[0]->get());
1341                     }
1342                     else
1343                     {
1344                         Scierror(999, _("%s: Wrong size for input argument #%d : A scalar or a vector of size %d expected.\n"), "grand", 4, 625);
1345                         return types::Function::Error;
1346                     }
1347
1348                     break;
1349                 }
1350                 case KISS:
1351                 case CLCG4:
1352                 {
1353                     for (int i = 0; i < 4; i++)
1354                     {
1355                         if (vectpDblInput[i]->isScalar() == false)
1356                         {
1357                             Scierror(999, _("%s: Wrong type for input argument #%d : A scalar expected.\n"), "grand", i + 2);
1358                             return types::Function::Error;
1359                         }
1360                     }
1361
1362                     if (current_base_gen == KISS)
1363                     {
1364                         ierr = set_state_kiss(  vectpDblInput[0]->get(0),
1365                                                 vectpDblInput[1]->get(0),
1366                                                 vectpDblInput[2]->get(0),
1367                                                 vectpDblInput[3]->get(0));
1368                     }
1369                     else // CLCG4
1370                     {
1371                         ierr = set_seed_clcg4(  ConfigVariable::getCurrentClcg4(),
1372                                                 vectpDblInput[0]->get(0),
1373                                                 vectpDblInput[1]->get(0),
1374                                                 vectpDblInput[2]->get(0),
1375                                                 vectpDblInput[3]->get(0));
1376                     }
1377
1378                     break;
1379                 }
1380                 case CLCG2:
1381                 {
1382                     for (int i = 0; i < 2; i++)
1383                     {
1384                         if (vectpDblInput[i]->isScalar() == false)
1385                         {
1386                             Scierror(999, _("%s: Wrong type for input argument #%d : A scalar expected.\n"), "grand", i + 2);
1387                             return types::Function::Error;
1388                         }
1389                     }
1390
1391                     ierr = set_state_clcg2(vectpDblInput[0]->get(0), vectpDblInput[1]->get(0));
1392                     break;
1393                 }
1394                 case URAND:
1395                 {
1396                     if (vectpDblInput[0]->isScalar() == false)
1397                     {
1398                         Scierror(999, _("%s: Wrong type for input argument #%d : A scalar expected.\n"), "grand", 2);
1399                         return types::Function::Error;
1400                     }
1401
1402                     ierr = set_state_urand(vectpDblInput[0]->get(0));
1403                     break;
1404                 }
1405             }
1406
1407             if (ierr == 0)
1408             {
1409                 Scierror(999, _("%s: Wrong value for the last %d input argument(s).\n"), "grand", in.size() - 1);
1410                 return types::Function::Error;
1411             }
1412
1413             break;
1414         }
1415         case 25: // phr2sd
1416         {
1417             delete pDblOut;
1418             if (pStrGenOrPhr->isScalar() == false)
1419             {
1420                 Scierror(999, _("%s: Wrong type for input argument #%d : One string expected.\n"), "grand", 2);
1421                 return types::Function::Error;
1422             }
1423
1424             types::Double* pDblOut = new types::Double(1, 2);
1425             int size = (int)wcslen(pStrGenOrPhr->get(0));
1426             int piOut[2];
1427             char* strPhr = wide_string_to_UTF8(pStrGenOrPhr->get(0));
1428
1429             C2F(phrtsd)(strPhr, &size, piOut, piOut + 1, size);
1430
1431             pDblOut->set(0, static_cast<double>(piOut[0]));
1432             pDblOut->set(1, static_cast<double>(piOut[1]));
1433
1434             FREE(strPhr);
1435             out.push_back(pDblOut);
1436             break;
1437         }
1438         case 26: // setcgn
1439         {
1440             delete pDblOut;
1441             if (current_base_gen != CLCG4)
1442             {
1443                 sciprint(_("The %s option affects only the %s generator\n"), "setcgn", "clcg4");
1444             }
1445
1446             if (vectpDblInput[0]->isScalar() == false)
1447             {
1448                 Scierror(999, _("%s: Wrong type for input argument #%d : A scalar expected.\n"), "grand", 2);
1449                 return types::Function::Error;
1450             }
1451
1452             if (vectpDblInput[0]->get(0) < 0 || vectpDblInput[0]->get(0) > Maxgen)
1453             {
1454                 Scierror(999, _("%s: Wrong value for input argument #%d : Must be between %d and %d.\n"), "grand", 0, Maxgen);
1455                 return types::Function::Error;
1456             }
1457
1458             ConfigVariable::setCurrentClcg4(static_cast<int>(vectpDblInput[0]->get(0)));
1459             double dOut = static_cast<double>(ConfigVariable::getCurrentClcg4());
1460             out.push_back(new types::Double(dOut));
1461             break;
1462         }
1463         case 27: // getcgn
1464         {
1465             delete pDblOut;
1466             double dOut = static_cast<double>(ConfigVariable::getCurrentClcg4());
1467             out.push_back(new types::Double(dOut));
1468             break;
1469         }
1470         case 28: // initgn
1471         {
1472             delete pDblOut;
1473             SeedType where;
1474             if (current_base_gen != CLCG4)
1475             {
1476                 sciprint(_("The %s option affects only the %s generator\n"), "initgn", "clcg4");
1477             }
1478
1479             if (vectpDblInput[0]->isScalar() == false)
1480             {
1481                 Scierror(999, _("%s: Wrong type for input argument #%d : A scalar expected.\n"), "grand", 2);
1482                 return types::Function::Error;
1483             }
1484
1485             if ( vectpDblInput[0]->get(0) != 0 &&
1486                     vectpDblInput[0]->get(0) != -1 &&
1487                     vectpDblInput[0]->get(0) != 1)
1488             {
1489                 Scierror(999, _("%s: Wrong value for input argument #%d : Must be between %d, %d or %d.\n"), "grand", 2, -1, 0, 1);
1490                 return types::Function::Error;
1491             }
1492
1493             where = (SeedType)(int)(vectpDblInput[0]->get(0) + 1);
1494             init_generator_clcg4(ConfigVariable::getCurrentClcg4(), where);
1495             out.push_back(vectpDblInput[0]);
1496             break;
1497         }
1498         case 29: // setall
1499         {
1500             delete pDblOut;
1501             if (current_base_gen != CLCG4)
1502             {
1503                 sciprint(_("The %s option affects only the %s generator\n"), "setall", "clcg4");
1504             }
1505
1506             for (int i = 0; i < 4; i++)
1507             {
1508                 if (vectpDblInput[i]->isScalar() == false)
1509                 {
1510                     Scierror(999, _("%s: Wrong type for input argument #%d : A scalar expected.\n"), "grand", i + 2);
1511                     return types::Function::Error;
1512                 }
1513             }
1514
1515             int ierr = set_initial_seed_clcg4(  vectpDblInput[0]->get(0),
1516                                                 vectpDblInput[1]->get(0),
1517                                                 vectpDblInput[2]->get(0),
1518                                                 vectpDblInput[3]->get(0));
1519             if (ierr == 0)
1520             {
1521                 Scierror(999, _("%s: Wrong value for the last %d input argument(s).\n"), "grand", 4);
1522                 return types::Function::Error;
1523             }
1524
1525             out.push_back(pStrMethod);
1526             break;
1527         }
1528         case 30: // advnst
1529         {
1530             delete pDblOut;
1531             if (current_base_gen != CLCG4)
1532             {
1533                 sciprint(_("The %s option affects only the %s generator\n"), "advnst", "clcg4");
1534             }
1535
1536             if (vectpDblInput[0]->isScalar() == false)
1537             {
1538                 Scierror(999, _("%s: Wrong type for input argument #%d : A scalar expected.\n"), "grand", 2);
1539                 return types::Function::Error;
1540             }
1541
1542             int k = static_cast<int>(vectpDblInput[0]->get(0));
1543
1544             if (k < 1)
1545             {
1546                 Scierror(999, _("%s: Wrong value for input argument #%d : Must be > %d.\n"), "grand", 2, 0);
1547                 return types::Function::Error;
1548             }
1549
1550             advance_state_clcg4(ConfigVariable::getCurrentClcg4(), k);
1551             out.push_back(new types::Double(static_cast<double>(k)));
1552             break;
1553         }
1554     }
1555
1556     return types::Function::OK;
1557 }
1558 /*--------------------------------------------------------------------------*/
1559
1560 template<class U>
1561 void sci_grand_prm(int iNumIter, U *pIn, types::InternalType** pOut)
1562 {
1563     U* pUTempo = NULL;
1564     types::InternalType* pITTempo = NULL;
1565     int* piDimsArray = NULL;
1566     int Dims = 0;
1567
1568     if ((pIn->getCols() == 1) && (pIn->getDims() == 2))
1569     {
1570         pOut[0]->getAs<U>()->resize(pIn->getRows(), iNumIter);
1571         pUTempo = pIn;
1572     }
1573     else if ((pIn->getRows() == 1) && (pIn->getDims() == 2))
1574     {
1575         pIn->transpose(pITTempo);
1576         pOut[0]->getAs<U>()->resize(iNumIter, pIn->getCols());
1577         pUTempo = pITTempo->getAs<U>();
1578     }
1579     else
1580     {
1581         piDimsArray = pOut[0]->getAs<U>()->getDimsArray();
1582         Dims = pOut[0]->getAs<U>()->getDims();
1583         piDimsArray[Dims] = iNumIter;
1584         pOut[0]->getAs<U>()->resize(piDimsArray, Dims + 1);
1585         pUTempo = pIn;
1586     }
1587
1588     int isize = pUTempo->getSize();
1589
1590     types::Double* pDblOut = new types::Double(isize, iNumIter, pUTempo->isComplex());
1591
1592     for (int i = 0; i < iNumIter; i++)
1593     {
1594         for (int j = 0; j < isize; j++)
1595         {
1596             pDblOut->set(isize * i + j, j);
1597         }
1598         C2F(genprm)(pDblOut->get() + (isize * i), &isize);
1599     }
1600
1601     if ((pIn->getCols() != 1) && (pIn->getRows() == 1) && (pIn->getDims() == 2))
1602     {
1603         pDblOut->transpose(pITTempo);
1604         delete pDblOut;
1605         pDblOut = pITTempo->getAs<types::Double>();
1606     }
1607
1608     if (pUTempo->isComplex() && pUTempo->isPoly() == false)
1609     {
1610         for (int i = 0; i < pOut[0]->getAs<U>()->getSize(); i++)
1611         {
1612             pOut[0]->getAs<U>()->set(i, pIn->get(static_cast<int>(pDblOut->get(i))));
1613             pOut[0]->getAs<U>()->setImg(i, pIn->getImg(static_cast<int>(pDblOut->get(i))));
1614         }
1615     }
1616     else
1617     {
1618         for (int i = 0; i < pOut[0]->getAs<U>()->getSize(); i++)
1619         {
1620             pOut[0]->getAs<U>()->set(i, pIn->get(static_cast<int>(pDblOut->get(i))));
1621         }
1622     }
1623
1624     if ((pIn->getCols() != 1) && (pIn->getRows() == 1) && (pIn->getDims() == 2))
1625     {
1626         delete pUTempo;
1627     }
1628
1629     delete pDblOut;
1630 }
1631
1632 /*--------------------------------------------------------------------------*/