Mercurial > hg > octave-lyh
changeset 14916:0b0569667939
maint: Periodic merge of default to jit
author | Jordi Gutiérrez Hermoso <jordigh@octave.org> |
---|---|
date | Wed, 23 May 2012 16:22:05 -0400 |
parents | 757f729fd41d (diff) cba58541954c (current diff) |
children | 232d8ab07932 |
files | configure.ac src/Makefile.am |
diffstat | 32 files changed, 966 insertions(+), 601 deletions(-) [+] |
line wrap: on
line diff
--- a/Makefile.am +++ b/Makefile.am @@ -86,7 +86,9 @@ CLEANFILES = $(BUILT_SOURCES) -DISTCLEANFILES = .gdbinit +DISTCLEANFILES = \ + .gdbinit \ + ChangeLog MAINTAINERCLEANFILES = $(BUILT_DISTFILES)
--- a/NEWS +++ b/NEWS @@ -59,6 +59,12 @@ Octave:array-as-scalar => Octave:array-to-scalar Octave:array-as-vector => Octave:array-to-vector + ** The colormap function now provides new options--"list", "register", + and "unregister"--to list all available colormap functions, and to + add or remove a function name from the list of known colormap + functions. Packages that implement extra colormaps should use these + commands with PKG_ADD and PKG_DEL statements. + ** Other new functions added in 3.8.0: colorcube splinefit @@ -89,12 +95,6 @@ static - ** The colormap function now provides new options "list", "register", - and "unregister" to list all available colormap functions, and to - add or remove a function name from teh list of known colormap - functions. Packages that implement extra colormaps should use these - commands with PKG_ADD and PKG_DEL statements. - Summary of important user-visible changes for version 3.6: ---------------------------------------------------------
--- a/configure.ac +++ b/configure.ac @@ -2430,8 +2430,9 @@ if $warn_msg_printed; then AC_MSG_NOTICE([]) - AC_MSG_NOTICE([NOTE: libraries may be skipped if a library is not found OR]) - AC_MSG_NOTICE([NOTE: if the library on your system is missing required features.]) + AC_MSG_NOTICE([NOTE: Libraries or auxiliary programs may be skipped if they are]) + AC_MSG_NOTICE([NOTE: not found OR if they are missing required features on your]) + AC_MSG_NOTICE([NOTE: system. ]) fi ### End of configure.
--- a/doc/interpreter/Makefile.am +++ b/doc/interpreter/Makefile.am @@ -198,6 +198,10 @@ ../../AUTHORS: preface.texi contributors.texi rm -f AUTHORS + if [ "x$(srcdir)" != "x." ] && [ -f $(srcdir)/contributors.texi ] && [ ! -f contributors.texi ]; then \ + cp $(srcdir)/contributors.texi contributors.texi; \ + touch -r $(srcdir)/contributors.texi contributors.texi; \ + fi -$(MAKEINFO) -D AUTHORSONLY \ --no-validate --no-headers --no-split --output AUTHORS $< mv AUTHORS ../../AUTHORS @@ -245,6 +249,7 @@ images \ images.awk \ images.mk \ + macros.texi \ mk_doc_cache.m \ mkcontrib.awk \ munge-texi.pl \
deleted file mode 100644 --- a/etc/README.devel +++ /dev/null @@ -1,82 +0,0 @@ -This directory contains development releases of Octave. - -If you want a stable, well-tested version of Octave, you should be -looking at ftp://ftp.gnu.org/gnu/octave. - -Development releases are provided for people who want to help test, -debug, and improve Octave. Very little testing is done before making -the development releases and they may even be made when Octave is in -an inconsistent state. It is possible that you will encounter a -very obvious bug, such as a failure to compile on *any* machine. It is -likely that such bugs will be fixed by the next development release, -so it really isn't necessary to report them unless they persist over -more than one release. - -Please DO report other bugs in the development releases as soon as you -find them. Bugs should be reported to the bug tracker at -'http://bugs.octave.org'. Please read read the bug reporting -guidelines (http://www.gnu.org/software/octave/bugs.html) before -submitting an item. - -If you have a fix for a bug, or an enhancement to submit, send your -patch to maintainers@octave.org or submit it to the patch -tracker at 'http://savannah.gnu.org/patch/?group=octave'. - -By adhering to the following guidelines you can minimize the work that -Octave maintainers need to do to apply your patch. Maintaining Octave -is a lot of work in the best of circumstances, and we can't keep up -unless you do your best to help. - - * Send an explanation with your changes of what problem they fix or - what improvement they bring about. For a bug fix, just include a - copy of the bug report, and explain why the change fixes the bug. - - * Always include a proper bug report for the problem you think you - have fixed. We need to convince ourselves that the change is - right before installing it. Even if it is right, we might have - trouble judging it if we don't have a way to reproduce the problem. - - * Include all the comments that are appropriate to help people - reading the source in the future understand why this change was - needed. - - * Don't mix together changes made for different reasons. Send them - _individually_. - - If you make two changes for separate reasons, then we might not - want to install them both. We might want to install just one. - - * Use `diff -c' to make your diffs. Diffs without context are hard - for us to install reliably. More than that, they make it hard for - us to study the diffs to decide whether we want to install them. - Unified diff format is better than contextless diffs, but not as - easy to read as `-c' format. - - If you have GNU diff, use `diff -cp', which shows the name of the - function that each change occurs in. - - * Write the change log entries for your changes. - - Read the `ChangeLog' file to see what sorts of information to put - in, and to learn the style that we use. The purpose of the - change log is to show people where to find what was changed. So - you need to be specific about what functions you changed; in - large functions, it's often helpful to indicate where within the - function the change was made. - - On the other hand, once you have shown people where to find the - change, you need not explain its purpose. Thus, if you add a new - function, all you need to say about it is that it is new. If you - feel that the purpose needs explaining, it probably does--but the - explanation will be much more useful if you put it in comments in - the code. - - If you would like your name to appear in the header line for who - made the change, send us the header line. - -If you would like to be on the very sharpest part of the bleeding -edge, you can now use Mercurial to access Octave's current development -sources. Instructions for checking out a copy are available on the -web at http://www.gnu.org/software/octave/download.html. - -Last updated: Sat Jan 22 21:26:18 PST 2011
deleted file mode 100644 --- a/etc/README.ftp +++ /dev/null @@ -1,38 +0,0 @@ -This directory contains the source for Octave, a high-level interactive -language for solving numerical problems. See the file README for more general -information, and the file NEWS for a list of recent changes. - -Binary distributions: --------------------- - -Octave binaries are not distributed from this site. - -Packaged versions of Octave for various GNU/Linux systems are available -with the major GNU/Linux distributions (Debian, Ubuntu, Red Hat, SUSE, etc.). - -Binary distributions of Octave for Mac OS X are available from both Fink -and MacPorts. - - http://www.finkproject.com - http://www.macports.org/ - -The file README.Windows provides instructions for installing Octave on -Windows systems. - -A note about .gz files: ----------------------- - -Files with names ending in `.gz' have been compressed with `gzip'. - -Unlike the compress utility, gzip is free of any known software -patents and tends to compress better anyway. Gzip can uncompress -`compress'-compressed files too, so you can install it as "uncompress" -and use it to handle both types of files. - -The gzip program is available in the directory /pub/gnu in shar, tar, -or gzipped tar format (for those who already have a prior version of -gzip and want faster data transmission). It works on virtually every -unix system, MSDOS, OS/2, and VMS. - - -Last updated: Thu Jan 20 10:14:49 PST 2011
deleted file mode 100644 --- a/etc/README.mirrors +++ /dev/null @@ -1,12 +0,0 @@ -If you have trouble transferring Octave from this site, you might try -one of those listed below. They mirror the /pub/octave directory on -ftp.octave.org - -FTP: - - site directory - ---- --------- - ftp.u-aizu.ac.jp /pub/SciEng/numanal/Octave - mirrors.fe.up.pt /pub/octave - -Last updated: Mon Jan 10 21:26:17 PST 2011
deleted file mode 100644 --- a/etc/README.snapshots +++ /dev/null @@ -1,161 +0,0 @@ -Octave Snapshots -- general info - -Snapshots are an "image" of the main Octave development tree, captured -at a particular random instant in time. When you use the snapshots, -you should be able to maintain a local copy of Octave that is -reasonably close to the official source tree used by the Octave -maintainers. - -The primary purpose of providing snapshots is to widen the group of -motivated developers that would like to help test, debug, and enhance -Octave, by providing you with access to the "latest and greatest" -source. This has several advantages, and several disadvantages. - - First the advantages: - - o Once we have a large base of motivated testers using the - snapshots, this should provide good coverage across all - currently supported Octave hosts and targets. If a new bug is - introduced in Octave due to fixing another bug or ongoing - development, it should become obvious much more quickly and - get fixed before the next general net release. This should - help to reduce the chances of Octave being released to the - general public with a major bug that went unnoticed during the - release cycle testing because they are machine dependent. We - hope to greatly improve Octave's stability and reliability by - involving more people and more execution environments in the - prerelease testing. - - o With access to the latest source, any diffs that you send to fix - bugs or add new features should be much easier for the Octave - maintainers to merge into the official source base (after - suitable review of course). This encourages us to merge your - changes quicker, while they are still "fresh". - - o Once your diffs are merged, you can obtain a new copy of Octave - containing your changes almost immediately. Thus you do not - have to maintain local copies of your changes for any longer - than it takes to get them merged into the official source - base. This encourages you to send in changes quicker. - - And the disadvantages: - - o The snapshot you get will be largely untested and of unknown - quality. It may fail to configure or compile. It may have - serious bugs. You should always keep a copy of the last known - working version before updating to the current snapshot, or at - least be able to regenerate a working version if the latest - snapshot is unusable in your environment for some reason. - - If a production version of Octave has a bug and a snapshot has - the fix, and you care about stability, you should put only the - fix for that particular problem into your production version. - Of course, if you are eager to test Octave, you can use the - snapshot versions in your daily work, but users who have not - been consulted about whether they feel like testing Octave should - generally have something which is at least as bug free as the - last released version. - - o Providing timely response to your questions, bug reports, and - submitted patches will require the Octave developers to - allocate time from an already thin time budget. Please try to - help us make this time as productive as possible. See the - section below about how to submit changes. - - -How to get the snapshots ------------------------- - -The current plan is to provide a full snapshot every week or so. For -now, diffs from previous versions will not be available. The files -will be available via anonymous ftp from ftp.octave.org, in the -directory /private/octave in the form of a tar files compressed with -GNU gzip. You can ftp gzip from ftp.octave.org in the directory -/pub/gnu. - -Even though the snapshots are available in a public place, we ask that -recipients not widely publicize the availability of the snapshots. -The motivation for this request is not to hoard them, but to avoid the -situation where the general Octave user base naively attempts to use -the snapshots, has trouble with them, complains publicly, and the -reputation of Octave declines because of a perception of instability -or lack of quality control. - - -Octave test suite ------------------ - -A test suite is distributed as an integral part of the snapshots. -However, to use it you will need to get a copy of the dejagnu testing -framework. Snapshots of dejagnu are available alongside the Octave -snapshots, using the same naming conventions as the Octave snapshots. -Once you have installed the dejagnu framework, a simple "make check" -in the Octave directory should be sufficient to run the tests. - -Note that the test suite is still quite limited. The test framework -itself might not install on your system if you have an environment -that is not similar to one that the Octave developers already use. -The tests themselves only cover a small portion of Octave features, -and what tests do exist for a feature are not exhaustive. New tests -are welcomed. - - -Bug reports ------------ - -Send bug reports to maintainers@octave.org. - -Note that since no testing is done on the snapshots, and snapshots may -even be made when Octave is in an inconsistent state, it may not be -unusual for an occasional snapshot to have a very obvious bug, such as -failure to compile on *any* machine. It is likely that such bugs will -be fixed by the next snapshot, so it really isn't necessary to report -them unless they persist over more than one snapshot. - -Missing files should always be reported, since they usually mean there -is a problem with the snapshot-generating process and we won't know -about them unless someone tells us. - -Bugs which are non-obvious, such as failure to compile on only a -specific machine, a new machine dependent or obscure bug (particularly -one not detected by the testsuite), etc. should be reported when you -discover them, or have a suggested patch to fix them. - - -FORMAT FOR PATCHES ------------------- - -If you have a fix for a bug, or an enhancement to submit, send your -patch to maintainers@octave.org. Here are some simple guidelines for -submitting patches: - - o Use "context diffs" for patches. A typical command for - generating context diffs is "diff -rc octave-old octave-new". - - o Use the "minimalist approach" for patches. That is, each patch - should address only one particular bug, new feature, etc. Do - not save up many unrelated changes and submit them all in one - big patch, since in general, the larger the patch the more - difficult it is for us to decide if the patch is either - correct or desirable. And if we find something about the - patch that needs to be corrected before it can be installed, - we would have to reject the entire patch, which might contain - changes which otherwise would be accepted if submitted - separately. - - o Submit a sample ChangeLog entry with your patch. See the - existing Octave ChangeLog for examples of what a ChangeLog - entry should look like. The emacs command ^X4A will create a - ChangeLog entry header for you. - - -Thanks, - -John W. Eaton -jwe@octave.org - -Wed, 31 Oct 2007 16:31:54 EDT - -This file was adapted from a similar document written by Fred Fish and -used by the GDB developers. -
--- a/libcruft/misc/lo-error.c +++ b/libcruft/misc/lo-error.c @@ -110,6 +110,7 @@ void liboctave_fatal_with_id (const char *id, const char *fmt, ...) { + (void) id; /*unused*/ va_list args; va_start (args, fmt); verror ("fatal", fmt, args); @@ -130,6 +131,7 @@ void liboctave_warning_with_id (const char *id, const char *fmt, ...) { + (void) id; /*unused*/ va_list args; va_start (args, fmt); verror ("warning", fmt, args);
--- a/liboctave/DiagArray2.h +++ b/liboctave/DiagArray2.h @@ -127,7 +127,7 @@ T operator () (octave_idx_type r, octave_idx_type c) const { #if defined (BOUNDS_CHECKING) - checkelem (r, c); + return checkelem (r, c); #else return elem (r, c); #endif
--- a/liboctave/Sparse.cc +++ b/liboctave/Sparse.cc @@ -235,7 +235,7 @@ (*current_liboctave_error_handler) ("Sparse::Sparse (const dim_vector&): dimension mismatch"); else - rep = new typename Sparse<T>::SparseRep (dv(0), dv(1)); + rep = new typename Sparse<T>::SparseRep (dv(0), dv(1), 0); } template <class T> @@ -301,11 +301,10 @@ (*current_liboctave_error_handler) ("sparse: column index %d out of bound %d", r.extent (nc), nc); - rep = new SparseRep (nr, nc); + rep = new typename Sparse<T>::SparseRep (nr, nc, (nzm > 0 ? nzm : 0)); dimensions = dim_vector (nr, nc); - octave_idx_type n = a.numel (), rl = r.length (nr), cl = c.length (nc); bool a_scalar = n == 1; if (a_scalar) @@ -324,6 +323,7 @@ if (n == 1 && a(0) != T ()) { change_capacity (nzm > 1 ? nzm : 1); + xcidx(0) = 0; xridx(0) = r(0); xdata(0) = a(0); for (octave_idx_type j = 0; j < nc; j++) @@ -352,6 +352,7 @@ new_nz += rd[i-1] != rd[i]; // Allocate result. change_capacity (nzm > new_nz ? nzm : new_nz); + xcidx (0) = 0; xcidx (1) = new_nz; octave_idx_type *rri = ridx (); T *rrd = data (); @@ -494,6 +495,7 @@ new_nz += rd[i-1] != rd[i]; // Allocate result. change_capacity (nzm > new_nz ? nzm : new_nz); + xcidx(0) = 0; xcidx(1) = new_nz; octave_idx_type *rri = ridx (); T *rrd = data ();
--- a/liboctave/oct-rand.cc +++ b/liboctave/oct-rand.cc @@ -419,19 +419,120 @@ return retval; } -Matrix -octave_rand::do_matrix (octave_idx_type n, octave_idx_type m, double a) +float +octave_rand::do_float_scalar (float a) { - Matrix retval; + float retval = 0.0; + + if (use_old_generators) + { + double da = a; + double dretval; + switch (current_distribution) + { + case uniform_dist: + F77_FUNC (dgenunf, DGENUNF) (0.0, 1.0, dretval); + break; + + case normal_dist: + F77_FUNC (dgennor, DGENNOR) (0.0, 1.0, dretval); + break; + + case expon_dist: + F77_FUNC (dgenexp, DGENEXP) (1.0, dretval); + break; - if (n >= 0 && m >= 0) - { - retval.clear (n, m); + case poisson_dist: + if (da < 0.0 || xisnan(da) || xisinf(da)) + dretval = octave_NaN; + else + { + // workaround bug in ignpoi, by calling with different Mu + F77_FUNC (dignpoi, DIGNPOI) (da + 1, dretval); + F77_FUNC (dignpoi, DIGNPOI) (da, dretval); + } + break; - if (n > 0 && m > 0) - fill (retval.capacity(), retval.fortran_vec(), a); + case gamma_dist: + if (da <= 0.0 || xisnan(da) || xisinf(da)) + retval = octave_NaN; + else + F77_FUNC (dgengam, DGENGAM) (1.0, da, dretval); + break; + + default: + (*current_liboctave_error_handler) + ("rand: invalid distribution ID = %d", current_distribution); + break; + } + retval = dretval; } else + { + switch (current_distribution) + { + case uniform_dist: + retval = oct_float_randu (); + break; + + case normal_dist: + retval = oct_float_randn (); + break; + + case expon_dist: + retval = oct_float_rande (); + break; + + case poisson_dist: + // Keep poisson distribution in double precision for accuracy + retval = oct_randp (a); + break; + + case gamma_dist: + retval = oct_float_randg (a); + break; + + default: + (*current_liboctave_error_handler) + ("rand: invalid distribution ID = %d", current_distribution); + break; + } + + save_state (); + } + + return retval; +} + +Array<double> +octave_rand::do_vector (octave_idx_type n, double a) +{ + Array<double> retval; + + if (n > 0) + { + retval.clear (n, 1); + + fill (retval.capacity(), retval.fortran_vec(), a); + } + else if (n < 0) + (*current_liboctave_error_handler) ("rand: invalid negative argument"); + + return retval; +} + +Array<float> +octave_rand::do_float_vector (octave_idx_type n, float a) +{ + Array<float> retval; + + if (n > 0) + { + retval.clear (n, 1); + + fill (retval.capacity(), retval.fortran_vec(), a); + } + else if (n < 0) (*current_liboctave_error_handler) ("rand: invalid negative argument"); return retval; @@ -452,19 +553,17 @@ return retval; } -Array<double> -octave_rand::do_vector (octave_idx_type n, double a) +FloatNDArray +octave_rand::do_float_nd_array (const dim_vector& dims, float a) { - Array<double> retval; + FloatNDArray retval; - if (n > 0) + if (! dims.all_zero ()) { - retval.clear (n, 1); + retval.clear (dims); - fill (retval.capacity (), retval.fortran_vec (), a); + fill (retval.capacity(), retval.fortran_vec(), a); } - else if (n < 0) - (*current_liboctave_error_handler) ("rand: invalid negative argument"); return retval; } @@ -693,3 +792,94 @@ return; } + +void +octave_rand::fill (octave_idx_type len, float *v, float a) +{ + if (len < 1) + return; + + switch (current_distribution) + { + case uniform_dist: + if (use_old_generators) + { +#define RAND_FUNC(x) F77_FUNC (dgenunf, DGENUNF) (0.0, 1.0, x) + MAKE_RAND (len); +#undef RAND_FUNC + } + else + oct_fill_float_randu (len, v); + break; + + case normal_dist: + if (use_old_generators) + { +#define RAND_FUNC(x) F77_FUNC (dgennor, DGENNOR) (0.0, 1.0, x) + MAKE_RAND (len); +#undef RAND_FUNC + } + else + oct_fill_float_randn (len, v); + break; + + case expon_dist: + if (use_old_generators) + { +#define RAND_FUNC(x) F77_FUNC (dgenexp, DGENEXP) (1.0, x) + MAKE_RAND (len); +#undef RAND_FUNC + } + else + oct_fill_float_rande (len, v); + break; + + case poisson_dist: + if (use_old_generators) + { + double da = a; + if (da < 0.0 || xisnan(da) || xisinf(da)) +#define RAND_FUNC(x) x = octave_NaN; + MAKE_RAND (len); +#undef RAND_FUNC + else + { + // workaround bug in ignpoi, by calling with different Mu + double tmp; + F77_FUNC (dignpoi, DIGNPOI) (da + 1, tmp); +#define RAND_FUNC(x) F77_FUNC (dignpoi, DIGNPOI) (da, x) + MAKE_RAND (len); +#undef RAND_FUNC + } + } + else + oct_fill_float_randp (a, len, v); + break; + + case gamma_dist: + if (use_old_generators) + { + double da = a; + if (da <= 0.0 || xisnan(da) || xisinf(da)) +#define RAND_FUNC(x) x = octave_NaN; + MAKE_RAND (len); +#undef RAND_FUNC + else +#define RAND_FUNC(x) F77_FUNC (dgengam, DGENGAM) (1.0, da, x) + MAKE_RAND (len); +#undef RAND_FUNC + } + else + oct_fill_float_randg (a, len, v); + break; + + default: + (*current_liboctave_error_handler) + ("rand: invalid distribution ID = %d", current_distribution); + break; + } + + save_state (); + + return; +}
--- a/liboctave/oct-rand.h +++ b/liboctave/oct-rand.h @@ -27,8 +27,8 @@ #include <string> #include "dColVector.h" -#include "dMatrix.h" #include "dNDArray.h" +#include "fNDArray.h" #include "lo-ieee.h" class @@ -136,13 +136,24 @@ return instance_ok () ? instance->do_scalar (a) : octave_NaN; } - // Return a matrix of numbers from the sequence, filled in column - // major order. - static Matrix matrix (octave_idx_type r, octave_idx_type c, double a = 1.0) + // Return the next number from the sequence. + static float float_scalar (float a = 1.0) { - return instance_ok () ? instance->do_matrix (r, c, a) : Matrix (); + return instance_ok () ? instance->do_float_scalar (a) : octave_Float_NaN; } + // Return an array of numbers from the sequence. + static Array<double> vector (octave_idx_type n, double a = 1.0) + { + return instance_ok () ? instance->do_vector (n, a) : Array<double> (); + } + + // Return an array of numbers from the sequence. + static Array<float> float_vector (octave_idx_type n, float a = 1.0) + { + return instance_ok () ? instance->do_float_vector (n, a) : Array<float> (); + } + // Return an N-dimensional array of numbers from the sequence, // filled in column major order. static NDArray nd_array (const dim_vector& dims, double a = 1.0) @@ -150,10 +161,12 @@ return instance_ok () ? instance->do_nd_array (dims, a) : NDArray (); } - // Return an array of numbers from the sequence. - static Array<double> vector (octave_idx_type n, double a = 1.0) + + // Return an N-dimensional array of numbers from the sequence, + // filled in column major order. + static FloatNDArray float_nd_array (const dim_vector& dims, float a = 1.0) { - return instance_ok () ? instance->do_vector (n, a) : Array<double> (); + return instance_ok () ? instance->do_float_nd_array (dims, a) : FloatNDArray (); } private: @@ -220,16 +233,22 @@ // Return the next number from the sequence. double do_scalar (double a = 1.); - // Return a matrix of numbers from the sequence, filled in column - // major order. - Matrix do_matrix (octave_idx_type r, octave_idx_type c, double a = 1.); + // Return the next number from the sequence. + float do_float_scalar (float a = 1.); + + // Return an array of numbers from the sequence. + Array<double> do_vector (octave_idx_type n, double a = 1.); + + // Return an array of numbers from the sequence. + Array<float> do_float_vector (octave_idx_type n, float a = 1.); // Return an N-dimensional array of numbers from the sequence, // filled in column major order. NDArray do_nd_array (const dim_vector& dims, double a = 1.); - // Return an array of numbers from the sequence. - Array<double> do_vector (octave_idx_type n, double a = 1.); + // Return an N-dimensional array of numbers from the sequence, + // filled in column major order. + FloatNDArray do_float_nd_array (const dim_vector& dims, float a = 1.); // Some helper functions. @@ -248,6 +267,8 @@ void switch_to_generator (int dist); void fill (octave_idx_type len, double *v, double a); + + void fill (octave_idx_type len, float *v, float a); }; #endif
--- a/liboctave/randgamma.c +++ b/liboctave/randgamma.c @@ -141,3 +141,59 @@ oct_fill_randg(a,1,&ret); return ret; } + +#undef NAN +#undef RUNI +#undef RNOR +#undef REXP +#define NAN octave_Float_NaN +#define RUNI oct_float_randu() +#define RNOR oct_float_randn() +#define REXP oct_float_rande() + +void +oct_fill_float_randg (float a, octave_idx_type n, float *r) +{ + octave_idx_type i; + /* If a < 1, start by generating gamma(1+a) */ + const float d = (a < 1. ? 1.+a : a) - 1./3.; + const float c = 1./sqrt(9.*d); + + /* Handle invalid cases */ + if (a <= 0 || INFINITE(a)) + { + for (i=0; i < n; i++) + r[i] = NAN; + return; + } + + for (i=0; i < n; i++) + { + float x, xsq, v, u; + frestart: + x = RNOR; + v = (1+c*x); + v *= v*v; + if (v <= 0) + goto frestart; /* rare, so don't bother moving up */ + u = RUNI; + xsq = x*x; + if (u >= 1.-0.0331*xsq*xsq && log(u) >= 0.5*xsq + d*(1-v+log(v))) + goto frestart; + r[i] = d*v; + } + if (a < 1) + { /* Use gamma(a) = gamma(1+a)*U^(1/a) */ + /* Given REXP = -log(U) then U^(1/a) = exp(-REXP/a) */ + for (i = 0; i < n; i++) + r[i] *= exp(-REXP/a); + } +} + +float +oct_float_randg (float a) +{ + float ret; + oct_fill_float_randg(a,1,&ret); + return ret; +}
--- a/liboctave/randgamma.h +++ b/liboctave/randgamma.h @@ -32,6 +32,9 @@ extern OCTAVE_API double oct_randg (double a); extern OCTAVE_API void oct_fill_randg (double a, octave_idx_type n, double *p); +extern OCTAVE_API float oct_float_randg (float a); +extern OCTAVE_API void oct_fill_float_randg (float a, octave_idx_type n, float *p); + #ifdef __cplusplus } #endif
--- a/liboctave/randmtzig.c +++ b/liboctave/randmtzig.c @@ -64,10 +64,6 @@ http://www.math.keio.ac.jp/matumoto/emt.html email: matumoto@math.keio.ac.jp - * 2006-04-01 David Bateman - * * convert for use in octave, declaring static functions only used - * here and adding oct_ to functions visible externally - * * inverse sense of ALLBITS * 2004-01-19 Paul Kienzle * * comment out main * add init_by_entropy, get_state, set_state @@ -87,6 +83,15 @@ * * 2005-02-23 Paul Kienzle * * fix -DHAVE_X86_32 flag and add -DUSE_X86_32=0|1 for explicit control + * + * 2006-04-01 David Bateman + * * convert for use in octave, declaring static functions only used + * here and adding oct_ to functions visible externally + * * inverse sense of ALLBITS + * + * 2012-05-18 David Bateman + * * Remove randu64 and ALLBIT option + * * Add the single precision generators */ /* @@ -96,9 +101,6 @@ available. This is not necessary if your architecture has /dev/urandom defined. - Compile with -DALLBITS to disable 53-bit random numbers. This is about - 50% slower than using 32-bit random numbers. - Uses implicit -Di386 or explicit -DHAVE_X86_32 to determine if CPU=x86. You can force X86 behaviour with -DUSE_X86_32=1, or suppress it with -DUSE_X86_32=0. You should also consider -march=i686 or similar for @@ -126,21 +128,28 @@ static uint32_t randi32(void) returns 32-bit unsigned int static uint64_t randi53(void) returns 53-bit unsigned int static uint64_t randi54(void) returns 54-bit unsigned int - static uint64_t randi64(void) returns 64-bit unsigned int - static double randu32(void) returns 32-bit uniform in (0,1) + static float randu32(void) returns 32-bit uniform in (0,1) static double randu53(void) returns 53-bit uniform in (0,1) double oct_randu(void) returns M-bit uniform in (0,1) double oct_randn(void) returns M-bit standard normal double oct_rande(void) returns N-bit standard exponential + float oct_float_randu(void) returns M-bit uniform in (0,1) + float oct_float_randn(void) returns M-bit standard normal + float oct_float_rande(void) returns N-bit standard exponential + === Array generators === void oct_fill_randi32(octave_idx_type, uint32_t []) void oct_fill_randi64(octave_idx_type, uint64_t []) + void oct_fill_randu(octave_idx_type, double []) void oct_fill_randn(octave_idx_type, double []) void oct_fill_rande(octave_idx_type, double []) + void oct_fill_float_randu(octave_idx_type, float []) + void oct_fill_float_randn(octave_idx_type, float []) + void oct_fill_float_rande(octave_idx_type, float []) */ #if defined (HAVE_CONFIG_H) @@ -180,6 +189,7 @@ static int left = 1; static int initf = 0; static int initt = 1; +static int inittf = 1; /* initializes state[MT_N] with a seed */ void @@ -378,34 +388,14 @@ #endif } -#if 0 -// FIXME -- this doesn't seem to be used anywhere; should it be removed? -static uint64_t -randi64 (void) -{ - const uint32_t lo = randi32(); - const uint32_t hi = randi32(); -#if HAVE_X86_32 - uint64_t u; - uint32_t *p = (uint32_t *)&u; - p[0] = lo; - p[1] = hi; - return u; -#else - return (((uint64_t)hi<<32)|lo); -#endif -} -#endif - -#ifdef ALLBITS /* generates a random number on (0,1)-real-interval */ -static double +static float randu32 (void) { - return ((double)randi32() + 0.5) * (1.0/4294967296.0); + return ((float)randi32() + 0.5) * (1.0/4294967296.0); /* divided by 2^32 */ } -#else + /* generates a random number on (0,1) with 53-bit resolution */ static double randu53 (void) @@ -414,35 +404,22 @@ const uint32_t b=randi32()>>6; return (a*67108864.0+b+0.4) * (1.0/9007199254740992.0); } -#endif /* Determine mantissa for uniform doubles */ double oct_randu (void) { -#ifdef ALLBITS + return randu53 (); +} + +/* Determine mantissa for uniform floats */ +float +oct_float_randu (void) +{ return randu32 (); -#else - return randu53 (); -#endif } /* ===== Ziggurat normal and exponential generators ===== */ -#ifdef ALLBITS -# define ZIGINT uint32_t -# define EMANTISSA 4294967296.0 /* 32 bit mantissa */ -# define ERANDI randi32() /* 32 bits for mantissa */ -# define NMANTISSA 2147483648.0 /* 31 bit mantissa */ -# define NRANDI randi32() /* 31 bits for mantissa + 1 bit sign */ -# define RANDU randu32() -#else -# define ZIGINT uint64_t -# define EMANTISSA 9007199254740992.0 /* 53 bit mantissa */ -# define ERANDI randi53() /* 53 bits for mantissa */ -# define NMANTISSA EMANTISSA -# define NRANDI randi54() /* 53 bits for mantissa + 1 bit sign */ -# define RANDU randu53() -#endif #define ZIGGURAT_TABLE_SIZE 256 @@ -454,6 +431,13 @@ #define ZIGGURAT_EXP_INV_R 0.129918765548341586 #define EXP_SECTION_AREA 0.0039496598225815571993 +#define ZIGINT uint64_t +#define EMANTISSA 9007199254740992.0 /* 53 bit mantissa */ +#define ERANDI randi53() /* 53 bits for mantissa */ +#define NMANTISSA EMANTISSA +#define NRANDI randi54() /* 53 bits for mantissa + 1 bit sign */ +#define RANDU randu53() + static ZIGINT ki[ZIGGURAT_TABLE_SIZE]; static double wi[ZIGGURAT_TABLE_SIZE], fi[ZIGGURAT_TABLE_SIZE]; static ZIGINT ke[ZIGGURAT_TABLE_SIZE]; @@ -592,7 +576,6 @@ * Of course, different compilers and operating systems may * have something to do with this. */ -#if !defined(ALLBITS) # if HAVE_X86_32 /* 53-bit mantissa, 1-bit sign, x86 32-bit architecture */ double x; @@ -615,14 +598,6 @@ const double x = ( r&1 ? -rabs : rabs) * wi[idx]; # endif /* !HAVE_X86_32 */ if (rabs < (int64_t)ki[idx]) -#else /* ALLBITS */ - /* 32-bit mantissa */ - const uint32_t r = randi32(); - const uint32_t rabs = r&LMASK; - const int idx = (int)(r&0xFF); - const double x = ((int32_t)r) * wi[idx]; - if (rabs < ki[idx]) -#endif /* ALLBITS */ return x; /* 99.3% of the time we return here 1st try */ else if (idx == 0) { @@ -677,6 +652,173 @@ } } +#undef ZIGINT +#undef EMANTISSA +#undef ERANDI +#undef NMANTISSA +#undef NRANDI +#undef RANDU + +#define ZIGINT uint32_t +#define EMANTISSA 4294967296.0 /* 32 bit mantissa */ +#define ERANDI randi32() /* 32 bits for mantissa */ +#define NMANTISSA 2147483648.0 /* 31 bit mantissa */ +#define NRANDI randi32() /* 31 bits for mantissa + 1 bit sign */ +#define RANDU randu32() + +static ZIGINT fki[ZIGGURAT_TABLE_SIZE]; +static float fwi[ZIGGURAT_TABLE_SIZE], ffi[ZIGGURAT_TABLE_SIZE]; +static ZIGINT fke[ZIGGURAT_TABLE_SIZE]; +static float fwe[ZIGGURAT_TABLE_SIZE], ffe[ZIGGURAT_TABLE_SIZE]; + + +static void +create_ziggurat_float_tables (void) +{ + int i; + float x, x1; + + /* Ziggurat tables for the normal distribution */ + x1 = ZIGGURAT_NOR_R; + fwi[255] = x1 / NMANTISSA; + ffi[255] = exp (-0.5 * x1 * x1); + + /* Index zero is special for tail strip, where Marsaglia and Tsang + * defines this as + * k_0 = 2^31 * r * f(r) / v, w_0 = 0.5^31 * v / f(r), f_0 = 1, + * where v is the area of each strip of the ziggurat. + */ + fki[0] = (ZIGINT) (x1 * ffi[255] / NOR_SECTION_AREA * NMANTISSA); + fwi[0] = NOR_SECTION_AREA / ffi[255] / NMANTISSA; + ffi[0] = 1.; + + for (i = 254; i > 0; i--) + { + /* New x is given by x = f^{-1}(v/x_{i+1} + f(x_{i+1})), thus + * need inverse operator of y = exp(-0.5*x*x) -> x = sqrt(-2*ln(y)) + */ + x = sqrt(-2. * log(NOR_SECTION_AREA / x1 + ffi[i+1])); + fki[i+1] = (ZIGINT)(x / x1 * NMANTISSA); + fwi[i] = x / NMANTISSA; + ffi[i] = exp (-0.5 * x * x); + x1 = x; + } + + fki[1] = 0; + + /* Zigurrat tables for the exponential distribution */ + x1 = ZIGGURAT_EXP_R; + fwe[255] = x1 / EMANTISSA; + ffe[255] = exp (-x1); + + /* Index zero is special for tail strip, where Marsaglia and Tsang + * defines this as + * k_0 = 2^32 * r * f(r) / v, w_0 = 0.5^32 * v / f(r), f_0 = 1, + * where v is the area of each strip of the ziggurat. + */ + fke[0] = (ZIGINT) (x1 * ffe[255] / EXP_SECTION_AREA * EMANTISSA); + fwe[0] = EXP_SECTION_AREA / ffe[255] / EMANTISSA; + ffe[0] = 1.; + + for (i = 254; i > 0; i--) + { + /* New x is given by x = f^{-1}(v/x_{i+1} + f(x_{i+1})), thus + * need inverse operator of y = exp(-x) -> x = -ln(y) + */ + x = - log(EXP_SECTION_AREA / x1 + ffe[i+1]); + fke[i+1] = (ZIGINT)(x / x1 * EMANTISSA); + fwe[i] = x / EMANTISSA; + ffe[i] = exp (-x); + x1 = x; + } + fke[1] = 0; + + inittf = 0; +} + +/* + * Here is the guts of the algorithm. As Marsaglia and Tsang state the + * algorithm in their paper + * + * 1) Calculate a random signed integer j and let i be the index + * provided by the rightmost 8-bits of j + * 2) Set x = j * w_i. If j < k_i return x + * 3) If i = 0, then return x from the tail + * 4) If [f(x_{i-1}) - f(x_i)] * U < f(x) - f(x_i), return x + * 5) goto step 1 + * + * Where f is the functional form of the distribution, which for a normal + * distribution is exp(-0.5*x*x) + */ + +float +oct_float_randn (void) +{ + if (inittf) + create_ziggurat_float_tables(); + + while (1) + { + /* 32-bit mantissa */ + const uint32_t r = randi32(); + const uint32_t rabs = r&LMASK; + const int idx = (int)(r&0xFF); + const float x = ((int32_t)r) * fwi[idx]; + if (rabs < fki[idx]) + return x; /* 99.3% of the time we return here 1st try */ + else if (idx == 0) + { + /* As stated in Marsaglia and Tsang + * + * For the normal tail, the method of Marsaglia[5] provides: + * generate x = -ln(U_1)/r, y = -ln(U_2), until y+y > x*x, + * then return r+x. Except that r+x is always in the positive + * tail!!!! Any thing random might be used to determine the + * sign, but as we already have r we might as well use it + * + * [PAK] but not the bottom 8 bits, since they are all 0 here! + */ + float xx, yy; + do + { + xx = - ZIGGURAT_NOR_INV_R * log (RANDU); + yy = - log (RANDU); + } + while ( yy+yy <= xx*xx); + return (rabs&0x100 ? -ZIGGURAT_NOR_R-xx : ZIGGURAT_NOR_R+xx); + } + else if ((ffi[idx-1] - ffi[idx]) * RANDU + ffi[idx] < exp(-0.5*x*x)) + return x; + } +} + +float +oct_float_rande (void) +{ + if (inittf) + create_ziggurat_float_tables(); + + while (1) + { + ZIGINT ri = ERANDI; + const int idx = (int)(ri & 0xFF); + const float x = ri * fwe[idx]; + if (ri < fke[idx]) + return x; // 98.9% of the time we return here 1st try + else if (idx == 0) + { + /* As stated in Marsaglia and Tsang + * + * For the exponential tail, the method of Marsaglia[5] provides: + * x = r - ln(U); + */ + return ZIGGURAT_EXP_R - log(RANDU); + } + else if ((ffe[idx-1] - ffe[idx]) * RANDU + ffe[idx] < exp(-x)) + return x; + } +} + /* Array generators */ void oct_fill_randu (octave_idx_type n, double *p) @@ -701,3 +843,27 @@ for (i = 0; i < n; i++) p[i] = oct_rande(); } + +void +oct_fill_float_randu (octave_idx_type n, float *p) +{ + octave_idx_type i; + for (i = 0; i < n; i++) + p[i] = oct_float_randu(); +} + +void +oct_fill_float_randn (octave_idx_type n, float *p) +{ + octave_idx_type i; + for (i = 0; i < n; i++) + p[i] = oct_float_randn(); +} + +void +oct_fill_float_rande (octave_idx_type n, float *p) +{ + octave_idx_type i; + for (i = 0; i < n; i++) + p[i] = oct_float_rande(); +}
--- a/liboctave/randmtzig.h +++ b/liboctave/randmtzig.h @@ -82,11 +82,19 @@ extern OCTAVE_API double oct_randn (void); extern OCTAVE_API double oct_rande (void); +extern OCTAVE_API float oct_float_randu (void); +extern OCTAVE_API float oct_float_randn (void); +extern OCTAVE_API float oct_float_rande (void); + /* === Array generators === */ extern OCTAVE_API void oct_fill_randu (octave_idx_type n, double *p); extern OCTAVE_API void oct_fill_randn (octave_idx_type n, double *p); extern OCTAVE_API void oct_fill_rande (octave_idx_type n, double *p); +extern OCTAVE_API void oct_fill_float_randu (octave_idx_type n, float *p); +extern OCTAVE_API void oct_fill_float_randn (octave_idx_type n, float *p); +extern OCTAVE_API void oct_fill_float_rande (octave_idx_type n, float *p); + #ifdef __cplusplus } #endif
--- a/liboctave/randpoisson.c +++ b/liboctave/randpoisson.c @@ -368,6 +368,46 @@ } } +static void +poisson_cdf_lookup_float(double lambda, float *p, size_t n) +{ + double t[TABLESIZE]; + + /* Precompute the table for the u up to and including 0.458. + * We will almost certainly need it. */ + int intlambda = (int)floor(lambda); + double P; + int tableidx; + size_t i = n; + + t[0] = P = exp(-lambda); + for (tableidx = 1; tableidx <= intlambda; tableidx++) { + P = P*lambda/(double)tableidx; + t[tableidx] = t[tableidx-1] + P; + } + + while (i-- > 0) { + double u = RUNI; + int k = (u > 0.458 ? intlambda : 0); + nextk: + if ( u <= t[k] ) { + p[i] = (float) k; + continue; + } + if (++k < tableidx) goto nextk; + + while (tableidx < TABLESIZE) { + P = P*lambda/(double)tableidx; + t[tableidx] = t[tableidx-1] + P; + if (t[tableidx] == t[tableidx-1]) t[tableidx] = 1.0; + tableidx++; + if (u <= t[tableidx-1]) break; + } + + p[i] = (float)(tableidx-1); + } +} + /* From Press, et al., Numerical Recipes */ static void poisson_rejection (double lambda, double *p, size_t n) @@ -392,6 +432,30 @@ } } +/* From Press, et al., Numerical Recipes */ +static void +poisson_rejection_float (double lambda, float *p, size_t n) +{ + double sq = sqrt(2.0*lambda); + double alxm = log(lambda); + double g = lambda*alxm - LGAMMA(lambda+1.0); + size_t i; + + for (i = 0; i < n; i++) + { + double y, em, t; + do { + do { + y = tan(M_PI*RUNI); + em = sq * y + lambda; + } while (em < 0.0); + em = floor(em); + t = 0.9*(1.0+y*y)*exp(em*alxm-flogfak(em)-g); + } while (RUNI > t); + p[i] = em; + } +} + /* The cutoff of L <= 1e8 in the following two functions before using * the normal approximation is based on: * > L=1e8; x=floor(linspace(0,2*L,1000)); @@ -463,3 +527,68 @@ } return ret; } + +/* Generate a set of poisson numbers with the same distribution */ +void +oct_fill_float_randp (float FL, octave_idx_type n, float *p) +{ + double L = FL; + octave_idx_type i; + if (L < 0.0 || INFINITE(L)) + { + for (i=0; i<n; i++) + p[i] = NAN; + } + else if (L <= 10.0) + { + poisson_cdf_lookup_float(L, p, n); + } + else if (L <= 1e8) + { + for (i=0; i<n; i++) + p[i] = pprsc(L); + } + else + { + /* normal approximation: from Phys. Rev. D (1994) v50 p1284 */ + const double sqrtL = sqrt(L); + for (i = 0; i < n; i++) + { + p[i] = floor(RNOR*sqrtL + L + 0.5); + if (p[i] < 0.0) + p[i] = 0.0; /* will probably never happen */ + } + } +} + +/* Generate one poisson variate */ +float +oct_float_randp (float FL) +{ + double L = FL; + float ret; + if (L < 0.0) ret = NAN; + else if (L <= 12.0) { + /* From Press, et al. Numerical recipes */ + double g = exp(-L); + int em = -1; + double t = 1.0; + do { + ++em; + t *= RUNI; + } while (t > g); + ret = em; + } else if (L <= 1e8) { + /* numerical recipes */ + poisson_rejection_float(L, &ret, 1); + } else if (INFINITE(L)) { + /* FIXME R uses NaN, but the normal approx. suggests that as + * limit should be inf. Which is correct? */ + ret = NAN; + } else { + /* normal approximation: from Phys. Rev. D (1994) v50 p1284 */ + ret = floor(RNOR*sqrt(L) + L + 0.5); + if (ret < 0.0) ret = 0.0; /* will probably never happen */ + } + return ret; +}
--- a/liboctave/randpoisson.h +++ b/liboctave/randpoisson.h @@ -32,6 +32,9 @@ extern OCTAVE_API double oct_randp (double L); extern OCTAVE_API void oct_fill_randp (double L, octave_idx_type n, double *p); +extern OCTAVE_API float oct_float_randp (float L); +extern OCTAVE_API void oct_fill_float_randp (float L, octave_idx_type n, float *p); + #ifdef __cplusplus } #endif
--- a/m4/acinclude.m4 +++ b/m4/acinclude.m4 @@ -522,8 +522,12 @@ ;; *) LEX='$(top_srcdir)/build-aux/missing flex' - warn_flex="I didn't find flex, but it's only a problem if you need to reconstruct lex.cc" - AC_MSG_WARN([$warn_flex]) + warn_flex=" + +I didn't find flex, but it's only a problem if you need to reconstruct +lex.cc, which is the case if you're building from VCS sources. +" + OCTAVE_CONFIGURE_WARNING([warn_flex]) ;; esac AC_SUBST(LFLAGS) @@ -538,8 +542,14 @@ ;; *) YACC='$(top_srcdir)/build-aux/missing bison' - warn_bison="I didn't find bison, but it's only a problem if you need to reconstruct parse.cc" - AC_MSG_WARN([$warn_bison]) + warn_bison=" + +I didn't find bison, but it's only a problem if you need to +reconstruct parse.cc, which is the case if you're building from VCS +sources. + +" + OCTAVE_CONFIGURE_WARNING([warn_bison]) ;; esac ]) @@ -562,7 +572,7 @@ AC_CHECK_PROGS(DEFAULT_PAGER, $octave_possible_pagers, []) if test -z "$DEFAULT_PAGER"; then warn_less="I couldn't find \`less', \`more', \`page', or \`pg'" - AC_MSG_WARN([$warn_less]) + OCTAVE_CONFIGURE_WARNING([warn_less]) fi fi ]) @@ -581,10 +591,14 @@ warn_gnuplot=yes GNUPLOT="$gp_default" + warn_gnuplot = " - AC_MSG_WARN([gnuplot not found. It isn't necessary to have gnuplot]) - AC_MSG_WARN([installed, but without native graphics or gnuplot]) - AC_MSG_WARN([you won't be able to use any of Octave's plotting commands.]) +gnuplot not found. It isn't necessary to have gnuplot installed, but +without native graphics or gnuplot you won't be able to use any of +Octave's plotting commands. + +" + OCTAVE_CONFIGURE_WARNING([warn_gnuplot]) fi fi AC_SUBST(GNUPLOT) @@ -597,8 +611,12 @@ AC_CHECK_PROG(GPERF, gperf, gperf, []) if test -z "$GPERF"; then GPERF='$(top_srcdir)/build-aux/missing gperf' - warn_gperf="I didn't find gperf, but it's only a problem if you need to reconstruct oct-gperf.h" - AC_MSG_WARN([$warn_gperf]) + warn_gperf=" + +I didn't find gperf, but it's only a problem if you need to +reconstruct oct-gperf.h +" + OCTAVE_CONFIGURE_WARNING([warn_gperf]) fi AC_SUBST(GPERF) ]) @@ -618,8 +636,14 @@ AC_CHECK_PROGS(GHOSTSCRIPT, [$gs_names]) if test -z "$GHOSTSCRIPT"; then GHOSTSCRIPT='$(top_srcdir)/build-aux/missing gs' - warn_ghostscript="I didn't find ghostscript, so reconstructing figures for the manual will fail, and saving graphics in some output formats will fail when using Octave" - AC_MSG_WARN([$warn_ghostscript]) + warn_ghostscript=" + +I didn't find ghostscript, so reconstructing figures for the manual +will fail, and saving graphics in some output formats will fail when +using Octave +" + + OCTAVE_CONFIGURE_WARNING([warn_ghostscript]) fi AC_SUBST(GHOSTSCRIPT) ]) @@ -631,8 +655,12 @@ AC_CHECK_PROG(TEXI2DVI, texi2dvi, texi2dvi, []) if test -z "$TEXI2DVI"; then TEXI2DVI='$(top_srcdir)/build-aux/missing texi2dvi' - warn_texi2dvi="I didn't find texi2dvi, but it's only a problem if you need to reconstruct the DVI version of the manual" - AC_MSG_WARN([$warn_texi2dvi]) + warn_texi2dvi=" + +I didn't find texi2dvi, but it's only a problem if you need to +reconstruct the DVI version of the manual +" + OCTAVE_CONFIGURE_WARNING([warn_texi2dvi]) fi AC_SUBST(TEXI2DVI) ]) @@ -654,8 +682,12 @@ fi if $missing; then TEXI2PDF='$(top_srcdir)/build-aux/missing texi2pdf' - warn_texi2pdf="I didn't find texi2pdf, but it's only a problem if you need to reconstruct the PDF version of the manual" - AC_MSG_WARN([$warn_texi2pdf]) + warn_texi2pdf=" + +I didn't find texi2pdf, but it's only a problem if you need to +reconstruct the PDF version of the manual +" + OCTAVE_CONFIGURE_WARNING([warn_texi2pdf]) fi AC_SUBST(TEXI2PDF) ]) @@ -1501,31 +1533,33 @@ else AC_MSG_RESULT([no]) if test "$FT2_CONFIG" = "no" ; then - AC_MSG_WARN([ + warn_ft2_config = " - The freetype-config script installed by FreeType 2 could not be found. - If FreeType 2 was installed in PREFIX, make sure PREFIX/bin is in - your path, or set the FT2_CONFIG environment variable to the - full path to freetype-config. - ]) +The freetype-config script installed by FreeType 2 could not be found. +If FreeType 2 was installed in PREFIX, make sure PREFIX/bin is in your +path, or set the FT2_CONFIG environment variable to the full path to +freetype-config. +" + OCTAVE_CONFIGURE_WARNING([warn_ft2_config]) else if test x$ft_config_is_lt = xyes ; then - AC_MSG_WARN([ + warn_ft2_too_old=" - Your installed version of the FreeType 2 library is too old. - If you have different versions of FreeType 2, make sure that - correct values for --with-ft-prefix or --with-ft-exec-prefix - are used, or set the FT2_CONFIG environment variable to the - full path to freetype-config. - ]) +Your installed version of the FreeType 2 library is too old. If you +have different versions of FreeType 2, make sure that correct values +for --with-ft-prefix or --with-ft-exec-prefix are used, or set the +FT2_CONFIG environment variable to the full path to freetype-config. +" + OCTAVE_CONFIGURE_WARNING([warn_ft2_too_old]) else - AC_MSG_WARN([ + warn_ft2_failed=" - The FreeType test program failed to run. If your system uses - shared libraries and they are installed outside the normal - system library path, make sure the variable LD_LIBRARY_PATH - (or whatever is appropiate for your system) is correctly set. - ]) +The FreeType test program failed to run. If your system uses shared +libraries and they are installed outside the normal system library +path, make sure the variable LD_LIBRARY_PATH (or whatever is +appropiate for your system) is correctly set. +" + OCTAVE_CONFIGURE_WARNING([warn_ft2_failed]) fi fi
--- a/scripts/pkg/private/load_package_dirs.m +++ b/scripts/pkg/private/load_package_dirs.m @@ -28,6 +28,11 @@ && installed_pkgs_lst{i}.loaded) continue; else + ## Insert this package at the front before recursing over dependencies. + if (! any (idx == i)) + idx = [i, idx]; + endif + if (handle_deps) deps = installed_pkgs_lst{i}.depends; if ((length (deps) > 1) @@ -36,8 +41,10 @@ for k = 1 : length (deps) for j = 1 : length (installed_pkgs_lst) if (strcmp (installed_pkgs_lst{j}.name, deps{k}.package)) - tmplidx (end + 1) = j; - break; + if (! any (idx == j)) + tmplidx (end + 1) = j; + break; + endif endif endfor endfor @@ -45,9 +52,6 @@ installed_pkgs_lst); endif endif - if (isempty (find(idx == i))) - idx (end + 1) = i; - endif endif endfor endfunction
--- a/scripts/plot/private/__go_draw_figure__.m +++ b/scripts/plot/private/__go_draw_figure__.m @@ -44,6 +44,7 @@ else bg_is_set = false; endif + fg_was_set = false; for i = nkids:-1:1 type = get (kids(i), "type"); @@ -142,6 +143,11 @@ if (isnumeric (fg) && strcmp (get (kids(i), "visible"), "on")) fprintf (plot_stream, "set obj 2 rectangle from graph 0,0 to graph 1,1 behind fc rgb \"#%02x%02x%02x\"\n", 255 * fg); fg_is_set = true; + fg_was_set = true; + elseif (fg_was_set) + fprintf (plot_stream, "unset obj 2\n"); + fg_is_set = false; + fg_was_set = false; else fg_is_set = false; endif
--- a/scripts/polynomial/module.mk +++ b/scripts/polynomial/module.mk @@ -1,7 +1,7 @@ FCN_FILE_DIRS += polynomial polynomial_PRIVATE_FCN_FILES = \ - polynomial/private/__splinefit__.m + polynomial/private/__splinefit__.m polynomial_FCN_FILES = \ polynomial/compan.m \ @@ -28,7 +28,8 @@ polynomial/roots.m \ polynomial/spline.m \ polynomial/splinefit.m \ - polynomial/unmkpp.m + polynomial/unmkpp.m \ + $(polynomial_PRIVATE_FCN_FILES) FCN_FILES += $(polynomial_FCN_FILES)
--- a/scripts/sparse/spdiags.m +++ b/scripts/sparse/spdiags.m @@ -79,8 +79,8 @@ ## Create new matrix of size mxn using v,c [j, i, v] = find (v); offset = max (min (c(:), n-m), 0); - j = j + offset(i); - i = j-c(:)(i); + j = j(:) + offset(i(:)); + i = j - c(:)(i(:)); idx = i > 0 & i <= m & j > 0 & j <= n; A = sparse (i(idx), j(idx), v(idx), m, n); endif @@ -90,4 +90,4 @@ %!assert (spdiags (zeros (1,0),1,1,1), sparse (0)) %!assert (spdiags (zeros (0,1),1,1,1), sparse (0)) - +%!assert (spdiags ([0.5 -1 0.5], 0:2, 1, 1), sparse(0.5))
--- a/scripts/time/datenum.m +++ b/scripts/time/datenum.m @@ -85,6 +85,7 @@ ## Days until start of month assuming year starts March 1. persistent monthstart = [306; 337; 0; 31; 61; 92; 122; 153; 184; 214; 245; 275]; + persistent monthlength = [31; 28; 31; 30; 31; 30; 31; 31; 30; 31; 30; 31]; if (nargin == 0 || nargin > 6 || (nargin > 2 && (ischar (year) || iscellstr (year)))) @@ -110,6 +111,19 @@ month(month<1) = 1; ## For compatibility. Otherwise allow negative months. + ## Treat fractional months, by converting the fraction to days + if (floor (month) != month) + fracmonth = month - floor (month); + month = floor (month); + if ((mod (month-1,12) + 1) == 2 && + (floor (year/4) - floor (year/100) + floor (year/400)) != 0) + ## leap year + day += fracmonth * 29; + else + day += fracmonth * monthlength ((mod (month-1,12) + 1)); + endif + endif + ## Set start of year to March by moving Jan. and Feb. to previous year. ## Correct for months > 12 by moving to subsequent years. year += fix ((month-14)/12);
--- a/scripts/time/datetick.m +++ b/scripts/time/datetick.m @@ -159,7 +159,7 @@ ## Day scale or less if (xmax - xmin < N / 24 / 60 / 60) scl = 1 / 24 / 60 / 60; - elseif (xmax - xmin < N / 24 / 6) + elseif (xmax - xmin < N / 24 / 60) scl = 1 / 24 / 60; else scl = 1 / 24; @@ -186,7 +186,7 @@ elseif (maxyear - minyear < N) sep = __calc_tick_sep__ (minmonth , maxmonth); xmin = datenum (ymin, sep * floor (minmonth / sep), 1); - xmax = datenum (ymin, sep * ceil (maxmonth / sep), 1); + xmax = datenum (ymax, sep * ceil (maxmonth / sep), 1); nticks = ceil (maxmonth / sep) - floor (minmonth / sep) + 1; else sep = __calc_tick_sep__ (minyear , maxyear);
--- a/src/DLD-FUNCTIONS/rand.cc +++ b/src/DLD-FUNCTIONS/rand.cc @@ -1,3 +1,4 @@ + /* Copyright (C) 1996-2012 John W. Eaton @@ -60,6 +61,7 @@ NDArray a; int idx = 0; dim_vector dims; + bool is_single = false; unwind_protect frame; // Restore current distribution on any exit. @@ -68,6 +70,19 @@ octave_rand::distribution (distribution); + if (nargin > 0 && args(nargin-1).is_string()) + { + std::string s_arg = args(nargin-1).string_value (); + + if (s_arg == "single") + { + is_single = true; + nargin--; + } + else if (s_arg == "double") + nargin--; + } + if (additional_arg) { if (nargin == 0) @@ -298,27 +313,54 @@ dims.chop_trailing_singletons (); - if (additional_arg) + if (is_single) { - if (a.length() == 1) - return octave_rand::nd_array (dims, a(0)); - else + if (additional_arg) { - if (a.dims() != dims) + if (a.length() == 1) + return octave_rand::float_nd_array (dims, a(0)); + else { - error ("%s: mismatch in argument size", fcn); - return retval; + if (a.dims() != dims) + { + error ("%s: mismatch in argument size", fcn); + return retval; + } + octave_idx_type len = a.length (); + FloatNDArray m (dims); + float *v = m.fortran_vec (); + for (octave_idx_type i = 0; i < len; i++) + v[i] = octave_rand::float_scalar (a(i)); + return m; } - octave_idx_type len = a.length (); - NDArray m (dims); - double *v = m.fortran_vec (); - for (octave_idx_type i = 0; i < len; i++) - v[i] = octave_rand::scalar (a(i)); - return m; } + else + return octave_rand::float_nd_array (dims); } else - return octave_rand::nd_array (dims); + { + if (additional_arg) + { + if (a.length() == 1) + return octave_rand::nd_array (dims, a(0)); + else + { + if (a.dims() != dims) + { + error ("%s: mismatch in argument size", fcn); + return retval; + } + octave_idx_type len = a.length (); + NDArray m (dims); + double *v = m.fortran_vec (); + for (octave_idx_type i = 0; i < len; i++) + v[i] = octave_rand::scalar (a(i)); + return m; + } + } + else + return octave_rand::nd_array (dims); + } } DEFUN_DLD (rand, args, , @@ -332,6 +374,8 @@ @deftypefnx {Loadable Function} {@var{v} =} rand (\"seed\")\n\ @deftypefnx {Loadable Function} {} rand (\"seed\", @var{v})\n\ @deftypefnx {Loadable Function} {} rand (\"seed\", \"reset\")\n\ +@deftypefnx {Loadable Function} {} rand (@dots{}, \"single\")\n\ +@deftypefnx {Loadable Function} {} rand (@dots{}, \"double\")\n\ Return a matrix with random elements uniformly distributed on the\n\ interval (0, 1). The arguments are handled the same as the arguments\n\ for @code{eye}.\n\ @@ -402,6 +446,9 @@ \n\ The state or seed of the generator can be reset to a new random value\n\ using the \"reset\" keyword.\n\ +\n\ +The class of the value returned can be controlled by a trailing \"double\"\n\ +or \"single\" argument. These are the only valid classes.\n\ @seealso{randn, rande, randg, randp}\n\ @end deftypefn") { @@ -499,6 +546,8 @@ @deftypefnx {Loadable Function} {@var{v} =} randn (\"seed\")\n\ @deftypefnx {Loadable Function} {} randn (\"seed\", @var{v})\n\ @deftypefnx {Loadable Function} {} randn (\"seed\", \"reset\")\n\ +@deftypefnx {Loadable Function} {} randn (@dots{}, \"single\")\n\ +@deftypefnx {Loadable Function} {} randn (@dots{}, \"double\")\n\ Return a matrix with normally distributed random\n\ elements having zero mean and variance one. The arguments are\n\ handled the same as the arguments for @code{rand}.\n\ @@ -506,6 +555,9 @@ By default, @code{randn} uses the Marsaglia and Tsang ``Ziggurat technique''\n\ to transform from a uniform to a normal distribution.\n\ \n\ +The class of the value returned can be controlled by a trailing \"double\"\n\ +or \"single\" argument. These are the only valid classes.\n\ +\n\ Reference: G. Marsaglia and W.W. Tsang,\n\ @cite{Ziggurat Method for Generating Random Variables},\n\ J. Statistical Software, vol 5, 2000,\n\ @@ -565,12 +617,17 @@ @deftypefnx {Loadable Function} {@var{v} =} rande (\"seed\")\n\ @deftypefnx {Loadable Function} {} rande (\"seed\", @var{v})\n\ @deftypefnx {Loadable Function} {} rande (\"seed\", \"reset\")\n\ +@deftypefnx {Loadable Function} {} rande (@dots{}, \"single\")\n\ +@deftypefnx {Loadable Function} {} rande (@dots{}, \"double\")\n\ Return a matrix with exponentially distributed random elements. The\n\ arguments are handled the same as the arguments for @code{rand}.\n\ \n\ By default, @code{randn} uses the Marsaglia and Tsang ``Ziggurat technique''\n\ to transform from a uniform to an exponential distribution.\n\ \n\ +The class of the value returned can be controlled by a trailing \"double\"\n\ +or \"single\" argument. These are the only valid classes.\n\ +\n\ Reference: G. Marsaglia and W.W. Tsang,\n\ @cite{Ziggurat Method for Generating Random Variables},\n\ J. Statistical Software, vol 5, 2000,\n\ @@ -632,6 +689,8 @@ @deftypefnx {Loadable Function} {@var{v} =} randg (\"seed\")\n\ @deftypefnx {Loadable Function} {} randg (\"seed\", @var{v})\n\ @deftypefnx {Loadable Function} {} randg (\"seed\", \"reset\")\n\ +@deftypefnx {Loadable Function} {} randg (@dots{}, \"single\")\n\ +@deftypefnx {Loadable Function} {} randg (@dots{}, \"double\")\n\ Return a matrix with @code{gamma(@var{a},1)} distributed random elements.\n\ The arguments are handled the same as the arguments for @code{rand},\n\ except for the argument @var{a}.\n\ @@ -711,6 +770,9 @@ @end example\n\ \n\ @end table\n\ +\n\ +The class of the value returned can be controlled by a trailing \"double\"\n\ +or \"single\" argument. These are the only valid classes.\n\ @seealso{rand, randn, rande, randp}\n\ @end deftypefn") { @@ -898,6 +960,8 @@ @deftypefnx {Loadable Function} {@var{v} =} randp (\"seed\")\n\ @deftypefnx {Loadable Function} {} randp (\"seed\", @var{v})\n\ @deftypefnx {Loadable Function} {} randp (\"seed\", \"reset\")\n\ +@deftypefnx {Loadable Function} {} randp (@dots{}, \"single\")\n\ +@deftypefnx {Loadable Function} {} randp (@dots{}, \"double\")\n\ Return a matrix with Poisson distributed random elements with mean value\n\ parameter given by the first argument, @var{l}. The arguments\n\ are handled the same as the arguments for @code{rand}, except for the\n\ @@ -928,6 +992,9 @@ L. Montanet, et al., @cite{Review of Particle Properties}, Physical Review\n\ D 50 p1284, 1994.\n\ @end table\n\ +\n\ +The class of the value returned can be controlled by a trailing \"double\"\n\ +or \"single\" argument. These are the only valid classes.\n\ @seealso{rand, randn, rande, randg}\n\ @end deftypefn") {
--- a/src/Makefile.am +++ b/src/Makefile.am @@ -100,6 +100,7 @@ defaults.in.h \ DOCSTRINGS \ find-defun-files.sh \ + gendoc.pl \ genprops.awk \ gl2ps.c \ graphics.in.h \ @@ -107,7 +108,6 @@ mk-pkg-add \ mkbuiltins \ mkdefs \ - mkgendoc \ mkoctfile.in.cc \ mkoctfile.in.sh \ mkops \ @@ -664,27 +664,16 @@ AM_CXXFLAGS := $(filter-out $(DLL_CXXDEFS), $(AM_CXXFLAGS) $(GRAPHICS_CFLAGS)) if AMCOND_BUILD_DOCS -.DOCSTRINGS: gendoc$(BUILD_EXEEXT) +.DOCSTRINGS: $(ALL_DEF_FILES) gendoc.pl if [ "x$(srcdir)" != "x." ] && [ -f $(srcdir)/DOCSTRINGS ] && [ ! -f DOCSTRINGS ]; then \ cp $(srcdir)/DOCSTRINGS DOCSTRINGS; \ touch -r $(srcdir)/DOCSTRINGS DOCSTRINGS; \ fi @echo "creating .DOCSTRINGS from .cc source files" - @./gendoc > $@ + @$(PERL) $(srcdir)/gendoc.pl $(ALL_DEF_FILES) > $@ $(top_srcdir)/build-aux/move-if-change $@ DOCSTRINGS touch $@ -doc-files: $(ALL_DEF_FILES) - echo $(ALL_DEF_FILES) > $@-t - mv $@-t $@ - -gendoc.cc: doc-files mkgendoc - $(srcdir)/mkgendoc doc-files > $@-t - mv $@-t $@ - -gendoc$(BUILD_EXEEXT): gendoc.cc - $(BUILD_CXX) $(BUILD_CXXFLAGS) -o $@ $^ $(BUILD_LDFLAGS) - all-local: $(OCT_STAMP_FILES) $(DLD_FUNCTIONS_PKG_ADD_FILE) .DOCSTRINGS else all-local: $(OCT_STAMP_FILES) $(DLD_FUNCTIONS_PKG_ADD_FILE) @@ -771,9 +760,6 @@ CLEANFILES = \ $(bin_SCRIPTS) \ $(DLD_FUNCTIONS_PKG_ADD_FILE) \ - doc-files \ - gendoc.cc \ - gendoc$(BUILD_EXEEXT) \ graphics-props.cc \ oct-parse.output
new file mode 100755 --- /dev/null +++ b/src/gendoc.pl @@ -0,0 +1,73 @@ +#! /usr/bin/perl -w +# +# Copyright (C) 2012 Rik Wehbring +# +# This file is part of Octave. +# +# Octave is free software; you can redistribute it and/or modify it +# under the terms of the GNU General Public License as published by the +# Free Software Foundation; either version 3 of the License, or (at +# your option) any later version. +# +# Octave is distributed in the hope that it will be useful, but WITHOUT +# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or +# FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License +# for more details. +# +# You should have received a copy of the GNU General Public License +# along with Octave; see the file COPYING. If not, see +# <http://www.gnu.org/licenses/>. + +unless (@ARGV > 1) { die "Usage: $0 df-file1 ..." } + +print <<__END_OF_MSG__; +### DO NOT EDIT! +### +### This file is generated automatically from Octave source files. +### Edit source files directly and run make to update this file. + +__END_OF_MSG__ + +DFFILE: foreach $df_fname (@ARGV) +{ + open (DF_FH, $df_fname) or die "Unable to open $df_fname"; + + $src_fname = ""; + @func_list = (); + @docstr = (); + + LINE: while (<DF_FH>) + { + if (/XDEFUN_FILE_NAME \("([^"]+)"/) + { + $src_fname = $1; + next LINE; + } + if (/XDEF/ and ! /XDEFALIAS/) + { + ## Decode 4 or 5 part macro definition. + ($func, $str) = /\("?(\w+)"?,[^,]+,[^,]+,(?:[^,]+,)?\s*"(.*)"\)\s*$/ ; + + unless ($func) { die "Unable to parse $df_fname at line $.\n" } + + push (@func_list, $func); + ## Do escape sequence expansion + $str =~ s/(?<!\\)\\n/\n/g; + $str =~ s/\\([^\\])/$1/g; + $str =~ s/\\\\/\\/g; + push (@docstr, $str); + } + } + close (DF_FH); + + ## Print results in DOCSTRING format + foreach $i (0 .. $#func_list) + { + $func = $func_list[$i]; + print "$func\n"; + print "\@c $func $src_fname\n"; + print $docstr[$i],"\n"; + } + +} +
deleted file mode 100755 --- a/src/mkgendoc +++ /dev/null @@ -1,124 +0,0 @@ -#! /bin/sh -# -# Copyright (C) 1999-2012 John W. Eaton -# -# This file is part of Octave. -# -# Octave is free software; you can redistribute it and/or modify it -# under the terms of the GNU General Public License as published by the -# Free Software Foundation; either version 3 of the License, or (at -# your option) any later version. -# -# Octave is distributed in the hope that it will be useful, but WITHOUT -# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or -# FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License -# for more details. -# -# You should have received a copy of the GNU General Public License -# along with Octave; see the file COPYING. If not, see -# <http://www.gnu.org/licenses/>. - -SED=${SED:-'sed'} - -if test $# -ne 1; then - echo "usage: mkgendoc f1" 1>&2 - exit 1 -fi - -DOC_FILES=`cat $1` - -if test -z "$DOC_FILES"; then - echo "mkgendoc: DOC_FILES is empty!" 1>&2 - exit 1 -fi - -cat << \EOF -// DO NOT EDIT! Generated automatically by mkgendoc - -#if defined (__DECCXX) -#define __USE_STD_IOSTREAM -#endif - -#include <iostream> -#include <string> - -#define XDEFUN_FILE_NAME(name) \ - std::string file_name = name; - -#define XDEFUN_DLD_INTERNAL(name, args_name, nargout_name, doc) \ - print_doc_string (#name, file_name, doc); - -#define XDEFUNX_DLD_INTERNAL(name, fname, args_name, nargout_name, doc) \ - print_doc_string (name, file_name, doc); - -#define XDEFUN_INTERNAL(name, args_name, nargout_name, doc) \ - print_doc_string (#name, file_name, doc); - -#define XDEFCONSTFUN_INTERNAL(name, args_name, nargout_name, doc) \ - print_doc_string (#name, file_name, doc); - -#define XDEFUNX_INTERNAL(name, fname, args_name, nargout_name, doc) \ - print_doc_string (name, file_name, doc); - -#define XDEFALIAS_INTERNAL(alias, name) - -#define XDEFVAR_INTERNAL(name, sname, defn, protect, chg_fcn, doc) \ - print_doc_string (#name, file_name, doc); - -#define XDEFCONST_INTERNAL(name, defn, doc) \ - print_doc_string (#name, file_name, doc); - -static void -print_doc_string (const std::string& name, const std::string& file_name, - const std::string& doc) -{ - std::cout << ""; - - size_t len = name.length (); - - if (name[0] == '"' && name[len-1] == '"') - std::cout << name.substr (1, len-2) << "\n"; - else - std::cout << name << "\n"; - - std::cout << "@c " << name << " " << file_name << "\n" << doc << "\n"; -} - -EOF - -for file in $DOC_FILES; do - fcn=`echo $file | $SED 's,.*/,,; s/\.df//; s/-/_/g;'` - echo "static void" - echo "print_${fcn}_doc_strings (void)" - echo "{" - cat $file - echo "}" - echo "" -done - -cat << \EOF - -int -main (void) -{ - std::cout - << "### DO NOT EDIT!\n" - << "###\n" - << "### This file is generated automatically from the Octave sources.\n" - << "### Edit those files instead and run make to update this file.\n" - << std::endl; - -EOF - -for file in $DOC_FILES; do - fcn=`echo $file | $SED 's,.*/,,; s/\.df//; s/-/_/g;'` - echo " print_${fcn}_doc_strings ();" -done - -cat << \EOF - - return 0; -} -EOF - -exit 0
--- a/src/oct-stream.cc +++ b/src/oct-stream.cc @@ -1107,9 +1107,12 @@ case 'i': { - int c1 = is.get (); - - if (! is.eof ()) + int c1 = EOF; + + while (is && (c1 = is.get ()) != EOF && isspace (c1)) + /* skip whitespace */; + + if (c1 != EOF) { if (c1 == '0') {
--- a/test/test_io.m +++ b/test/test_io.m @@ -250,6 +250,12 @@ %!assert (sscanf (['ab'; 'cd'], '%s'), 'acbd') +%!assert (sscanf ('02:08:30', '%i:%i:%i'), [2; 0]); +%!assert (sscanf ('02:08:30', '%d:%d:%d'), [2; 8; 30]); + +%!assert (sscanf ('0177 08', '%i'), [127; 0; 8]); +%!assert (sscanf ('0177 08', '%d'), [177; 8]); + %!test %! [val, count, msg, pos] = sscanf ("3I2", "%f"); %! assert (val, 3);