home *** CD-ROM | disk | FTP | other *** search
/ Geek Gadgets 1 / ADE-1.bin / ade-dist / octave-1.1.1p1-src.tgz / tar.out / fsf / octave / src / f-rand.cc < prev    next >
C/C++ Source or Header  |  1996-09-28  |  6KB  |  304 lines
  1. // f-rand.cc                                           -*- C++ -*-
  2. /*
  3.  
  4. Copyright (C) 1993, 1994, 1995 John W. Eaton
  5.  
  6. This file is part of Octave.
  7.  
  8. Octave is free software; you can redistribute it and/or modify it
  9. under the terms of the GNU General Public License as published by the
  10. Free Software Foundation; either version 2, or (at your option) any
  11. later version.
  12.  
  13. Octave is distributed in the hope that it will be useful, but WITHOUT
  14. ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
  15. FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
  16. for more details.
  17.  
  18. You should have received a copy of the GNU General Public License
  19. along with Octave; see the file COPYING.  If not, write to the Free
  20. Software Foundation, 675 Mass Ave, Cambridge, MA 02139, USA.
  21.  
  22. */
  23.  
  24. #ifdef HAVE_CONFIG_H
  25. #include "config.h"
  26. #endif
  27.  
  28. #include <time.h>
  29.  
  30. #include "tree-const.h"
  31. #include "f77-uscore.h"
  32. #include "error.h"
  33. #include "gripes.h"
  34. #include "utils.h"
  35. #include "help.h"
  36. #include "defun-dld.h"
  37.  
  38. // Possible distributions of random numbers.
  39. enum rand_dist { uniform, normal };
  40.  
  41. // Current distribution of random numbers.
  42. static rand_dist current_distribution = uniform;
  43.  
  44. extern "C"
  45. {
  46.   int *F77_FCN (dgennor) (double*, double*, double*);
  47.   int *F77_FCN (dgenunf) (double*, double*, double*);
  48.   int *F77_FCN (setall) (int*, int*);
  49.   int *F77_FCN (getsd) (int*, int*);
  50. }
  51.  
  52. static double
  53. curr_rand_seed (void)
  54. {
  55.   union d2i { double d; int i[2]; };
  56.   union d2i u;
  57.   F77_FCN (getsd) (&(u.i[0]), &(u.i[1]));
  58.   return u.d;
  59. }
  60.  
  61. static int
  62. force_to_fit_range (int i, int lo, int hi)
  63. {
  64.   assert (hi > lo && lo >= 0 && hi > lo);
  65.  
  66.   i = i > 0 ? i : -i;
  67.  
  68.   if (i < lo)
  69.     i = lo;
  70.   else if (i > hi)
  71.     i = i % hi;
  72.  
  73.   return i;
  74. }
  75.  
  76. static void
  77. set_rand_seed (double val)
  78. {
  79.   union d2i { double d; int i[2]; };
  80.   union d2i u;
  81.   u.d = val;
  82.   int i0 = force_to_fit_range (u.i[0], 1, 2147483563);
  83.   int i1 = force_to_fit_range (u.i[1], 1, 2147483399);
  84.   F77_FCN (setall) (&i0, &i1);
  85. }
  86.  
  87. static char *
  88. curr_rand_dist (void)
  89. {
  90.   if (current_distribution == uniform)
  91.     return "uniform";
  92.   else if (current_distribution == normal)
  93.     return "normal";
  94.   else
  95.     {
  96.       panic_impossible ();
  97.       return 0;
  98.     }
  99. }
  100.  
  101. DEFUN_DLD_BUILTIN ("rand", Frand, Srand, 2, 1,
  102.   "rand                  -- generate a random value\n\
  103. \n\
  104. rand (N)              -- generate N x N matrix\n\
  105. rand (A)              -- generate matrix the size of A\n\
  106. rand (N, M)           -- generate N x M matrix\n\
  107. rand (\"dist\")         -- get current distribution\n\
  108. rand (DISTRIBUTION)   -- set distribution type (\"normal\" or \"uniform\"\n\
  109. rand (SEED)           -- get current seed\n\
  110. rand (SEED, N)        -- set seed")
  111. {
  112.   Octave_object retval;
  113.  
  114.   int nargin = args.length ();
  115.  
  116.   if (nargin > 2 || nargout > 1)
  117.     {
  118.       print_usage ("rand");
  119.       return retval;
  120.     }
  121.  
  122.   static int initialized = 0;
  123.   if (! initialized)
  124.     {
  125. // Make the random number generator give us a different sequence every
  126. // time we start octave unless we specifically set the seed.  The
  127. // technique used below will cycle monthly, but it it does seem to
  128. // work ok to give fairly different seeds each time Octave starts.
  129.  
  130. #if 0
  131.       int s0 = 1234567890;
  132.       int s1 = 123456789;
  133. #else
  134.       time_t now;
  135.       struct tm *tm;
  136.  
  137.       time (&now);
  138.       tm = localtime (&now);
  139.  
  140.       int hour = tm->tm_hour + 1;
  141.       int minute = tm->tm_min + 1;
  142.       int second = tm->tm_sec + 1;
  143.  
  144.       int s0 = tm->tm_mday * hour * minute * second;
  145.       int s1 = hour * minute * second;
  146. #endif
  147.       s0 = force_to_fit_range (s0, 1, 2147483563);
  148.       s1 = force_to_fit_range (s1, 1, 2147483399);
  149.  
  150.       F77_FCN (setall) (&s0, &s1);
  151.       initialized = 1;
  152.     }
  153.  
  154.   int n = 0;
  155.   int m = 0;
  156.   if (nargin == 0)
  157.     {
  158.       n = 1;
  159.       m = 1;
  160.       goto gen_matrix;
  161.     }
  162.   else if (nargin == 1)
  163.     {
  164.       tree_constant tmp = args(0);
  165.  
  166.       if (tmp.is_string ())
  167.     {
  168.       char *s_arg = tmp.string_value ();
  169.       if (strcmp (s_arg, "dist") == 0)
  170.         {
  171.           char *s = curr_rand_dist ();
  172.           retval(0) = s;
  173.         }
  174.       else if (strcmp (s_arg, "seed") == 0)
  175.         {
  176.           double d = curr_rand_seed ();
  177.           retval(0) = d;
  178.         }
  179.       else if (strcmp (s_arg, "uniform") == 0)
  180.         current_distribution = uniform;
  181.       else if (strcmp (s_arg, "normal") == 0)
  182.         current_distribution = normal;
  183.       else
  184.         error ("rand: unrecognized string argument");
  185.     }
  186.       else if (tmp.is_scalar_type ())
  187.     {
  188.       double dval = tmp.double_value ();
  189.  
  190.       if (xisnan (dval))
  191.         {
  192.           error ("rand: NaN is invalid a matrix dimension");
  193.         }
  194.       else
  195.         {
  196.           m = n = NINT (tmp.double_value ());
  197.  
  198.           if (! error_state)
  199.         goto gen_matrix;
  200.         }
  201.     }
  202.       else if (tmp.is_range ())
  203.     {
  204.       Range r = tmp.range_value ();
  205.       n = 1;
  206.       m = r.nelem ();
  207.       goto gen_matrix;
  208.     }
  209.       else if (tmp.is_matrix_type ())
  210.     {
  211.       n = args(0).rows ();
  212.       m = args(0).columns ();
  213.       goto gen_matrix;
  214.     }
  215.       else
  216.     {
  217.       gripe_wrong_type_arg ("rand", tmp);
  218.       return retval;
  219.     }
  220.     }
  221.   else if (nargin == 2)
  222.     {
  223.       if (args(0).is_string ()
  224.       && strcmp (args(0).string_value (), "seed") == 0)
  225.     {
  226.       double d = args(1).double_value ();
  227.  
  228.       if (! error_state)
  229.         set_rand_seed (d);
  230.     }
  231.       else
  232.     {
  233.       double dval = args(0).double_value ();
  234.  
  235.       if (xisnan (dval))
  236.         {
  237.           error ("rand: NaN is invalid as a matrix dimension");
  238.         }
  239.       else
  240.         {
  241.           n = NINT (dval);
  242.  
  243.           if (! error_state)
  244.         {
  245.           m = NINT (args(1).double_value ());
  246.  
  247.           if (! error_state)
  248.             goto gen_matrix;
  249.         }
  250.         }
  251.     }
  252.     }
  253.  
  254.   return retval;
  255.  
  256.  gen_matrix:
  257.  
  258.   if (n == 0 || m == 0)
  259.     {
  260.       Matrix m;
  261.       retval.resize (1, m);
  262.     }
  263.   else if (n > 0 && m > 0)
  264.     {
  265.       Matrix rand_mat (n, m);
  266.       for (int j = 0; j < m; j++)
  267.     for (int i = 0; i < n; i++)
  268.       {
  269.         double d_zero = 0.0;
  270.         double d_one = 1.0;
  271.         double val;
  272.         switch (current_distribution)
  273.           {
  274.           case uniform:
  275.         F77_FCN (dgenunf) (&d_zero, &d_one, &val);
  276.         rand_mat.elem (i, j) = val;
  277.         break;
  278.  
  279.           case normal:
  280.         F77_FCN (dgennor) (&d_zero, &d_one, &val);
  281.         rand_mat.elem (i, j) = val;
  282.         break;
  283.  
  284.           default:
  285.         panic_impossible ();
  286.         break;
  287.           }
  288.       }
  289.  
  290.       retval(0) = rand_mat;
  291.     }
  292.   else
  293.     error ("rand: invalid negative argument");
  294.  
  295.   return retval;
  296. }
  297.  
  298. /*
  299. ;;; Local Variables: ***
  300. ;;; mode: C++ ***
  301. ;;; page-delimiter: "^/\\*" ***
  302. ;;; End: ***
  303. */
  304.