fb6ac96ce2bb318c7655dba19caaadfb79f28cdd
[scilab.git] / scilab / modules / randlib / sci_gateway / c / sci_grand.c
1 /*
2 * Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
3 * Copyright (C) ENPC
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 /*------------------------------------------------------------------------
14 *    Interface for grand
15 *    jpc@cermics.enpc.fr
16 *    stuff to deal with several generators added
17 *         by Bruno Pincon (12/11/2001)
18 *
19 --------------------------------------------------------------------------*/
20 #include <string.h>
21 #include <math.h>
22 #include "localization.h"
23 #include "stack-c.h"
24
25 /** external functions to be called through this interface **/
26
27 #include "grand.h"
28 #include "clcg4.h"
29 #include "others_generators.h"
30 #include "sciprint.h"
31 #include "Scierror.h"
32 #include "gw_randlib.h"
33
34 enum {MT, KISS, CLCG4, CLCG2, URAND, FSULTRA};
35
36 /* the current generator : */
37 static int current_gen = MT;
38
39 /* for clcg4 : the current virtual gen (current_clcg4 in [0, Maxgen]) */
40 static int current_clcg4 = 0;
41
42 /* clcg4 must be called with the virtual generator number */
43 static unsigned long int clcg4_with_gen(void)
44 {
45     return ( clcg4(current_clcg4) );
46 }
47
48 #define NbGenInScilab 6
49
50 /*  pointers onto the generators func */
51 unsigned long int (*gen[NbGenInScilab])(void) = { randmt, kiss,  clcg4_with_gen, clcg2 , urandc , fsultra};
52
53 /*  names at the scilab level */
54 static char *names_gen[NbGenInScilab] = { "mt",  "kiss", "clcg4", "clcg2", "urand", "fsultra" };
55
56 /* all the generators provided integers in [0, RngMaxInt] :        */
57 static
58 unsigned long RngMaxInt[NbGenInScilab] = { 4294967295ul,  /* mt    */
59         4294967295ul,  /* kiss  */
60         2147483646ul,  /* clcg4 */
61         2147483561ul,  /* clcg2 */
62         2147483647ul,  /* urand */
63         4294967295ul
64                                          }; /* fsultra*/
65 /* the factors (1/(RngMaxInt+1)) to get reals in [0,1) :           */
66 static
67 double factor[NbGenInScilab] = { 2.3283064365386963e-10,  /* mt    */
68                                  2.3283064365386963e-10,  /* kiss  */
69                                  4.6566128752457969e-10,  /* clcg4 */
70                                  4.6566130595601735e-10,  /* clcg2 */
71                                  4.6566128730773926e-10,  /* urand */
72                                  2.3283064365386963e-10
73                                }; /* fsultra*/
74
75 double C2F(ranf)(void)
76 {
77     /* random deviate from U[0,1) */
78     return ( (double) gen[current_gen]() * factor[current_gen] );
79 }
80
81 double ignlgi(void)
82 {
83     /* random deviate from Ui[0,RngMaxInt] (direct output of the current gen) */
84     return ( (double) gen[current_gen]() );
85 }
86
87 double C2F(ignuin)(double *a, double *b)
88 {
89     /*  random deviate from Ui[a,b]
90     *  it is assumed that : (i)  a and b are integers (stored in double)
91     *                       (ii) b-a+1 <= RngMaxInt[current_gen]
92     *  (these verif are done at the calling level)
93     *
94     *  We use the classic method with a minor difference : to choose
95     *  uniformly an int in [a,b] (ie d=b-a+1 numbers) with a generator
96     *  which provides uniformly integers in [0,RngMaxInt] (ie m=RngMaxInt+1
97     *  numbers) we do the Euclidian division :
98     *                                           m = q d + r,   r in [0,d-1]
99     *
100     *  and accept only numbers l in [0, qd-1], then the output is k = a + (l mod d)
101     *  (ie numbers falling in [qd , RngMaxInt] are rejected).
102     *  The problem is that RngMaxInt is 2^32-1 for mt and kiss so that RngMaxInt+1 = 0
103     *  with the 32 bits unsigned int arithmetic. So in place of rejected r
104     *  numbers we reject r+1 by using RngMaxInt in place of m. The constraint is
105     *  then that (b-a+1) <= RngMaxInt and if we doesn't want to deal we each generator
106     *  we take (b-a+1) <= Min RngMaxInt =  2147483561 (clcg2)
107     */
108     unsigned long k, d = (unsigned long)((*b - *a) + 1), qd;
109
110     if ( d == 1)
111     {
112         return (*a);
113     }
114
115     qd = RngMaxInt[current_gen] - RngMaxInt[current_gen] % d;
116     do
117     {
118         k = (unsigned long)ignlgi();
119     }
120     while ( k >= qd );
121     return ( *a + (double)(k % d) );
122 }
123
124 /**************************************************
125 *  hand written interface for the randlib
126 ***********************************************************************/
127
128 int sci_Rand(char *fname, unsigned long fname_len)
129 {
130     int minrhs = 1, maxrhs = 10, minlhs = 1, maxlhs = 2;
131     int ResL, ResC, suite, m2, n2, l2, m1, n1, l1, ls, ms, ns, la, lr, lb, lc;
132     int l3, l4;
133     int i;
134
135     Nbvars = 0;
136     CheckRhs(minrhs, maxrhs);
137     CheckLhs(minlhs, maxlhs);
138     if ( GetType(1) != sci_matrix)
139     {
140         int un = 1, deux = 2, dim_state_mt = 625, dim_state_fsultra = 40, dim_state_4 = 4;
141         GetRhsVar(1, STRING_DATATYPE, &ms, &ns, &ls);
142         if ( strcmp(cstk(ls), "getsd") == 0)
143         {
144             if ( Rhs != 1 )
145             {
146                 Scierror(999, _("%s: Wrong number of input argument: %d expected with option '%s'.\n"), fname, 1, "getsd");
147                 return 0;
148             }
149             if ( Lhs != 1 )
150             {
151                 Scierror(999, _("%s: Wrong number of output argument: %d expected the option '%s'.\n"), fname, 1, "getsd");
152                 return 0;
153             }
154
155             switch (current_gen)
156             {
157                 case(MT) :
158                     CreateVar(Rhs + 2, MATRIX_OF_DOUBLE_DATATYPE, &dim_state_mt, &un, &lr);
159                     get_state_mt(stk(lr));
160                     break;
161                 case(KISS) :
162                     CreateVar(Rhs + 2, MATRIX_OF_DOUBLE_DATATYPE, &dim_state_4, &un, &lr);
163                     get_state_kiss(stk(lr));
164                     break;
165                 case(CLCG4) :
166                     CreateVar(Rhs + 2, MATRIX_OF_DOUBLE_DATATYPE, &dim_state_4, &un, &lr);
167                     get_state_clcg4(current_clcg4, stk(lr));
168                     break;
169                 case(CLCG2) :
170                     CreateVar(Rhs + 2, MATRIX_OF_DOUBLE_DATATYPE, &deux, &un, &lr);
171                     get_state_clcg2(stk(lr));
172                     break;
173                 case(URAND) :
174                     CreateVar(Rhs + 2, MATRIX_OF_DOUBLE_DATATYPE, &un, &un, &lr);
175                     get_state_urand(stk(lr));
176                     break;
177                 case(FSULTRA) :
178                     CreateVar(Rhs + 2, MATRIX_OF_DOUBLE_DATATYPE, &dim_state_fsultra, &un, &lr);
179                     get_state_fsultra(stk(lr));
180                     break;
181             };
182             LhsVar(1) = Rhs + 2;
183             PutLhsVar();
184             return 0;
185         }
186         else if ( strcmp(cstk(ls), "setall") == 0 )
187         {
188             if ( current_gen != CLCG4 )
189             {
190                 sciprint(_("The %s option affects only the clcg4 generator\n"), "setall");
191             }
192             if ( Rhs != 5 )
193             {
194                 Scierror(999, _("%s: Wrong number of input arguments: %d expected with option '%s'.\n"), fname, 5, "setall");
195                 return 0;
196             }
197             GetRhsVar(2, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &l1);
198             if ( m1*n1 != 1)
199             {
200                 Scierror(999, _("%s: Wrong type for input argument #%d: Scalar expected.\n"), fname, 2);
201                 return 0;
202             }
203             GetRhsVar(3, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &l2);
204             if ( m1*n1 != 1)
205             {
206                 Scierror(999, _("%s: Wrong type for input argument #%d: Scalar expected.\n"), fname, 3);
207                 return 0;
208             }
209             GetRhsVar(4, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &l3);
210             if ( m1*n1 != 1)
211             {
212                 Scierror(999, _("%s: Wrong type for input argument #%d: Scalar expected.\n"), fname, 4);
213                 return 0;
214             }
215             GetRhsVar(5, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &l4);
216             if ( m1*n1 != 1)
217             {
218                 Scierror(999, _("%s: Wrong type for input argument #%d: Scalar expected.\n"), fname, 5);
219                 return 0;
220             }
221
222             if (! set_initial_seed_clcg4(*stk(l1), *stk(l2), *stk(l3), *stk(l4)) )
223             {
224                 /* => seeds were not good  (info is displayed by the function) */
225                 SciError(999);
226                 return 0;
227             }
228             LhsVar(1) = 1;
229             PutLhsVar();
230             return(0);
231         }
232         else if ( strcmp(cstk(ls), "setsd") == 0 )
233         {
234             switch (current_gen)
235             {
236                 case(MT) :
237                     if ( Rhs != 2 )
238                     {
239                         Scierror(999, _("%s: Wrong number of input arguments: %d expected for '%s' with the %s generator.\n"), fname, 2, "setsd", "mt");
240                         return 0;
241                     }
242                     GetRhsVar(2, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &l1);
243                     if ( m1*n1 == 1)          /* simple init of mt     */
244                     {
245                         if (! set_state_mt_simple(*stk(l1)) )
246                         {
247                             SciError(999);
248                             return(0);
249                         };
250                     }
251                     else if ( m1*n1 == 625 )  /* init of all the state */
252                     {
253                         if (! set_state_mt(stk(l1)))
254                         {
255                             SciError(999);
256                             return(0);
257                         };
258                     }
259                     else
260                     {
261                         Scierror(999, _("%s: Wrong values for input argument: Vector of %d or %d values for %s expected.\n"), fname, 1, 625, "mt");
262                         return 0;
263                     };
264                     break;
265
266                 case(FSULTRA) :
267                     if ( Rhs == 2 ) /* init via a "complete" state */
268                     {
269                         GetRhsVar(2, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &l1);
270                         if ( m1 != 40  ||  n1 != 1)
271                         {
272                             Scierror(999, _("%s: Wrong size for input argument #%d: %dx%d expected.\n"), fname, 2, 40, 1);
273                             return 0;
274                         };
275                         if (! set_state_fsultra(stk(l1)) )
276                         {
277                             SciError(999);
278                             return(0);
279                         }
280                         ;
281                     }
282                     else if ( Rhs == 3 ) /* init with 2 integers (like before) */
283                     {
284                         GetRhsVar(2, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &l1);
285                         if ( m1*n1 != 1)
286                         {
287                             Scierror(999, _("%s: Wrong type for input argument #%d: Scalar expected.\n"), fname, 2);
288                             return 0;
289                         };
290                         GetRhsVar(3, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &l2);
291                         if ( m1*n1 != 1)
292                         {
293                             Scierror(999, _("%s: Wrong type for input argument #%d: Scalar expected.\n"), fname, 3);
294                             return 0;
295                         };
296                         if (! set_state_fsultra_simple(*stk(l1), *stk(l2)) )
297                         {
298                             SciError(999);
299                             return(0);
300                         };
301                     }
302                     else
303                     {
304                         Scierror(999, _("%s: Wrong number of input arguments: %d or %d expected for '%s' option with the %s generator.\n"), fname, 2, 3, "setsd", "fsultra");
305                         return 0;
306                     }
307                     break;
308
309                 case(KISS) :
310                 case(CLCG4) :
311                     if ( Rhs != 5 )
312                     {
313                         Scierror(999, _("%s: Wrong number of input arguments: expected %d for '%s' option with the %s or %s generator.\n"), fname, 5, "setsd", "kiss", "clcg4");
314                         return 0;
315                     }
316                     GetRhsVar(2, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &l1);
317                     if ( m1*n1 != 1)
318                     {
319                         Scierror(999, _("%s: Wrong type for input argument #%d: Scalar expected.\n"), fname, 2);
320                         return 0;
321                     }
322                     GetRhsVar(3, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &l2);
323                     if ( m1*n1 != 1)
324                     {
325                         Scierror(999, _("%s: Wrong type for input argument #%d: Scalar expected.\n"), fname, 3);
326                         return 0;
327                     }
328                     GetRhsVar(4, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &l3);
329                     if ( m1*n1 != 1)
330                     {
331                         Scierror(999, _("%s: Wrong type for input argument #%d: Scalar expected.\n"), fname, 4);
332                         return 0;
333                     }
334                     GetRhsVar(5, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &l4);
335                     if ( m1*n1 != 1)
336                     {
337                         Scierror(999, _("%s: Wrong type for input argument #%d: Scalar expected.\n"), fname, 5);
338                         return 0;
339                     }
340                     if (current_gen == KISS)
341                     {
342                         if (! set_state_kiss(*stk(l1), *stk(l2), *stk(l3), *stk(l4)))
343                         {
344                             SciError(999);
345                             return 0;
346                         };
347                     }
348                     else
349                     {
350                         if (! set_seed_clcg4(current_clcg4, *stk(l1), *stk(l2), *stk(l3), *stk(l4)))
351                         {
352                             SciError(999);
353                             return 0;
354                         };
355                     }
356                     break;
357
358                 case(CLCG2) :
359                     if ( Rhs != 3 )
360                     {
361                         Scierror(999, _("%s: Wrong number of input arguments: %d expected for '%s' option with the %s generator.\n"), fname, 3, "setsd", "clcg2");
362                         return 0;
363                     }
364                     GetRhsVar(2, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &l1);
365                     if ( m1*n1 != 1)
366                     {
367                         Scierror(999, _("%s: Wrong type for input argument #%d: Scalar expected.\n"), fname, 2);
368                         return 0;
369                     };
370                     GetRhsVar(3, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &l2);
371                     if ( m1*n1 != 1)
372                     {
373                         Scierror(999, _("%s: Wrong type for input argument #%d: Scalar expected.\n"), fname, 3);
374                         return 0;
375                     };
376                     if (! set_state_clcg2(*stk(l1), *stk(l2)))
377                     {
378                         SciError(999);
379                         return 0;
380                     };
381                     break;
382
383                 case(URAND) :
384                     if ( Rhs != 2 )
385                     {
386                         Scierror(999, _("%s: Wrong number of input arguments: %d expected for '%s' option with the %s generator.\n"), fname, 2, "setsd", "urand");
387                         return 0;
388                     }
389                     GetRhsVar(2, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &l1);
390                     if ( m1*n1 != 1)
391                     {
392                         Scierror(999, _("%s: Wrong type for input argument #%d: Scalar expected.\n"), fname, 2);
393                         return 0;
394                     };
395                     if (! set_state_urand(*stk(l1)))
396                     {
397                         SciError(999);
398                         return 0;
399                     };
400                     break;
401             };
402             LhsVar(1) = 0;
403             PutLhsVar();
404             return 0;
405         }
406         else if (strcmp("phr2sd", cstk(ls)) == 0)
407         {
408             if ( Rhs != 2 )
409             {
410                 Scierror(999, _("%s: Wrong number of input arguments: %d expected with option '%s'.\n"), fname, 2, "phr2sd");
411                 return 0;
412             }
413             if ( Lhs > 1 )
414             {
415                 Scierror(999, _("%s: Wrong number of output argument: %d expected with option '%s'.\n"), fname, 1, "phr2sd");
416
417                 return 0;
418             }
419             GetRhsVar(2, STRING_DATATYPE, &m1, &n1, &l1);
420             CreateVar(3, MATRIX_OF_INTEGER_DATATYPE, &un, &deux, &l2);
421
422             C2F(phrtsd)(cstk(l1), &m1, istk(l2), istk(l2 + 1), m1);
423             LhsVar(1) = 3;
424             PutLhsVar();
425             return 0;
426         }
427
428         else if (strcmp("initgn", cstk(ls)) == 0)
429         {
430             SeedType Where;
431             if ( current_gen != CLCG4 )
432             {
433                 sciprint(_("%s: The %s option affects only the %s generator\n"), fname, "initgn", "clcg4");
434             }
435             if ( Rhs != 2)
436             {
437                 Scierror(999, _("%s: Wrong number of input arguments: %d expected with option '%s'.\n"), fname, 2, "initgn");
438                 return 0;
439             }
440             GetRhsVar(2, MATRIX_OF_INTEGER_DATATYPE, &m1, &n1, &l1);
441             if ( *istk(l1) != 0 && *istk(l1) != -1 && *istk(l1) != 1)
442             {
443                 Scierror(999, _("%s: Wrong value for input argument #%d: %d, %d or %d expected.\n"), fname, 2, -1, 0, 1);
444                 return 0;
445             }
446             Where = (SeedType) (*istk(l1) + 1);
447             init_generator_clcg4(current_clcg4, Where);
448             LhsVar(1) = 2;
449             PutLhsVar();
450             return 0;
451         }
452         else if (strcmp("setcgn", cstk(ls)) == 0)
453         {
454             if ( current_gen != CLCG4 )
455             {
456                 sciprint(_("The %s option affects only the %s generator\n"), "setcgn", "clcg4");
457             }
458             if ( Rhs != 2)
459             {
460                 Scierror(999, _("%s: Wrong number of input arguments: %d expected with option '%s'.\n"), fname, 2, "setcgn");
461                 return 0;
462             }
463             GetRhsVar(2, MATRIX_OF_INTEGER_DATATYPE, &m1, &n1, &l1);
464             if ( *istk(l1) < 0 || *istk(l1) > Maxgen )
465             {
466                 Scierror(999, _("%s: Wrong value for input argument #%d: Must be between %d and %d.\n"), fname, 2, 0, Maxgen);
467                 return 0;
468             }
469             current_clcg4 = *istk(l1);
470             LhsVar(1) = 2;
471             PutLhsVar();
472             return 0;
473         }
474         else if (strcmp("advnst", cstk(ls)) == 0)
475         {
476             int k;
477             if ( current_gen != CLCG4 )
478             {
479                 sciprint(_("The %s option affects only the %s generator\n"), "advnst", "clcg4");
480             }
481             if ( Rhs != 2)
482             {
483                 Scierror(999, _("%s: Wrong number of input arguments: %d expected with option '%s'.\n"), fname, 2, "advnst");
484                 return 0;
485             }
486             GetRhsVar(2, MATRIX_OF_INTEGER_DATATYPE, &m1, &n1, &l1);
487             k = *istk(l1);
488             if ( k < 1 )
489             {
490                 Scierror(999, _("%s: Wrong value for input argument #%d: Must be > %d.\n"), fname, 2, 0);
491                 return 0;
492             }
493             advance_state_clcg4(current_clcg4, k);
494             LhsVar(1) = 2;
495             PutLhsVar();
496             return 0;
497         }
498         else if (strcmp("getcgn", cstk(ls)) == 0)
499         {
500             if ( Rhs != 1)
501             {
502                 Scierror(999, _("%s: Wrong number of input argument: %d expected with option '%s'.\n"), fname, 1, "getcgn");
503                 return 0;
504             }
505             if ( current_gen != CLCG4 )
506             {
507                 sciprint(_("This information concerns only the clcg4 generator\n"));
508             }
509             CreateVar(2, MATRIX_OF_INTEGER_DATATYPE, &un, &un, &l1);
510             *istk(l1) = current_clcg4;
511             LhsVar(1) = 2;
512             PutLhsVar();
513             return 0;
514         }
515         else if (strcmp("setgen", cstk(ls)) == 0)
516         {
517             int msb, nsb, lsb;
518             if ( Rhs != 2)
519             {
520                 Scierror(999, _("%s: Wrong number of input arguments: %d expected with option '%s'.\n"), fname, 2, "setgen");
521                 return 0;
522             }
523             GetRhsVar(2, STRING_DATATYPE, &msb, &nsb, &lsb);
524             if (strcmp("mt", cstk(lsb)) == 0)
525             {
526                 current_gen = MT;
527             }
528             else if (strcmp("kiss", cstk(lsb)) == 0)
529             {
530                 current_gen = KISS;
531             }
532             else if (strcmp("clcg4", cstk(lsb)) == 0)
533             {
534                 current_gen = CLCG4;
535             }
536             else if (strcmp("clcg2", cstk(lsb)) == 0)
537             {
538                 current_gen = CLCG2;
539             }
540             else if (strcmp("urand", cstk(lsb)) == 0)
541             {
542                 current_gen = URAND;
543             }
544             else if (strcmp("fsultra", cstk(lsb)) == 0)
545             {
546                 current_gen = FSULTRA;
547             }
548             else
549             {
550                 Scierror(999, _("%s: Wrong value for input argument #%d: '%s', '%s', '%s', '%s', '%s' or '%s' expected.\n"), fname, 2, "mt", "kiss", "clcg4", "clcg2", "urand", "fsultra");
551                 return 0;
552             }
553             LhsVar(1) = 2;
554             PutLhsVar();
555             return 0;
556         }
557         else if (strcmp("getgen", cstk(ls)) == 0)
558         {
559             int l_un = 1;
560             if ( Rhs != 1)
561             {
562                 Scierror(999, _("%s: Wrong number of input argument: %d expected with option '%s'.\n"), fname, 1, "getgen");
563                 return 0;
564             }
565             CreateVarFromPtr( Rhs + 2, MATRIX_OF_STRING_DATATYPE, &l_un, &l_un, &names_gen[current_gen]);
566             LhsVar(1) = Rhs + 2;
567             PutLhsVar();
568             return 0;
569         }
570         else
571         {
572             Scierror(999, _("%s Wrong value for input argument #%d: %s.\n"), fname, 1, cstk(ls));
573
574             return 0;
575         }
576     }
577     minrhs = 2;
578     CheckRhs(minrhs, maxrhs);
579     if ( GetType(2) == sci_matrix ) /** m,n,'string' */
580     {
581         GetRhsVar(1, MATRIX_OF_INTEGER_DATATYPE, &m1, &n1, &l1);
582         if ( m1*n1 != 1)
583         {
584             Scierror(999, _("%s: Wrong type for input argument #%d: Scalar expected.\n"), fname, 1);
585             return 0;
586         }
587         ResL = *istk(l1);
588         GetRhsVar(2, MATRIX_OF_INTEGER_DATATYPE, &m2, &n2, &l2);
589         if ( m2*n2 != 1)
590         {
591             Scierror(999, _("%s: Wrong type for input argument #%d: Scalar expected.\n"), fname, 2);
592             return 0;
593         }
594         ResC = *istk(l2);
595         GetRhsVar(3, STRING_DATATYPE, &ms, &ns, &ls);
596         suite = 4;
597         if (ResL < 0 && (ResL != -1 || ResC != -1)) //ResL=-1 & ResC=-1 => eye
598         {
599             Scierror(999, _("%s: Wrong value for input argument #%d: Positive scalar expected.\n"), fname, 1);
600             return 0;
601         }
602
603         if (ResC < 0 && (ResL != -1 || ResC != -1)) //ResL=-1 & ResC=-1 => eye
604         {
605             Scierror(999, _("%s: Wrong value for input argument #%d: Positive scalar expected.\n"), fname, 2);
606             return 0;
607         }
608     }
609     else
610     {
611         GetRhsVar(1, MATRIX_OF_INTEGER_DATATYPE, &ResL, &ResC, &l1);
612         GetRhsVar(2, STRING_DATATYPE, &ms, &ns, &ls);
613         suite = 3;
614     }
615     if ( strcmp(cstk(ls), "bet") == 0)
616     {
617         double minlog = 1.e-37;
618         if ( Rhs != suite + 1)
619         {
620             Scierror(999, _("Missing A and B for beta law\n"));
621             return 0;
622         }
623         GetRhsVar(suite, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &la);
624         if ( m1*n1 != 1)
625         {
626             Scierror(999, _("%s: Wrong size for input argument: Scalar expected for %s.\n"), fname, "A");
627             return 0;
628         }
629         GetRhsVar(suite + 1, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &lb);
630         if ( m1*n1 != 1)
631         {
632             Scierror(999, _("%s: Wrong size for input argument: Scalar expected for %s.\n"), fname, "B");
633             return 0;
634         }
635         CreateVar(suite + 2, MATRIX_OF_DOUBLE_DATATYPE, &ResL, &ResC, &lr);
636         if ( *stk(la) < minlog || *stk(lb) < minlog)
637         {
638             Scierror(999, _("Rand(...,'bet',..): A or B < %f\n"), minlog);
639             return 0;
640         }
641         for ( i = 0 ; i < ResL * ResC ; i++)
642         {
643             *stk(lr + i) = C2F(genbet)(stk(la), stk(lb));
644         }
645         LhsVar(1) = suite + 2;
646         PutLhsVar();
647         return 0;
648     }
649     else if ( strcmp(cstk(ls), "f") == 0)
650     {
651         if ( Rhs != suite + 1)
652         {
653             Scierror(999, _("Missing Dfn and Dfd for F law\n"));
654             return 0;
655         }
656         GetRhsVar(suite, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &la);
657         if ( m1*n1 != 1)
658         {
659             Scierror(999, _("%s: Wrong size for input argument: Scalar expected for %s.\n"), fname, "Dfn");
660             return 0;
661         }
662         GetRhsVar(suite + 1, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &lb);
663         if ( m1*n1 != 1)
664         {
665             Scierror(999, _("%s: Wrong size for input argument: Scalar expected for %s.\n"), fname, "Dfd");
666             return 0;
667         }
668         CreateVar(suite + 2, MATRIX_OF_DOUBLE_DATATYPE, &ResL, &ResC, &lr);
669         if ( *stk(la) <= 0.0 || *stk(lb) <= 0.0)
670         {
671             Scierror(999, _("Degrees of freedom nonpositive\n"));
672             return 0;
673         }
674         for ( i = 0 ; i < ResL * ResC ; i++)
675         {
676             *stk(lr + i) = C2F(genf)(stk(la), stk(lb));
677         }
678         LhsVar(1) = suite + 2;
679         PutLhsVar();
680         return 0;
681     }
682     else if ( strcmp(cstk(ls), "mul") == 0)
683     {
684         int l_i, nn, ncat;
685         double ptot;
686         if ( suite != 3 || ResL*ResC != 1)
687         {
688             Scierror(999, _("%s: Wrong value for input argument #%d: Must be the number of random deviate.\n"), fname, 1);
689             return 0;
690         }
691         nn = *istk(l1);
692         if ( Rhs != suite + 1)
693         {
694             Scierror(999, _("Missing N and P for MULtinomial law\n"));
695             return 0;
696         }
697         GetRhsVar(suite, MATRIX_OF_INTEGER_DATATYPE, &m1, &n1, &la);
698         if ( m1*n1 != 1)
699         {
700             Scierror(999, _("%s: Wrong size for input argument: Scalar expected for %s.\n"), fname, "N");
701             return 0;
702         }
703         GetRhsVar(suite + 1, MATRIX_OF_DOUBLE_DATATYPE, &m2, &n2, &lb);
704         if ( n2 != 1 )
705         {
706             Scierror(999, _("%s: Wrong size for input argument: Column vector expected.\n"), fname);
707             return 0;
708         }
709         ncat = m2 + 1;
710         CreateVar(suite + 2, MATRIX_OF_INTEGER_DATATYPE, &ncat, &nn, &lr);
711         if ( *istk(la) < 0 )
712         {
713             Scierror(999, _("N < 0\n"));
714             return 0;
715         }
716         if ( ncat <= 1)
717         {
718             Scierror(999, _("Ncat <= 1\n"));
719             return 0;
720         }
721         ptot = 0.0;
722         for ( l_i = 0 ; l_i < ncat - 1 ; l_i++ )
723         {
724             if ( *stk(lb + l_i) < 0.0 )
725             {
726                 Scierror(999, _("P(%d) < 0\n"), l_i + 1);
727                 return 0;
728             }
729             if ( *stk(lb + l_i) > 1.0 )
730             {
731                 Scierror(999, _("P(%d) > 1\n"), l_i + 1);
732                 return 0;
733             }
734             ptot += *stk(lb + l_i);
735         }
736         if ( ptot > 1.0)
737         {
738             Scierror(999, _("Sum of P(i) > 1\n"));
739             return 0;
740         }
741         for ( l_i = 0 ; l_i < nn ; l_i++)
742         {
743             C2F(genmul)(istk(la), stk(lb), &ncat, istk(lr + ncat * l_i));
744         }
745         LhsVar(1) = suite + 2;
746         PutLhsVar();
747         return 0;
748     }
749     else if ( strcmp(cstk(ls), "gam") == 0)
750     {
751         if ( Rhs != suite + 1)
752
753             /*  ETRE PLUS CONSISTANT ICI : choisir entre shape , rate ou
754             bien A et R (idem pour le man)
755             */
756         {
757             Scierror(999, _("Missing shape and rate for Gamma law\n"));
758             return 0;
759         }
760         GetRhsVar(suite, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &la);
761         if ( m1*n1 != 1)
762         {
763             Scierror(999, _("%s: Wrong size for input argument: Scalar expected for %s.\n"), fname, "shape");
764             return 0;
765         }
766         GetRhsVar(suite + 1, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &lb);
767         if ( m1*n1 != 1)
768         {
769             Scierror(999, _("%s: Wrong size for input argument: Scalar expected for %s.\n"), fname, "rate");
770             return 0;
771         }
772         CreateVar(suite + 2, MATRIX_OF_DOUBLE_DATATYPE, &ResL, &ResC, &lr);
773         if ( (*stk(la)) <= 0.0 ||  (*stk(lb)) <= 0.0 )
774         {
775             Scierror(999, _("grand(..'gam',A,R) : A <= 0.0 or R <= 0.0\n"));
776             return 0;
777         }
778         for ( i = 0 ; i < ResL * ResC ; i++)
779         {
780             /** WARNING : order is changed in parameters for
781             compatibility between Rand(...'gam',..) and cdfgam
782             **/
783             *stk(lr + i) = C2F(gengam)(stk(lb), stk(la));
784         }
785         LhsVar(1) = suite + 2;
786         PutLhsVar();
787         return 0;
788     }
789
790     else if ( strcmp(cstk(ls), "nor") == 0)
791     {
792         if ( Rhs != suite + 1)
793         {
794             Scierror(999, _("Missing Av and Sd for Normal law\n"));
795             return 0;
796         }
797         GetRhsVar(suite, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &la);
798         if ( m1*n1 != 1)
799         {
800             Scierror(999, _("%s: Wrong size for input argument: Scalar expected for %s.\n"), fname, "Av");
801             return 0;
802         }
803         GetRhsVar(suite + 1, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &lb);
804         if ( m1*n1 != 1)
805         {
806             Scierror(999, _("%s: Wrong size for input argument: Scalar expected for %s.\n"), fname, "Sd");
807             return 0;
808         }
809         CreateVar(suite + 2, MATRIX_OF_DOUBLE_DATATYPE, &ResL, &ResC, &lr);
810         if ( *stk(lb) < 0 )
811         {
812             Scierror(999, _("SD < 0.0\n"));
813             return 0;
814         }
815         for ( i = 0 ; i < ResL * ResC ; i++)
816         {
817             *stk(lr + i) = C2F(gennor)(stk(la), stk(lb));
818         }
819         LhsVar(1) = suite + 2;
820         PutLhsVar();
821         return 0;
822     }
823     else if ( strcmp(cstk(ls), "unf") == 0)
824     {
825         double low = 0, high = 0;
826         if ( Rhs != suite + 1)
827         {
828             Scierror(999, _("Missing Low and High for Uniform Real law\n"));
829             return 0;
830         }
831
832         GetRhsVar(suite, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &la);
833         if ( m1*n1 != 1)
834         {
835             Scierror(999, _("%s: Wrong type for input argument: Scalar expected.\n"), fname);
836             return 0;
837         }
838
839         GetRhsVar(suite + 1, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &lb);
840         if ( m1*n1 != 1)
841         {
842             Scierror(999, _("%s: Wrong type for input argument: Scalar expected.\n"), fname);
843             return 0;
844         }
845
846         CreateVar(suite + 2, MATRIX_OF_DOUBLE_DATATYPE, &ResL, &ResC, &lr);
847         low = *stk(la);
848         high =  *stk(lb);
849         if ( low > high )
850         {
851             Scierror(999, _("%s: Wrong type for input argument. Low < High expected.\n"), fname);
852             return 0;
853         }
854         for ( i = 0 ; i < ResL * ResC ; i++)
855         {
856             *stk(lr + i) = low + (high - low) * C2F(ranf)();
857         }
858         LhsVar(1) = suite + 2;
859         PutLhsVar();
860         return 0;
861     }
862     else if ( strcmp(cstk(ls), "uin") == 0)
863     {
864         double a = 0, b = 0;
865         if ( Rhs != suite + 1)
866         {
867             Scierror(999, _("Missing Low and High for Uniform int law\n"));
868             return 0;
869         }
870
871         GetRhsVar(suite, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &la);
872
873         if ( m1*n1 != 1)
874         {
875             Scierror(999, _("%s: Wrong type for input argument: Scalar expected.\n"), fname);
876             return 0;
877         }
878
879         GetRhsVar(suite + 1, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &lb);
880         if ( m1*n1 != 1)
881         {
882             Scierror(999, _("%s: Wrong type for input argument: Scalar expected.\n"), fname);
883             return 0;
884         }
885         a = *stk(la);
886         b = *stk(lb);
887
888         if ( a > b )
889         {
890             Scierror(999, _("%s: Wrong type for input argument. Low < High expected.\n"), fname);
891             return 0;
892         }
893
894         if ( a != floor(a) || b != floor(b) || (b - a + 1) > 2147483561 )
895         {
896             Scierror(999, _("a and b must integers with (b-a+1) <= 2147483561"));
897             return 0;
898         }
899
900         CreateVar(suite + 2, MATRIX_OF_DOUBLE_DATATYPE, &ResL, &ResC, &lr);
901         for ( i = 0 ; i < ResL * ResC ; i++)
902         {
903             *stk(lr + i) = C2F(ignuin)(stk(la), stk(lb));
904         }
905         LhsVar(1) = suite + 2;
906         PutLhsVar();
907         return 0;
908     }
909     else if ( strcmp(cstk(ls), "lgi") == 0)
910     {
911         if ( Rhs != suite - 1 )
912         {
913             Scierror(999, _("%s: Wrong number of input argument: %d expected with option '%s'.\n"), fname, suite - 1, "lgi");
914             return 0;
915         }
916         CreateVar(suite, MATRIX_OF_DOUBLE_DATATYPE, &ResL, &ResC, &lr);
917         for ( i = 0 ; i < ResL * ResC ; i++)
918         {
919             *stk(lr + i) = ignlgi();
920         }
921         LhsVar(1) = suite;
922         PutLhsVar();
923         return 0;
924     }
925     else if ( strcmp(cstk(ls), "prm") == 0)
926     {
927         int nn;
928         if ( suite != 3 || ResL*ResC != 1)
929         {
930             Scierror(999, _("%s: Wrong value for input argument: Number of random simulation expected.\n"), fname);
931             return 0;
932         }
933         nn = *istk(l1);
934         if ( Rhs != suite)
935         {
936             Scierror(999, _("Missing vect for random permutation\n"));
937             return 0;
938         }
939         GetRhsVar(suite, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &la);
940         if ( n1 != 1)
941         {
942             Scierror(999, _("%s: Wrong type for input argument: Column vector expected.\n"), fname);
943             return 0;
944         }
945         CreateVar(suite + 1, MATRIX_OF_DOUBLE_DATATYPE, &m1, &nn, &lr);
946         for ( i = 0 ; i < nn ; i++)
947         {
948             int j ;
949             for (j = 0; j < m1 ; j++ )
950             {
951                 *stk(lr + (m1)*i + j) = *stk(la + j);
952             }
953             C2F(genprm)(stk(lr + (m1)*i), &m1);
954         }
955         LhsVar(1) = suite + 1;
956         PutLhsVar();
957         return 0;
958     }
959     else if ( strcmp(cstk(ls), "nbn") == 0)
960     {
961         if ( Rhs != suite + 1)
962         {
963             Scierror(999, _("Missing N and P for Negative Binomial law\n"));
964             return 0;
965         }
966         GetRhsVar(suite, MATRIX_OF_INTEGER_DATATYPE, &m1, &n1, &la);
967         if ( m1*n1 != 1)
968         {
969             Scierror(999, _("%s: Wrong type for input argument: Scalar expected.\n"), fname);
970             return 0;
971         }
972         GetRhsVar(suite + 1, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &lb);
973         if ( m1*n1 != 1)
974         {
975             Scierror(999, _("%s: Wrong type for input argument: Scalar expected.\n"), fname);
976             return 0;
977         }
978         CreateVar(suite + 2, MATRIX_OF_DOUBLE_DATATYPE, &ResL, &ResC, &lr);
979         if ( *stk(lb) < 0.0 || *stk(lb) > 1.0 )
980         {
981             Scierror(999, _("P is not in [0,1]\n"));
982             return 0;
983         }
984         if ( *istk(la) < 0 )
985         {
986             Scierror(999, _("N < 0\n"));
987             return 0;
988         }
989         for ( i = 0 ; i < ResL * ResC ; i++)
990         {
991             *stk(lr + i) = (double) C2F(ignnbn)(istk(la), stk(lb));
992         }
993         LhsVar(1) = suite + 2;
994         PutLhsVar();
995         return 0;
996     }
997     else if ( strcmp(cstk(ls), "bin") == 0)
998     {
999         if ( Rhs != suite + 1)
1000         {
1001             Scierror(999, _("Missing N and P for Binomial law\n"));
1002             return 0;
1003         }
1004         GetRhsVar(suite, MATRIX_OF_INTEGER_DATATYPE, &m1, &n1, &la);
1005         if ( m1*n1 != 1)
1006         {
1007             Scierror(999, _("%s: Wrong type for input argument: Scalar expected.\n"), fname);
1008             return 0;
1009         }
1010         GetRhsVar(suite + 1, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &lb);
1011         if ( m1*n1 != 1)
1012         {
1013             Scierror(999, _("%s: Wrong type for input argument: Scalar expected.\n"), fname);
1014             return 0;
1015         }
1016         CreateVar(suite + 2, MATRIX_OF_DOUBLE_DATATYPE, &ResL, &ResC, &lr);
1017         if ( *stk(lb) < 0.0 || *stk(lb) > 1.0 )
1018         {
1019             Scierror(999, _("P is not in [0,1]\n"));
1020             return 0;
1021         }
1022         if ( *istk(la) < 0 )
1023         {
1024             Scierror(999, _("N < 0\n"));
1025             return 0;
1026         }
1027         for ( i = 0 ; i < ResL * ResC ; i++)
1028         {
1029             *stk(lr + i) = (double) C2F(ignbin)(istk(la), stk(lb));
1030         }
1031         LhsVar(1) = suite + 2;
1032         PutLhsVar();
1033         return 0;
1034     }
1035
1036     else if ( strcmp(cstk(ls), "mn") == 0)
1037     {
1038         int nn, un = 1, work, mp, parm, ierr;
1039         if ( suite != 3 || ResL*ResC != 1)
1040         {
1041             Scierror(999, _("%s: Wrong value for input argument #%d: Must be the number of random simulation.\n"), fname, 1);
1042             return 0;
1043         }
1044         nn = *istk(l1);
1045         if ( Rhs != suite + 1)
1046         {
1047             Scierror(999, _("Missing Mean and Cov for Multivariate Normal law\n"));
1048             return 0;
1049         }
1050         GetRhsVar(suite, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &la);
1051         if ( n1 != 1)
1052         {
1053             Scierror(999, _("%s: Wrong type for input argument: Column vector expected.\n"), fname);
1054             return 0;
1055         }
1056         GetRhsVar(suite + 1, MATRIX_OF_DOUBLE_DATATYPE, &m2, &n2, &lb);
1057         if ( m2 != n2 )
1058         {
1059             Scierror(999, _("%s: Wrong type for input argument: Square matrix expected.\n"), fname);
1060             return 0;
1061         }
1062         if ( m2 != m1 )
1063         {
1064             Scierror(999, _("%s: Wrong type for input arguments: Mean and Cov have incompatible dimensions\n"), fname);
1065             return 0;
1066         }
1067
1068         CreateVar(suite + 2, MATRIX_OF_DOUBLE_DATATYPE, &m1, &nn, &lr);
1069         CreateVar(suite + 3, MATRIX_OF_DOUBLE_DATATYPE, &m1, &un, &work);
1070         mp = m1 * (m1 + 3) / 2 + 1;
1071         CreateVar(suite + 4, MATRIX_OF_DOUBLE_DATATYPE, &mp, &un, &parm);
1072         if ( m1 <= 0 )
1073         {
1074             Scierror(999, _("%s: Wrong size for input arguments: Mean and Cov are of null size.\n"), fname);
1075             return 0;
1076         }
1077         C2F(setgmn)(stk(la), stk(lb), &m2, &m1, stk(parm), &ierr);
1078         if ( ierr == 1)
1079         {
1080             SciError(999);
1081             return 0;
1082         }
1083         for ( i = 0 ; i < nn ; i++)
1084         {
1085             C2F(genmn)(stk(parm), stk(lr + (m1)*i), stk(work));
1086         }
1087         LhsVar(1) = suite + 2;
1088         PutLhsVar();
1089         return 0;
1090     }
1091     else if ( strcmp(cstk(ls), "markov") == 0)
1092     {
1093         int nn, n1p1, lr1, j, icur, mm, jj;
1094         if ( suite != 3 || ResL*ResC != 1)
1095         {
1096             Scierror(999, _("%s: Wrong value for input argument #%d: Must be the number of random simulation.\n"), fname, 1);
1097             return 0;
1098         }
1099         nn = *istk(l1);
1100         if ( Rhs != suite + 1 )
1101         {
1102             Scierror(999, _("%s: Missing P matrix and X0 for Markov chain\n"), "fname");
1103             return 0;
1104         }
1105         GetRhsVar(suite, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &la);
1106         GetRhsVar(suite + 1, MATRIX_OF_INTEGER_DATATYPE, &m2, &n2, &lb);
1107         if ( m1 != n1 && m1 != 1 )
1108         {
1109             Scierror(999, _("%s: Wrong type for input argument #%d: Square matrix or row vector expected.\n"), fname, 2);
1110             return 0;
1111         }
1112
1113         if ( m2*n2 == 0 )
1114         {
1115             Scierror(999, _("X0 is empty\n"));
1116             return 0;
1117         }
1118
1119         for ( i = 0 ; i < m2 * n2 ; i++)
1120             if ( *istk(lb + i) - 1 < 0 || *istk(lb + i) - 1 >= n1 )
1121             {
1122                 Scierror(999, _("%s: X0(%d) must be in the range [1,%d[\n"), fname, i + 1, n1 + 1);
1123                 return 0;
1124             }
1125         mm = m2 * n2;
1126         CreateVar(suite + 2, MATRIX_OF_INTEGER_DATATYPE, &mm, &nn, &lr);
1127
1128         n1p1 = n1 + 1;
1129         CreateVar(suite + 3, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1p1, &lr1);
1130         for ( i = 0 ; i < m1 ; i++ )
1131         {
1132             double ptot = 0.0;
1133             for ( j = 0 ; j < n1 ; j++ )
1134             {
1135                 if ( *stk(la + i + m1 * j) < 0 )
1136                 {
1137                     Scierror(999, _("P(%d,%d) < 0\n"), i + 1, j + 1);
1138                     return 0;
1139                 }
1140                 if ( *stk(la + i + m1 * j) > 1 )
1141                 {
1142                     Scierror(999, _("P(%d,%d) > 1\n"), i + 1, j + 1);
1143                     return 0;
1144                 }
1145                 ptot += *stk(la + i + m1 * j) ;
1146             }
1147             if ( fabs(ptot - 1.0) > 1e-8 )
1148             {
1149                 Scierror(999, _("Sum of P(%d,1:%d)=%f ~= 1\n"), i + 1, n1, ptot);
1150                 return 0;
1151             }
1152         }
1153         /** Computing the cumulative sum of the P matrix **/
1154         for ( i = 0 ; i < m1 ; i++)
1155         {
1156             double cumsum = 0.0;
1157             *stk(lr1 + i) = cumsum;
1158             for ( j = 1; j < n1p1 ; j++ )
1159             {
1160                 cumsum += *stk(la + i + m1 * (j - 1));
1161                 *stk(lr1 + i + m1 * j) = cumsum;
1162             }
1163         }
1164         for ( jj = 0 ; jj < mm ; jj++)
1165         {
1166             icur = *istk(lb + jj) - 1;
1167             for ( i = 0 ; i < nn ; i++)
1168             {
1169                 int niv = 0;
1170                 double rr = C2F(ranf)();
1171                 if ( m1 == 1 )
1172                 {
1173                     icur = 0;
1174                 }
1175                 while ( rr >= *stk(lr1 + icur + m1 * niv) && niv < n1p1 )
1176                 {
1177                     niv++;
1178                 }
1179                 /** projection to avoid boundaries **/
1180                 niv = Max(Min(niv, n1), 1);
1181                 *istk(lr + jj + mm * i) = niv ;
1182                 icur = niv - 1;
1183             }
1184         }
1185         LhsVar(1) = suite + 2;
1186         PutLhsVar();
1187         return 0;
1188     }
1189     else if ( strcmp(cstk(ls), "def") == 0)
1190     {
1191         if ( Rhs != suite - 1 )
1192         {
1193             Scierror(999, _("%s: Wrong number of input argument.\n"), fname);
1194             return 0;
1195         }
1196         CreateVar(suite, MATRIX_OF_DOUBLE_DATATYPE, &ResL, &ResC, &lr);
1197         for ( i = 0 ; i < ResL * ResC ; i++)
1198         {
1199             *stk(lr + i) = C2F(ranf)();
1200         }
1201         LhsVar(1) = suite;
1202         PutLhsVar();
1203         return 0;
1204     }
1205
1206     else if ( strcmp(cstk(ls), "nch") == 0)
1207     {
1208         if ( Rhs != suite + 1)
1209         {
1210             Scierror(999, _("Missing Df and Xnonc for non-central chi-square law\n"));
1211             return 0;
1212         }
1213         GetRhsVar(suite, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &la);
1214         if ( m1*n1 != 1)
1215         {
1216             Scierror(999, _("%s: Wrong type for input argument: Scalar expected.\n"), fname);
1217             return 0;
1218         }
1219         GetRhsVar(suite + 1, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &lb);
1220         if ( m1*n1 != 1)
1221         {
1222             Scierror(999, _("%s: Wrong type for input argument: Scalar expected.\n"), fname);
1223             return 0;
1224         }
1225         CreateVar(suite + 2, MATRIX_OF_DOUBLE_DATATYPE, &ResL, &ResC, &lr);
1226         if ( *stk(la) < 1.0 || *stk(lb) < 0.0 )
1227         {
1228             Scierror(999, _("DF < 1 or XNONC < 0\n"));
1229             return 0;
1230         }
1231         for ( i = 0 ; i < ResL * ResC ; i++)
1232         {
1233             *stk(lr + i) = C2F(gennch)(stk(la), stk(lb));
1234         }
1235         LhsVar(1) = suite + 2;
1236         PutLhsVar();
1237         return 0;
1238     }
1239     else if ( strcmp(cstk(ls), "nf") == 0)
1240     {
1241         if ( Rhs != suite + 2)
1242         {
1243             Scierror(999, _("Missing Dfn, Dfd and Xnonc for non-central F law\n"));
1244             return 0;
1245         }
1246         GetRhsVar(suite, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &la);
1247         if ( m1*n1 != 1)
1248         {
1249             Scierror(999, _("%s: Wrong type for input argument: Scalar expected.\n"), fname);
1250             return 0;
1251         }
1252         GetRhsVar(suite + 1, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &lb);
1253         if ( m1*n1 != 1)
1254         {
1255             Scierror(999, _("%s: Wrong type for input argument: Scalar expected.\n"), fname);
1256             return 0;
1257         }
1258         GetRhsVar(suite + 2, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &lc);
1259         if ( m1*n1 != 1)
1260         {
1261             Scierror(999, _("%s: Wrong type for input argument: Scalar expected.\n"), fname);
1262             return 0;
1263         }
1264         CreateVar(suite + 3, MATRIX_OF_DOUBLE_DATATYPE, &ResL, &ResC, &lr);
1265         if ( *stk(la) < 1.0 || *stk(lb) < 0.0 || *stk(lc) < 0.0 )
1266         {
1267             Scierror(999, _("DF < 1.0 or DF <= 0.0 or Xnonc < 0.0\n"));
1268             return 0;
1269         }
1270         for ( i = 0 ; i < ResL * ResC ; i++)
1271         {
1272             *stk(lr + i) = C2F(gennf)(stk(la), stk(lb), stk(lc));
1273         }
1274         LhsVar(1) = suite + 3;
1275         PutLhsVar();
1276         return 0;
1277     }
1278
1279     else if ( strcmp(cstk(ls), "chi") == 0)
1280     {
1281         if ( Rhs != suite )
1282         {
1283             Scierror(999, _("Missing Df for chi-square law\n"));
1284             return 0;
1285         }
1286         GetRhsVar(suite, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &la);
1287         if ( m1*n1 != 1)
1288         {
1289             Scierror(999, _("%s: Wrong type for input argument: Scalar expected.\n"), fname);
1290             return 0;
1291         }
1292         CreateVar(suite + 1, MATRIX_OF_DOUBLE_DATATYPE, &ResL, &ResC, &lr);
1293         if  ( *stk(la) <= 0.0)
1294         {
1295             Scierror(999, _("Rand: DF <= 0\n"));
1296             return 0;
1297         }
1298         for ( i = 0 ; i < ResL * ResC ; i++)
1299         {
1300             *stk(lr + i) = C2F(genchi)(stk(la));
1301         }
1302         LhsVar(1) = suite + 1;
1303         PutLhsVar();
1304         return 0;
1305     }
1306     else if ( strcmp(cstk(ls), "poi") == 0)
1307     {
1308         if ( Rhs != suite )
1309         {
1310             Scierror(999, _("Missing Av for Poisson law\n"));
1311             return 0;
1312         }
1313         GetRhsVar(suite, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &la);
1314         if ( m1*n1 != 1)
1315         {
1316             Scierror(999, _("%s: Wrong type for input argument: Scalar expected.\n"), fname);
1317             return 0;
1318         }
1319         CreateVar(suite + 1, MATRIX_OF_DOUBLE_DATATYPE, &ResL, &ResC, &lr);
1320         if ( *stk(la) < 0.0 )
1321         {
1322             Scierror(999, _("Av < 0\n"));
1323             return 0;
1324         }
1325         for ( i = 0 ; i < ResL * ResC ; i++)
1326         {
1327             *stk(lr + i) = (double) C2F(ignpoi)(stk(la));
1328         }
1329         LhsVar(1) = suite + 1;
1330         PutLhsVar();
1331         return 0;
1332     }
1333     else if ( strcmp(cstk(ls), "geom") == 0)
1334     {
1335         double p;
1336         if ( Rhs != suite )
1337         {
1338             Scierror(999, _("Missing p for Geometric law\n"));
1339             return 0;
1340         }
1341         GetRhsVar(suite, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &la);
1342         if ( m1*n1 != 1)
1343         {
1344             Scierror(999, _("%s: Wrong type for input argument: Scalar expected.\n"), fname);
1345             return 0;
1346         }
1347         p = *stk(la);
1348         if ( p < 1.3e-307 || p > 1 )
1349         {
1350             Scierror(999, _("%s: Wrong value for input argument: Must be between '%s' and %d.\n"), fname, "pmin", 1);
1351             return 0;
1352         }
1353
1354         CreateVar(suite + 1, MATRIX_OF_DOUBLE_DATATYPE, &ResL, &ResC, &lr);
1355         for ( i = 0 ; i < ResL * ResC ; i++)
1356         {
1357             *stk(lr + i) = igngeom(p);
1358         }
1359         LhsVar(1) = suite + 1;
1360         PutLhsVar();
1361         return 0;
1362     }
1363
1364     else if ( strcmp(cstk(ls), "exp") == 0)
1365     {
1366         if ( Rhs != suite )
1367         {
1368             Scierror(999, _("Missing Av for exponential law\n"));
1369             return 0;
1370         }
1371         GetRhsVar(suite, MATRIX_OF_DOUBLE_DATATYPE, &m1, &n1, &la);
1372         if ( m1*n1 != 1)
1373         {
1374             Scierror(999, _("%s: Wrong type for input argument: Scalar expected.\n"), fname);
1375             return 0;
1376         }
1377         CreateVar(suite + 1, MATRIX_OF_DOUBLE_DATATYPE, &ResL, &ResC, &lr);
1378         if ( *stk(la) < 0.0 )
1379         {
1380             Scierror(999, _("Av < 0.0\n"));
1381             return 0;
1382         }
1383         for ( i = 0 ; i < ResL * ResC ; i++)
1384         {
1385             *stk(lr + i) = C2F(genexp)(stk(la));
1386         }
1387         LhsVar(1) = suite + 1;
1388         PutLhsVar();
1389         return 0;
1390     }
1391
1392     else
1393     {
1394         Scierror(999, _("%s: Wrong value for input argument %s.\n"), fname, cstk(ls));
1395         return 0;
1396     }
1397 }