Mercurial > hg > octave-max
changeset 9042:97aa01a85ea4
Merge documentation cleanup changes to main branch.
author | Rik <rdrider0-list@yahoo.com> |
---|---|
date | Wed, 25 Mar 2009 18:13:08 -0700 |
parents | 12ca81f1fa99 (diff) 853f96e8008f (current diff) |
children | cfad8f9a77fa |
files | NEWS scripts/plot/print.m src/DLD-FUNCTIONS/cellfun.cc src/DLD-FUNCTIONS/find.cc src/parse.y |
diffstat | 94 files changed, 1801 insertions(+), 517 deletions(-) [+] |
line wrap: on
line diff
--- a/ChangeLog +++ b/ChangeLog @@ -1,3 +1,12 @@ +2009-03-23 Jaroslav Hajek <highegg@gmail.com> + + * NEWS: More updates. + +2009-03-20 Jaroslav Hajek <highegg@gmail.com> + + * aclocal.m4 (OCTAVE_CMATH_FUNC): New macro. + * configure.in: Use it. + 2009-03-09 John W. Eaton <jwe@octave.org> * run-octave.in: Use doc-cache instead of DOC for doc cache file.
--- a/NEWS +++ b/NEWS @@ -202,6 +202,39 @@ The performance of the sum, prod, sumsq, cumsum, cumprod, any, all, max and min functions has been significantly improved. + ** Sorting and searching. + + The performance of sort has been improved, especially when sorting + indices are requested. An efficient built-in issorted implementation + was added. sortrows now uses a more efficient algorithm, especially + in the homegeneous case. lookup is now a built-in function performing + a binary search, optimized for long runs of close elements. Lookup + also works with cell arrays of strings. + + ** Range arithmetics + + For some operations on ranges, Octave will attempt to keep the result as a + range. These include negation, adding a scalar, subtracting a scalar, and + multiplying by a scalar. Ranges with zero increment are allowed and can be + constructed using the built-in function `ones'. + + ** Various performance improvements. + + Performance of a number of other built-in operations and functions was + improved, including: + + * logical operations + * comparison operators + * element-wise power + * accumarray + * cellfun + * isnan + * isinf + * isfinite + * nchoosek + * repmat + * strcmp + ** 64-bit integer arithmetic. Arithmetic with 64-bit integers (int64 and uint64 types) is fully @@ -262,21 +295,22 @@ addtodate hypot reallog bicgstab idivide realpow - cgs info realsqrt - command_line_path interp1q rectint - contrast isdebugmode regexptranslate - convn isfloat restoredefaultpath - cummin isstrprop roundb - cummax log1p rundemos - datetick lsqnonneg runlength - display matlabroot saveobj - expm1 namelengthmax spaugment - filemarker nargoutchk strchr - fstat pathdef strvcat - full perl subspace - fzero prctile symvar - genvarname quantile treelayout - histc re_read_readline_init_file validatestring + cellslices info realsqrt + cgs interp1q rectint + command_line_path isdebugmode regexptranslate + contrast isfloat restoredefaultpath + convn isstrprop roundb + cummin log1p rundemos + cummax lsqnonneg runlength + datetick matlabroot saveobj + display namelengthmax spaugment + expm1 nargoutchk strchr + filemarker pathdef strvcat + fstat perl subspace + full prctile symvar + fzero quantile treelayout + genvarname re_read_readline_init_file validatestring + histc ** Changes to strcat.
--- a/aclocal.m4 +++ b/aclocal.m4 @@ -1285,6 +1285,41 @@ AC_SUBST([FT2_LIBS])]) dnl end of freetype2.m4 +dnl Check whether a math mapper function is available in <cmath>. +dnl Will define HAVE_CMATH_FUNC if there is a double variant and +dnl HAVE_CMATH_FUNCF if there is a float variant. +dnl Currently capable of checking for functions with single +dnl argument and returning bool/int/real. +AC_DEFUN([OCTAVE_CMATH_FUNC],[ +AC_MSG_CHECKING([for std::$1 in <cmath>]) +AC_LANG_PUSH(C++) +AC_COMPILE_IFELSE([AC_LANG_PROGRAM([[ +#include <cmath> +void take_func (bool (*func) (double x)); +void take_func (int (*func) (double x)); +void take_func (double (*func) (double x)); +]], +[[ +take_func(std::$1); +]])], +[AC_MSG_RESULT([yes]) + AC_DEFINE(HAVE_CMATH_[]AS_TR_CPP($1),1,[Define if <cmath> provides $1])], +[AC_MSG_RESULT([no])]) +AC_MSG_CHECKING([for std::$1 (float variant) in <cmath>]) +AC_COMPILE_IFELSE([AC_LANG_PROGRAM([[ +#include <cmath> +void take_func (bool (*func) (float x)); +void take_func (int (*func) (float x)); +void take_func (float (*func) (float x)); +]], +[[ +take_func(std::$1); +]])], +[AC_MSG_RESULT([yes]) + AC_DEFINE(HAVE_CMATH_[]AS_TR_CPP($1)F,1,[Define if <cmath> provides float variant of $1])], +[AC_MSG_RESULT([no])]) +]) + dnl Check whether fast signed integer arithmetics using bit tricks dnl can be used in oct-inttypes.h. Defines HAVE_FAST_INT_OPS if dnl the following conditions hold:
--- a/configure.in +++ b/configure.in @@ -1776,6 +1776,12 @@ [Define if your system has a single-arg prototype for gettimeofday.])]) fi +dnl Maybe <cmath> defines the IEEE functions we need. + +OCTAVE_CMATH_FUNC(isnan) +OCTAVE_CMATH_FUNC(isinf) +OCTAVE_CMATH_FUNC(isfinite) + dnl Would like to get rid of this cruft, and just have dnl dnl AC_CHECK_FUNCS(finite isnan isinf)
--- a/doc/ChangeLog +++ b/doc/ChangeLog @@ -1,3 +1,8 @@ +2009-03-25 John W. Eaton <jwe@octave.org> + + * interpreter/munge-texi.cc (process_texi_input_file): + Copy leading comment with file name info to output. + 2009-03-09 John W. Eaton <jwe@octave.org> * interpreter/Makefile.in (DISTFILES): Use doc-cache instead of
--- a/doc/interpreter/munge-texi.cc +++ b/doc/interpreter/munge-texi.cc @@ -245,7 +245,23 @@ { std::string doc_string = help_text[symbol_name]; + size_t len = doc_string.length (); + int j = 0; + + // If there is a leading comment with the file + // name, copy it to the output. + if (len > 1 + && doc_string[j] == '@' + && doc_string[j+1] == 'c') + { + j = 2; + while (doc_string[j++] != '\n') + /* find eol */; + + os << doc_string.substr (0, j); + } + while (doc_string[j] == ' ') j++;
--- a/liboctave/Array-d.cc +++ b/liboctave/Array-d.cc @@ -27,7 +27,7 @@ // Instantiate Arrays of double values. -#include "lo-ieee.h" +#include "lo-mappers.h" #include "Array.h" #include "Array.cc" #include "oct-locbuf.h" @@ -40,19 +40,19 @@ inline bool sort_isnan<double> (double x) { - return lo_ieee_isnan (x); + return xisnan (x); } static bool nan_ascending_compare (double x, double y) { - return lo_ieee_isnan (y) ? ! lo_ieee_isnan (x) : x < y; + return xisnan (y) ? ! xisnan (x) : x < y; } static bool nan_descending_compare (double x, double y) { - return lo_ieee_isnan (x) ? ! lo_ieee_isnan (y) : x > y; + return xisnan (x) ? ! xisnan (y) : x > y; } Array<double>::compare_fcn_type
--- a/liboctave/Array-f.cc +++ b/liboctave/Array-f.cc @@ -27,7 +27,7 @@ // Instantiate Arrays of float values. -#include "lo-ieee.h" +#include "lo-mappers.h" #include "Array.h" #include "Array.cc" #include "oct-locbuf.h" @@ -40,19 +40,19 @@ inline bool sort_isnan<float> (float x) { - return lo_ieee_isnan (x); + return xisnan (x); } static bool nan_ascending_compare (float x, float y) { - return lo_ieee_isnan (y) ? ! lo_ieee_isnan (x) : x < y; + return xisnan (y) ? ! xisnan (x) : x < y; } static bool nan_descending_compare (float x, float y) { - return lo_ieee_isnan (x) ? ! lo_ieee_isnan (y) : x > y; + return xisnan (x) ? ! xisnan (y) : x > y; } Array<float>::compare_fcn_type
--- a/liboctave/Array.cc +++ b/liboctave/Array.cc @@ -461,6 +461,7 @@ class rec_permute_helper { octave_idx_type *dim, *stride; + bool use_blk; int top; public: @@ -472,37 +473,78 @@ dim = new octave_idx_type [2*n]; // A hack to avoid double allocation stride = dim + n; - top = 0; // Get cumulative dimensions. OCTAVE_LOCAL_BUFFER (octave_idx_type, cdim, n+1); cdim[0] = 1; for (int i = 1; i < n+1; i++) cdim[i] = cdim[i-1] * dv(i-1); - int k = 0; - // Reduce leading identity - for (; k < n && perm(k) == k; k++) ; - if (k > 0) + // Setup the permuted strides. + for (int k = 0; k < n; k++) { - dim[0] = cdim[k]; - stride[0] = 1; + int kk = perm(k); + dim[k] = dv(kk); + stride[k] = cdim[kk]; } - else - top = -1; - - // Setup the permuted strides. - for (; k < n; k++) + + // Reduce contiguous runs. + top = 0; + for (int k = 1; k < n; k++) { - ++top; - int kk = perm(k); - dim[top] = dv(kk); - stride[top] = cdim[kk]; + if (stride[k] == stride[top]*dim[top]) + dim[top] *= dim[k]; + else + { + top++; + dim[top] = dim[k]; + stride[top] = stride[k]; + } } + // Determine whether we can use block transposes. + use_blk = top >= 1 && stride[1] == 1 && stride[0] == dim[1]; + } ~rec_permute_helper (void) { delete [] dim; } + // Helper method for fast blocked transpose. + template <class T> + T *blk_trans (const T *src, T *dest, octave_idx_type nr, octave_idx_type nc) const + { + static const octave_idx_type m = 8; + OCTAVE_LOCAL_BUFFER (T, blk, m*m); + for (octave_idx_type kr = 0; kr < nr; kr += m) + for (octave_idx_type kc = 0; kc < nc; kc += m) + { + octave_idx_type lr = std::min (m, nr - kr); + octave_idx_type lc = std::min (m, nc - kc); + if (lr == m && lc == m) + { + const T *ss = src + kc * nr + kr; + for (octave_idx_type j = 0; j < m; j++) + for (octave_idx_type i = 0; i < m; i++) + blk[j*m+i] = ss[j*nr + i]; + T *dd = dest + kr * nc + kc; + for (octave_idx_type j = 0; j < m; j++) + for (octave_idx_type i = 0; i < m; i++) + dd[j*nc+i] = blk[i*m+j]; + } + else + { + const T *ss = src + kc * nr + kr; + for (octave_idx_type j = 0; j < lc; j++) + for (octave_idx_type i = 0; i < lr; i++) + blk[j*m+i] = ss[j*nr + i]; + T *dd = dest + kr * nc + kc; + for (octave_idx_type j = 0; j < lr; j++) + for (octave_idx_type i = 0; i < lc; i++) + dd[j*nc+i] = blk[i*m+j]; + } + } + + return dest + nr*nc; + } private: // Recursive N-d generalized transpose @@ -521,8 +563,9 @@ dest += len; } - } + else if (use_blk && lev == 1) + dest = blk_trans (src, dest, dim[1], dim[0]); else { octave_idx_type step = stride[lev], len = dim[lev];
--- a/liboctave/Array.h +++ b/liboctave/Array.h @@ -601,6 +601,26 @@ return result; } + + // This is non-breakable map, suitable for fast functions. Efficiency + // relies on compiler's ability to inline a function pointer. This seems + // to be OK with recent GCC. + template <class U> + Array<U> + fastmap (U (*fcn) (typename ref_param<T>::type)) const + { + octave_idx_type len = length (); + + const T *m = data (); + + Array<U> result (dims ()); + U *p = result.fortran_vec (); + + std::transform (m, m + len, p, fcn); + + return result; + } + }; #endif
--- a/liboctave/CColVector.cc +++ b/liboctave/CColVector.cc @@ -524,9 +524,7 @@ { octave_idx_type len = a.length(); - if (len < 1) - is.clear (std::ios::badbit); - else + if (len > 0) { double tmp; for (octave_idx_type i = 0; i < len; i++)
--- a/liboctave/CDiagMatrix.cc +++ b/liboctave/CDiagMatrix.cc @@ -30,6 +30,7 @@ #include "Array-util.h" #include "lo-error.h" +#include "lo-ieee.h" #include "mx-base.h" #include "mx-inlines.cc" #include "oct-cmplx.h"
--- a/liboctave/CMatrix.cc +++ b/liboctave/CMatrix.cc @@ -3640,9 +3640,7 @@ octave_idx_type nr = a.rows (); octave_idx_type nc = a.cols (); - if (nr < 1 || nc < 1) - is.clear (std::ios::badbit); - else + if (nr > 0 && nc > 0) { Complex tmp; for (octave_idx_type i = 0; i < nr; i++)
--- a/liboctave/CNDArray.cc +++ b/liboctave/CNDArray.cc @@ -754,6 +754,24 @@ dims ()); } +boolNDArray +ComplexNDArray::isnan (void) const +{ + return ArrayN<bool> (fastmap<bool> (xisnan)); +} + +boolNDArray +ComplexNDArray::isinf (void) const +{ + return ArrayN<bool> (fastmap<bool> (xisinf)); +} + +boolNDArray +ComplexNDArray::isfinite (void) const +{ + return ArrayN<bool> (fastmap<bool> (xfinite)); +} + ComplexNDArray conj (const ComplexNDArray& a) { @@ -898,9 +916,7 @@ { octave_idx_type nel = a.nelem (); - if (nel < 1 ) - is.clear (std::ios::badbit); - else + if (nel > 0) { Complex tmp; for (octave_idx_type i = 0; i < nel; i++)
--- a/liboctave/CNDArray.h +++ b/liboctave/CNDArray.h @@ -98,6 +98,9 @@ ComplexNDArray& insert (const ComplexNDArray& a, const Array<octave_idx_type>& ra_idx); NDArray abs (void) const; + boolNDArray isnan (void) const; + boolNDArray isinf (void) const; + boolNDArray isfinite (void) const; friend ComplexNDArray conj (const ComplexNDArray& a);
--- a/liboctave/CRowVector.cc +++ b/liboctave/CRowVector.cc @@ -434,9 +434,7 @@ { octave_idx_type len = a.length(); - if (len < 1) - is.clear (std::ios::badbit); - else + if (len > 0) { Complex tmp; for (octave_idx_type i = 0; i < len; i++)
--- a/liboctave/CSparse.cc +++ b/liboctave/CSparse.cc @@ -646,28 +646,28 @@ octave_idx_type nz = nnz (); SparseComplexMatrix retval (nc, nr, nz); - OCTAVE_LOCAL_BUFFER (octave_idx_type, w, nr + 1); - for (octave_idx_type i = 0; i < nr; i++) - w[i] = 0; for (octave_idx_type i = 0; i < nz; i++) - w[ridx(i)]++; + retval.xcidx (ridx (i) + 1)++; + // retval.xcidx[1:nr] holds the row degrees for rows 0:(nr-1) nz = 0; - for (octave_idx_type i = 0; i < nr; i++) + for (octave_idx_type i = 1; i <= nr; i++) { - retval.xcidx(i) = nz; - nz += w[i]; - w[i] = retval.xcidx(i); + const octave_idx_type tmp = retval.xcidx (i); + retval.xcidx (i) = nz; + nz += tmp; } - retval.xcidx(nr) = nz; - w[nr] = nz; + // retval.xcidx[1:nr] holds row entry *start* offsets for rows 0:(nr-1) for (octave_idx_type j = 0; j < nc; j++) for (octave_idx_type k = cidx(j); k < cidx(j+1); k++) { - octave_idx_type q = w [ridx(k)]++; + octave_idx_type q = retval.xcidx (ridx (k) + 1)++; retval.xridx (q) = j; retval.xdata (q) = conj (data (k)); } + assert (nnz () == retval.xcidx (nr)); + // retval.xcidx[1:nr] holds row entry *end* offsets for rows 0:(nr-1) + // and retval.xcidx[0:(nr-1)] holds their row entry *start* offsets return retval; } @@ -7470,9 +7470,7 @@ octave_idx_type nc = a.cols (); octave_idx_type nz = a.nzmax (); - if (nr < 1 || nc < 1) - is.clear (std::ios::badbit); - else + if (nr > 0 && nc > 0) { octave_idx_type itmp, jtmp, jold = 0; Complex tmp;
--- a/liboctave/ChangeLog +++ b/liboctave/ChangeLog @@ -1,3 +1,73 @@ +2009-03-25 John W. Eaton <jwe@octave.org> + + * Makefile.in (MATRIX_INC): Add Sparse-diag-op-defs.h and + Sparse-perm-op-defs.h to the list. + +2009-03-25 Jaroslav Hajek <highegg@gmail.com> + + * oct-inttypes.cc (INT_DOUBLE_BINOP_DECL (*, uint64), + INT_DOUBLE_BINOP_DECL (*, int64)): x -> y where appropriate. + +2009-03-25 Jaroslav Hajek <highegg@gmail.com> + + * Array.cc (rec_permute_helper::use_blk): New field. + (rec_permute_helper::blk_trans): New method. + (rec_permute_helper::rec_permute_helper): Use smart reductions, + detect possibility of using blocked transpose. + (rec_permute_helper::do_permute): Use blocked transpose if possible. + +2009-03-23 Jaroslav Hajek <highegg@gmail.com> + + * idx-vector.cc (convert_index(double,...)): Simplify. + +2009-03-21 Jaroslav Hajek <highegg@gmail.com> + + * Array-d.cc: lo_ieee_isnan -> xisnan. + * Array-f.cc: Ditto. + * oct-inttypes.cc: Ditto. + * oct-inttypes.h: Ditto. + * CDiagMatrix.cc: Add missing include. + * fCDiagMatrix.cc: Ditto. + +2009-03-20 Jaroslav Hajek <highegg@gmail.com> + + * CColVector.cc, CMatrix.cc, CNDArray.cc, CRowVector.cc, CSparse.cc, + boolSparse.cc, dColVector.cc, dMatrix.cc, dNDArray.cc, dRowVector.cc, + dSparse.cc, fCColVector.cc, fCMatrix.cc, fCNDArray.cc, fCRowVector.cc, + fColVector.cc, fMatrix.cc, fNDArray.cc, fRowVector.cc, intNDArray.cc: + Allow empty arrays in stream input operators. + +2009-03-20 Jaroslav Hajek <highegg@gmail.com> + + * Array.h (Array<T>::fastmap): New method. + * dNDArray.cc (NDArray::isnan, NDArray::isinf, NDArray::isfinite): + New methods. + * dNDArray.h: Declare them. + * fNDArray.cc (FloatNDArray::isnan, FloatNDArray::isinf, + FloatNDArray::isfinite): New methods. + * fNDArray.h: Declare them. + * CNDArray.cc (ComplexNDArray::isnan, ComplexNDArray::isinf, + ComplexNDArray::isfinite): New methods. + * CNDArray.h: Declare them. + * fCNDArray.cc (FloatComplexNDArray::isnan, FloatComplexNDArray::isinf, + FloatComplexNDArray::isfinite): New methods. + * fCNDArray.h: Declare them. + * lo-mappers.h (xisnan, xisinf, xfinite): If possible, use definitions + from <cmath>. + +2009-03-18 Jaroslav Hajek <highegg@gmail.com> + + * oct-norm.cc (get_eps): Remove that hack. + (higham): Use std::numeric_limits instead. + Include OCTAVE_QUIT. + +2009-03-16 Jason Riedy <jason@acm.org> + + * Sparse.cc (transpose): Eliminate the workspace by computing in + retval.xcidx. + * CSparse.cc (hermitian): Eliminate the workspace by computing in + retval.xcidx. + 2009-03-14 Jaroslav Hajek <highegg@gmail.com> * mx-op-decl.h (NDS_BOOL_OP_DECLS, SND_BOOL_OP_DECLS, NDND_BOOL_OP_DECLS): Support compound binary ops.
--- a/liboctave/Makefile.in +++ b/liboctave/Makefile.in @@ -53,9 +53,9 @@ dbleGEPBAL.h dbleHESS.h dbleLU.h dbleQR.h dbleQRP.h dbleSCHUR.h \ dbleSVD.h boolSparse.h CSparse.h dSparse.h MSparse-defs.h MSparse.h \ Sparse.h sparse-base-lu.h SparseCmplxLU.h SparsedbleLU.h \ - sparse-base-chol.h SparseCmplxCHOL.h \ - SparsedbleCHOL.h SparseCmplxQR.h SparseQR.h Sparse-op-defs.h \ - MatrixType.h PermMatrix.h \ + sparse-base-chol.h SparseCmplxCHOL.h SparsedbleCHOL.h \ + SparseCmplxQR.h SparseQR.h Sparse-op-defs.h Sparse-diag-op-defs.h \ + Sparse-perm-op-defs.h MatrixType.h PermMatrix.h \ int8NDArray.h uint8NDArray.h int16NDArray.h uint16NDArray.h \ int32NDArray.h uint32NDArray.h int64NDArray.h uint64NDArray.h \ intNDArray.h \
--- a/liboctave/Sparse.cc +++ b/liboctave/Sparse.cc @@ -1062,28 +1062,28 @@ octave_idx_type nz = nnz (); Sparse<T> retval (nc, nr, nz); - OCTAVE_LOCAL_BUFFER (octave_idx_type, w, nr + 1); - for (octave_idx_type i = 0; i < nr; i++) - w[i] = 0; for (octave_idx_type i = 0; i < nz; i++) - w[ridx(i)]++; + retval.xcidx (ridx (i) + 1)++; + // retval.xcidx[1:nr] holds the row degrees for rows 0:(nr-1) nz = 0; - for (octave_idx_type i = 0; i < nr; i++) + for (octave_idx_type i = 1; i <= nr; i++) { - retval.xcidx(i) = nz; - nz += w[i]; - w[i] = retval.xcidx(i); + const octave_idx_type tmp = retval.xcidx (i); + retval.xcidx (i) = nz; + nz += tmp; } - retval.xcidx(nr) = nz; - w[nr] = nz; + // retval.xcidx[1:nr] holds row entry *start* offsets for rows 0:(nr-1) for (octave_idx_type j = 0; j < nc; j++) for (octave_idx_type k = cidx(j); k < cidx(j+1); k++) { - octave_idx_type q = w [ridx(k)]++; + octave_idx_type q = retval.xcidx (ridx (k) + 1)++; retval.xridx (q) = j; retval.xdata (q) = data (k); } + assert (nnz () == retval.xcidx (nr)); + // retval.xcidx[1:nr] holds row entry *end* offsets for rows 0:(nr-1) + // and retval.xcidx[0:(nr-1)] holds their row entry *start* offsets return retval; }
--- a/liboctave/boolSparse.cc +++ b/liboctave/boolSparse.cc @@ -184,9 +184,7 @@ octave_idx_type nc = a.cols (); octave_idx_type nz = a.nzmax (); - if (nr < 1 || nc < 1) - is.clear (std::ios::badbit); - else + if (nr > 0 && nc > 0) { octave_idx_type itmp, jtmp, jold = 0; bool tmp;
--- a/liboctave/dColVector.cc +++ b/liboctave/dColVector.cc @@ -321,9 +321,7 @@ { octave_idx_type len = a.length(); - if (len < 1) - is.clear (std::ios::badbit); - else + if (len > 0) { double tmp; for (octave_idx_type i = 0; i < len; i++)
--- a/liboctave/dMatrix.cc +++ b/liboctave/dMatrix.cc @@ -3058,9 +3058,7 @@ octave_idx_type nr = a.rows (); octave_idx_type nc = a.cols (); - if (nr < 1 || nc < 1) - is.clear (std::ios::badbit); - else + if (nr > 0 && nc > 0) { double tmp; for (octave_idx_type i = 0; i < nr; i++)
--- a/liboctave/dNDArray.cc +++ b/liboctave/dNDArray.cc @@ -874,6 +874,24 @@ dims ()); } +boolNDArray +NDArray::isnan (void) const +{ + return ArrayN<bool> (fastmap<bool> (xisnan)); +} + +boolNDArray +NDArray::isinf (void) const +{ + return ArrayN<bool> (fastmap<bool> (xisinf)); +} + +boolNDArray +NDArray::isfinite (void) const +{ + return ArrayN<bool> (fastmap<bool> (xfinite)); +} + Matrix NDArray::matrix_value (void) const { @@ -947,9 +965,7 @@ { octave_idx_type nel = a.nelem (); - if (nel < 1 ) - is.clear (std::ios::badbit); - else + if (nel > 0) { double tmp; for (octave_idx_type i = 0; i < nel; i++)
--- a/liboctave/dNDArray.h +++ b/liboctave/dNDArray.h @@ -109,6 +109,9 @@ NDArray& insert (const NDArray& a, const Array<octave_idx_type>& ra_idx); NDArray abs (void) const; + boolNDArray isnan (void) const; + boolNDArray isinf (void) const; + boolNDArray isfinite (void) const; ComplexNDArray fourier (int dim = 1) const; ComplexNDArray ifourier (int dim = 1) const;
--- a/liboctave/dRowVector.cc +++ b/liboctave/dRowVector.cc @@ -291,9 +291,7 @@ { octave_idx_type len = a.length(); - if (len < 1) - is.clear (std::ios::badbit); - else + if (len > 0) { double tmp; for (octave_idx_type i = 0; i < len; i++)
--- a/liboctave/dSparse.cc +++ b/liboctave/dSparse.cc @@ -7588,9 +7588,7 @@ octave_idx_type nc = a.cols (); octave_idx_type nz = a.nzmax (); - if (nr < 1 || nc < 1) - is.clear (std::ios::badbit); - else + if (nr > 0 && nc > 0) { octave_idx_type itmp, jtmp, jold = 0; double tmp;
--- a/liboctave/fCColVector.cc +++ b/liboctave/fCColVector.cc @@ -524,9 +524,7 @@ { octave_idx_type len = a.length(); - if (len < 1) - is.clear (std::ios::badbit); - else + if (len > 0) { float tmp; for (octave_idx_type i = 0; i < len; i++)
--- a/liboctave/fCDiagMatrix.cc +++ b/liboctave/fCDiagMatrix.cc @@ -30,6 +30,7 @@ #include "Array-util.h" #include "lo-error.h" +#include "lo-ieee.h" #include "mx-base.h" #include "mx-inlines.cc" #include "oct-cmplx.h"
--- a/liboctave/fCMatrix.cc +++ b/liboctave/fCMatrix.cc @@ -3633,9 +3633,7 @@ octave_idx_type nr = a.rows (); octave_idx_type nc = a.cols (); - if (nr < 1 || nc < 1) - is.clear (std::ios::badbit); - else + if (nr > 0 && nc > 0) { FloatComplex tmp; for (octave_idx_type i = 0; i < nr; i++)
--- a/liboctave/fCNDArray.cc +++ b/liboctave/fCNDArray.cc @@ -749,6 +749,24 @@ dims ()); } +boolNDArray +FloatComplexNDArray::isnan (void) const +{ + return ArrayN<bool> (fastmap<bool> (xisnan)); +} + +boolNDArray +FloatComplexNDArray::isinf (void) const +{ + return ArrayN<bool> (fastmap<bool> (xisinf)); +} + +boolNDArray +FloatComplexNDArray::isfinite (void) const +{ + return ArrayN<bool> (fastmap<bool> (xfinite)); +} + FloatComplexNDArray conj (const FloatComplexNDArray& a) { @@ -893,9 +911,7 @@ { octave_idx_type nel = a.nelem (); - if (nel < 1 ) - is.clear (std::ios::badbit); - else + if (nel > 0) { FloatComplex tmp; for (octave_idx_type i = 0; i < nel; i++)
--- a/liboctave/fCNDArray.h +++ b/liboctave/fCNDArray.h @@ -98,6 +98,9 @@ FloatComplexNDArray& insert (const FloatComplexNDArray& a, const Array<octave_idx_type>& ra_idx); FloatNDArray abs (void) const; + boolNDArray isnan (void) const; + boolNDArray isinf (void) const; + boolNDArray isfinite (void) const; friend FloatComplexNDArray conj (const FloatComplexNDArray& a);
--- a/liboctave/fCRowVector.cc +++ b/liboctave/fCRowVector.cc @@ -434,9 +434,7 @@ { octave_idx_type len = a.length(); - if (len < 1) - is.clear (std::ios::badbit); - else + if (len > 0) { FloatComplex tmp; for (octave_idx_type i = 0; i < len; i++)
--- a/liboctave/fColVector.cc +++ b/liboctave/fColVector.cc @@ -321,9 +321,7 @@ { octave_idx_type len = a.length(); - if (len < 1) - is.clear (std::ios::badbit); - else + if (len > 0) { float tmp; for (octave_idx_type i = 0; i < len; i++)
--- a/liboctave/fMatrix.cc +++ b/liboctave/fMatrix.cc @@ -3057,9 +3057,7 @@ octave_idx_type nr = a.rows (); octave_idx_type nc = a.cols (); - if (nr < 1 || nc < 1) - is.clear (std::ios::badbit); - else + if (nr > 0 && nc > 0) { float tmp; for (octave_idx_type i = 0; i < nr; i++)
--- a/liboctave/fNDArray.cc +++ b/liboctave/fNDArray.cc @@ -829,6 +829,24 @@ dims ()); } +boolNDArray +FloatNDArray::isnan (void) const +{ + return ArrayN<bool> (fastmap<bool> (xisnan)); +} + +boolNDArray +FloatNDArray::isinf (void) const +{ + return ArrayN<bool> (fastmap<bool> (xisinf)); +} + +boolNDArray +FloatNDArray::isfinite (void) const +{ + return ArrayN<bool> (fastmap<bool> (xfinite)); +} + FloatMatrix FloatNDArray::matrix_value (void) const { @@ -902,9 +920,7 @@ { octave_idx_type nel = a.nelem (); - if (nel < 1 ) - is.clear (std::ios::badbit); - else + if (nel > 0) { float tmp; for (octave_idx_type i = 0; i < nel; i++)
--- a/liboctave/fNDArray.h +++ b/liboctave/fNDArray.h @@ -106,6 +106,9 @@ FloatNDArray& insert (const FloatNDArray& a, const Array<octave_idx_type>& ra_idx); FloatNDArray abs (void) const; + boolNDArray isnan (void) const; + boolNDArray isinf (void) const; + boolNDArray isfinite (void) const; FloatComplexNDArray fourier (int dim = 1) const; FloatComplexNDArray ifourier (int dim = 1) const;
--- a/liboctave/fRowVector.cc +++ b/liboctave/fRowVector.cc @@ -291,9 +291,7 @@ { octave_idx_type len = a.length(); - if (len < 1) - is.clear (std::ios::badbit); - else + if (len > 0) { float tmp; for (octave_idx_type i = 0; i < len; i++)
--- a/liboctave/idx-vector.cc +++ b/liboctave/idx-vector.cc @@ -182,18 +182,9 @@ inline octave_idx_type convert_index (double x, bool& conv_error, octave_idx_type& ext) { - octave_idx_type i; - if (xisnan (x) || xisinf (x)) - { - i = 0; - conv_error = true; - } - else - { - i = static_cast<octave_idx_type> (x); - if (static_cast<double> (i) != x) - conv_error = true; - } + octave_idx_type i = static_cast<octave_idx_type> (x); + if (static_cast<double> (i) != x) + conv_error = true; return convert_index (i, conv_error, ext); }
--- a/liboctave/intNDArray.cc +++ b/liboctave/intNDArray.cc @@ -146,9 +146,7 @@ { octave_idx_type nel = a.nelem (); - if (nel < 1 ) - is.clear (std::ios::badbit); - else + if (nel > 0) { T tmp;
--- a/liboctave/lo-mappers.cc +++ b/liboctave/lo-mappers.cc @@ -190,23 +190,29 @@ // double -> bool mappers. +#if ! defined(HAVE_CMATH_ISNAN) bool xisnan (double x) { return lo_ieee_isnan (x); } +#endif +#if ! defined(HAVE_CMATH_ISFINITE) bool xfinite (double x) { return lo_ieee_finite (x); } +#endif +#if ! defined(HAVE_CMATH_ISINF) bool xisinf (double x) { return lo_ieee_isinf (x); } +#endif bool octave_is_NA (double x) @@ -321,28 +327,6 @@ // complex -> bool mappers. bool -xisnan (const Complex& x) -{ - return (xisnan (real (x)) || xisnan (imag (x))); -} - -bool -xfinite (const Complex& x) -{ - double rx = real (x); - double ix = imag (x); - - return (xfinite (rx) && ! xisnan (rx) - && xfinite (ix) && ! xisnan (ix)); -} - -bool -xisinf (const Complex& x) -{ - return (xisinf (real (x)) || xisinf (imag (x))); -} - -bool octave_is_NA (const Complex& x) { return (octave_is_NA (real (x)) || octave_is_NA (imag (x))); @@ -524,23 +508,29 @@ // float -> bool mappers. +#if ! defined(HAVE_CMATH_ISNANF) bool xisnan (float x) { return lo_ieee_isnan (x); } +#endif +#if ! defined(HAVE_CMATH_ISFINITEF) bool xfinite (float x) { return lo_ieee_finite (x); } +#endif +#if ! defined(HAVE_CMATH_ISINFF) bool xisinf (float x) { return lo_ieee_isinf (x); } +#endif bool octave_is_NA (float x) @@ -655,28 +645,6 @@ // complex -> bool mappers. bool -xisnan (const FloatComplex& x) -{ - return (xisnan (real (x)) || xisnan (imag (x))); -} - -bool -xfinite (const FloatComplex& x) -{ - float rx = real (x); - float ix = imag (x); - - return (xfinite (rx) && ! xisnan (rx) - && xfinite (ix) && ! xisnan (ix)); -} - -bool -xisinf (const FloatComplex& x) -{ - return (xisinf (real (x)) || xisinf (imag (x))); -} - -bool octave_is_NA (const FloatComplex& x) { return (octave_is_NA (real (x)) || octave_is_NA (imag (x)));
--- a/liboctave/lo-mappers.h +++ b/liboctave/lo-mappers.h @@ -25,6 +25,7 @@ #define octave_liboctave_mappers_h 1 #include "oct-cmplx.h" +#include "lo-math.h" // Double Precision extern OCTAVE_API double arg (double x); @@ -46,9 +47,24 @@ inline bool xisnan (bool) { return false; } inline bool xisnan (char) { return false; } +#if defined (HAVE_CMATH_ISNAN) +inline bool xisnan (double x) +{ return std::isnan (x); } +#else extern OCTAVE_API bool xisnan (double x); +#endif +#if defined (HAVE_CMATH_ISFINITE) +inline bool xfinite (double x) +{ return std::isfinite (x); } +#else extern OCTAVE_API bool xfinite (double x); +#endif +#if defined (HAVE_CMATH_ISINF) +inline bool xisinf (double x) +{ return std::isinf (x); } +#else extern OCTAVE_API bool xisinf (double x); +#endif extern OCTAVE_API bool octave_is_NA (double x); extern OCTAVE_API bool octave_is_NaN_or_NA (double x) GCC_ATTR_DEPRECATED; @@ -70,9 +86,15 @@ extern OCTAVE_API Complex xroundb (const Complex& x); extern OCTAVE_API Complex signum (const Complex& x); -extern OCTAVE_API bool xisnan (const Complex& x); -extern OCTAVE_API bool xfinite (const Complex& x); -extern OCTAVE_API bool xisinf (const Complex& x); +inline bool +xisnan (const Complex& x) +{ return (xisnan (real (x)) || xisnan (imag (x))); } +inline bool +xfinite (const Complex& x) +{ return (xfinite (real (x)) && xfinite (imag (x))); } +inline bool +xisinf (const Complex& x) +{ return (xisinf (real (x)) || xisinf (imag (x))); } extern OCTAVE_API bool octave_is_NA (const Complex& x); extern OCTAVE_API bool octave_is_NaN_or_NA (const Complex& x); @@ -96,9 +118,25 @@ extern OCTAVE_API FloatComplex xlog2 (const FloatComplex& x, int& exp); extern OCTAVE_API float xexp2 (float x); +#if defined (HAVE_CMATH_ISNANF) +inline bool xisnan (float x) +{ return std::isnan (x); } +#else extern OCTAVE_API bool xisnan (float x); +#endif +#if defined (HAVE_CMATH_ISFINITEF) +inline bool xfinite (float x) +{ return std::isfinite (x); } +#else extern OCTAVE_API bool xfinite (float x); +#endif +#if defined (HAVE_CMATH_ISINFF) +inline bool xisinf (float x) +{ return std::isinf (x); } +#else extern OCTAVE_API bool xisinf (float x); +#endif + extern OCTAVE_API bool octave_is_NA (float x); extern OCTAVE_API bool octave_is_NaN_or_NA (float x) GCC_ATTR_DEPRECATED; @@ -120,9 +158,15 @@ extern OCTAVE_API FloatComplex xroundb (const FloatComplex& x); extern OCTAVE_API FloatComplex signum (const FloatComplex& x); -extern OCTAVE_API bool xisnan (const FloatComplex& x); -extern OCTAVE_API bool xfinite (const FloatComplex& x); -extern OCTAVE_API bool xisinf (const FloatComplex& x); +inline bool +xisnan (const FloatComplex& x) +{ return (xisnan (real (x)) || xisnan (imag (x))); } +inline bool +xfinite (const FloatComplex& x) +{ return (xfinite (real (x)) && xfinite (imag (x))); } +inline bool +xisinf (const FloatComplex& x) +{ return (xisinf (real (x)) || xisinf (imag (x))); } extern OCTAVE_API bool octave_is_NA (const FloatComplex& x); extern OCTAVE_API bool octave_is_NaN_or_NA (const FloatComplex& x);
--- a/liboctave/mx-op-defs.h +++ b/liboctave/mx-op-defs.h @@ -1120,6 +1120,13 @@ #define MPM_BIN_OPS(R, M, PM) \ MPM_MULTIPLY_OP(M, PM); +#define NDND_MAPPER_BODY(R, NAME) \ + R retval (dims ()); \ + octave_idx_type n = numel (); \ + for (octave_idx_type i = 0; i < n; i++) \ + retval.xelem (i) = NAME (elem (i)); \ + return retval; + #endif
--- a/liboctave/oct-inttypes.cc +++ b/liboctave/oct-inttypes.cc @@ -407,7 +407,7 @@ { return x / octave_uint64 (static_cast<uint64_t> (2)); } - else if (y < 0 || lo_ieee_isnan (x) || lo_ieee_isinf (x)) + else if (y < 0 || xisnan (y) || xisinf (y)) { return octave_uint64 (x.value () * y); } @@ -442,7 +442,7 @@ { return x / octave_int64 (static_cast<uint64_t> (4*y)); } - else if (lo_ieee_isnan (x) || lo_ieee_isinf (x)) + else if (xisnan (y) || xisinf (y)) { return octave_int64 (x.value () * y); }
--- a/liboctave/oct-inttypes.h +++ b/liboctave/oct-inttypes.h @@ -33,11 +33,11 @@ #include "lo-traits.h" #include "lo-math.h" #include "oct-types.h" -#include "lo-ieee.h" #include "lo-mappers.h" #ifdef OCTAVE_INT_USE_LONG_DOUBLE inline long double xround (long double x) { return roundl (x); } +inline long double xisnan (long double x) { return xisnan (static_cast<double> (x)); } #endif // Undefine min/max if needed (this may happen under Windows) @@ -289,7 +289,7 @@ // Compute proper thresholds. static const S thmin = compute_threshold (static_cast<S> (min_val ()), min_val ()); static const S thmax = compute_threshold (static_cast<S> (max_val ()), max_val ()); - if (lo_ieee_isnan (value)) + if (xisnan (value)) { fnan = true; return static_cast<T> (0);
--- a/liboctave/oct-norm.cc +++ b/liboctave/oct-norm.cc @@ -413,6 +413,7 @@ RR lambda = 0, mu = 0; for (octave_idx_type k = 0; k < m.columns (); k++) { + OCTAVE_QUIT; VectorT col (m.column (k)); if (k > 0) higham_subp (y, col, 4*k, p, lambda, mu); @@ -430,6 +431,7 @@ int iter = 0; while (iter < maxiter) { + OCTAVE_QUIT; y = m*x; gamma1 = gamma; gamma = vector_norm (y, p); @@ -449,9 +451,6 @@ return gamma; } -// TODO: is there a better way? -inline float get_eps (float) { return FLT_EPSILON; } -inline double get_eps (double) { return DBL_EPSILON; } // derive column vector and SVD types static const char *p_less1_gripe = "xnorm: p must be at least 1"; @@ -477,7 +476,8 @@ else if (p > 1) { VectorT x; - res = higham (m, p, std::sqrt (get_eps (p)), max_norm_iter, x); + const R sqrteps = std::sqrt (std::numeric_limits<R>::epsilon ()); + res = higham (m, p, sqrteps, max_norm_iter, x); } else (*current_liboctave_error_handler) (p_less1_gripe); @@ -497,7 +497,8 @@ else if (p > 1) { VectorT x; - res = higham (m, p, std::sqrt (get_eps (p)), max_norm_iter, x); + const R sqrteps = std::sqrt (std::numeric_limits<R>::epsilon ()); + res = higham (m, p, sqrteps, max_norm_iter, x); } else (*current_liboctave_error_handler) (p_less1_gripe);
--- a/scripts/ChangeLog +++ b/scripts/ChangeLog @@ -1,3 +1,52 @@ +2009-03-25 Kai Habel <kai.habel@gmx.de> + + * general/gradient.m: Fix calculation for more than two + dimensions. Change interpretation of vector arguments from + spacing to coordinates. New tests. + +2009-03-25 John W. Eaton <jwe@octave.org> + + * mkdoc: Pass full file name to gethelp. + * gethelp.cc (main): Handle second argument. Write comment with + full file name to output. + +2009-03-24 Ben Abbott <bpabbott@mac.com> + + * plot/gnuplot_drawnow.m: When printing, pass scalar plot_stream + to __gnuplot_draw_figure__, and close all plot streams when done. + +2009-03-24 John W. Eaton <jwe@octave.org> + + * general/isa.m: Handle parent classes. + +2009-03-23 Ben Abbott <bpabbott@mac.com> + + * plot/gnuplot_drawnow.m: Check that gnuplot has internal variable + "GPVAL_TERMINALS". + * plot/__gnuplot_has_feature__.m: Add "variable_GPVAL_TERMINALS". + +2009-03-21 Ben Abbott <bpabbott@mac.com> + + * plot/gnuplot_drawnow.m: Verify the gnuplot terminal is supported. + * plot/__gnuplot_get_var__.m: Add function to get gnuplot variables. + * plot/print.m: Restore the behavior for option -S<num>,<num>. + +2009-03-19 Jaroslav Hajek <highegg@gmail.com> + + * optimization/fsolve.m (guarded_eval): Simplify & fix missing + semicolon. + +2009-03-17 Jaroslav Hajek <highegg@gmail.com> + + * optimization/__fdjac__.m: Pass in fvec to save one evaluation. + * optimization/fsolve.m: Avoid redundant reevaluation when using + FD jacobians. Document how it can be done with user jacobians. Make + first iteration special and call outputfcn after it. Skip updates + unless two successful iterations have occured. + * optimization/__dogleg__.m: Add missing alpha in the zero-gradient + case. + * optimization/fsolve.m: Remove autodg (not used), simplify. + 2009-03-14 Jaroslav Hajek <highegg@gmail.com> * statistics/base/var.m: a -> x.
--- a/scripts/general/gradient.m +++ b/scripts/general/gradient.m @@ -17,36 +17,36 @@ ## <http://www.gnu.org/licenses/>. ## -*- texinfo -*- -## @deftypefn {Function File} {@var{x} =} gradient (@var{m}) -## @deftypefnx {Function File} {[@var{x}, @var{y}, @dots{}] =} gradient (@var{m}) +## @deftypefn {Function File} {@var{dx} =} gradient (@var{m}) +## @deftypefnx {Function File} {[@var{dx}, @var{dy}, @var{dz}, @dots{}] =} gradient (@var{m}) ## @deftypefnx {Function File} {[@dots{}] =} gradient (@var{m}, @var{s}) -## @deftypefnx {Function File} {[@dots{}] =} gradient (@var{m}, @var{dx}, @var{dy}, @dots{}) +## @deftypefnx {Function File} {[@dots{}] =} gradient (@var{m}, @var{x}, @var{y}, @var{z}, @dots{}) ## @deftypefnx {Function File} {[@dots{}] =} gradient (@var{f}, @var{x0}) ## @deftypefnx {Function File} {[@dots{}] =} gradient (@var{f}, @var{x0}, @var{s}) -## @deftypefnx {Function File} {[@dots{}] =} gradient (@var{f}, @var{x0}, @var{dx}, @var{dy}, @dots{}) +## @deftypefnx {Function File} {[@dots{}] =} gradient (@var{f}, @var{x0}, @var{x}, @var{y}, @dots{}) ## ## Calculate the gradient of sampled data, or of a function. If @var{m} ## is a vector, calculate the one dimensional gradient of @var{m}. If -## @var{m} is a matrix the gradient is calculated for each row. +## @var{m} is a matrix the gradient is calculated for each dimension. ## -## @code{[@var{x}, @var{y}] = gradient (@var{m})} calculates the one -## dimensional gradient for each direction if @var{m} if @var{m} is a +## @code{[@var{dx}, @var{dy}] = gradient (@var{m})} calculates the one +## dimensional gradient for @var{x} and @var{y} direction if @var{m} is a ## matrix. Additional return arguments can be use for multi-dimensional ## matrices. ## -## Spacing values between two points can be provided by the -## @var{dx}, @var{dy} or @var{h} parameters. If @var{h} is supplied it -## is assumed to be the spacing in all directions. Otherwise, separate -## values of the spacing can be supplied by the @var{dx}, etc variables. -## A scalar value specifies an equidistant spacing, while a vector value -## can be used to specify a variable spacing. The length must match -## their respective dimension of @var{m}. +## A constant spacing between two points can be provided by the +## @var{s} parameter. If @var{s} is a scalar, it is assumed to be the spacing +## for all dimensions. +## Otherwise, separate values of the spacing can be supplied by +## the @var{x}, @dots{} arguments. Scalar values specify an equidistant spacing. +## Vector values for the @var{x}, @dots{} arguments specify the coordinate for that +## dimension. The length must match their respective dimension of @var{m}. ## ## At boundary points a linear extrapolation is applied. Interior points ## are calculated with the first approximation of the numerical gradient ## ## @example -## y'(i) = 1/(x(i+1)-x(i-1)) *(y(i-1)-y(i+1)). +## y'(i) = 1/(x(i+1)-x(i-1)) * (y(i-1)-y(i+1)). ## @end example ## ## If the first argument @var{f} is a function handle, the gradient of the @@ -83,81 +83,85 @@ function varargout = matrix_gradient (m, varargin) transposed = false; if (isvector (m)) - ## make a column vector. + ## make a row vector. transposed = (size (m, 2) == 1); m = m(:)'; endif nd = ndims (m); sz = size (m); + if (length(sz) > 1) + tmp = sz(1); sz(1) = sz(2); sz(2) = tmp; + endif + if (nargin > 2 && nargin != nd + 1) print_usage () endif - + + ## cell d stores a spacing vector for each dimension d = cell (1, nd); if (nargin == 1) + ## no spacing given - assume 1.0 for all dimensions for i = 1:nd - d{i} = ones (sz(i), 1); + d{i} = ones (sz(i) - 1, 1); endfor elseif (nargin == 2) if (isscalar (varargin{1})) + ## single scalar value for all dimensions for i = 1:nd - d{i} = varargin{1} * ones (sz(i), 1); + d{i} = varargin{1} * ones (sz(i) - 1, 1); endfor else - for i = 1:nd - d{i} = varargin{1}; - endfor + ## vector for one-dimensional derivative + d{1} = diff (varargin{1}(:)); endif else + ## have spacing value for each dimension + if (length(varargin) != nd) + error ("dimensions and number of spacing values do not match."); + end for i = 1:nd if (isscalar (varargin{i})) - ## Why the hell did Matlab decide to swap these two values? - if (i == 1) - d{2} = varargin{1} * ones (sz(2), 1); - elseif (i == 2) - d{1} = varargin{2} * ones (sz(1), 1); - else - d{i} = varargin{i} * ones (sz(i), 1); - endif + d{i} = varargin{i} * ones (sz(i) - 1, 1); else - ## Why the hell did Matlab decide to swap these two values? - if (i == 1) - d{2} = varargin{1}; - elseif (i == 2) - d{1} = varargin{2}; - else - d{i} = varargin{i}; - endif + d{i} = diff (varargin{i}(:)); endif endfor endif - for i = 1:max (2, min (nd, nargout)) - mr = sz(i); - mc = prod ([sz(1:i-1), sz(i+1:nd)]); + m = shiftdim (m, 1); + for i = 1:min (nd, nargout) + mr = rows (m); + mc = numel (m) / mr; Y = zeros (size (m), class (m)); - + if (mr > 1) ## Top and bottom boundary. - Y(1,:) = diff (m(1:2,:)) / d{i}(1); - Y(mr,:) = diff (m(mr-1:mr,:)) / d{i}(mr-1); + Y(1,:) = diff (m(1:2, :)) / d{i}(1); + Y(mr,:) = diff (m(mr-1:mr, :) / d{i}(mr - 1)); endif if (mr > 2) ## Interior points. Y(2:mr-1,:) = ((m(3:mr,:) - m(1:mr-2,:)) - ./ kron (d{i}(1:mr-2) + d{i}(2:mr-1), ones (1, mc))); + ./ kron (d{i}(1:mr-2) + d{i}(2:mr-1), ones (1, mc))); endif - varargout{i} = ipermute (Y, [i:nd,1:i-1]); - m = permute (m, [2:nd,1]); + + ## turn multi-dimensional matrix in a way, that gradient + ## along x-direction is calculated first then y, z, ... + + if (i == 1) + varargout{i} = shiftdim (Y, nd - 1); + m = shiftdim (m, nd - 1); + elseif (i == 2) + varargout{i} = Y; + m = shiftdim (m, 2); + else + varargout{i} = shiftdim (Y, nd - i + 1); + m = shiftdim (m, 1); + endif endfor - ## Why the hell did Matlab decide to swap these two values? - tmp = varargout{1}; - varargout{1} = varargout{2}; - varargout{2} = tmp; - if (transposed) varargout{1} = varargout{1}.'; endif @@ -200,7 +204,7 @@ ## Calculate the gradient p0 = mat2cell (p0, num_points, ones (1, dim)); - varargout = cell (1,dim); + varargout = cell (1, dim); for d = 1:dim s = delta (d); df_dx = (f (p0{1:d-1}, p0{d}+s, p0{d+1:end}) @@ -216,7 +220,40 @@ %!test %! data = [1, 2, 4, 2]; %! dx = gradient (data); +%! dx2 = gradient (data, 0.25); +%! dx3 = gradient (data, [0.25, 0.5, 1, 3]); %! assert (dx, [1, 3/2, 0, -2]); +%! assert (dx2, [4, 6, 0, -8]); +%! assert (dx3, [4, 4, 0, -1]); +%! assert (size_equal(data, dx)); + +%!test +%! [Y,X,Z,U] = ndgrid (2:2:8,1:5,4:4:12,3:5:30); +%! [dX,dY,dZ,dU] = gradient (X); +%! assert (all(dX(:)==1)); +%! assert (all(dY(:)==0)); +%! assert (all(dZ(:)==0)); +%! assert (all(dU(:)==0)); +%! [dX,dY,dZ,dU] = gradient (Y); +%! assert (all(dX(:)==0)); +%! assert (all(dY(:)==2)); +%! assert (all(dZ(:)==0)); +%! assert (all(dU(:)==0)); +%! [dX,dY,dZ,dU] = gradient (Z); +%! assert (all(dX(:)==0)); +%! assert (all(dY(:)==0)); +%! assert (all(dZ(:)==4)); +%! assert (all(dU(:)==0)); +%! [dX,dY,dZ,dU] = gradient (U); +%! assert (all(dX(:)==0)); +%! assert (all(dY(:)==0)); +%! assert (all(dZ(:)==0)); +%! assert (all(dU(:)==5)); +%! assert (size_equal(dX, dY, dZ, dU, X, Y, Z, U)); +%! [dX,dY,dZ,dU] = gradient (U, 5.0); +%! assert (all(dU(:)==1)); +%! [dX,dY,dZ,dU] = gradient (U, 1.0, 2.0, 3.0, 2.5); +%! assert (all(dU(:)==2)); %!test %! x = 0:10;
--- a/scripts/general/isa.m +++ b/scripts/general/isa.m @@ -41,7 +41,14 @@ elseif (strcmp (cname, "numeric")) retval = any (strcmp (class (x), fnum_classes)); else - retval = strcmp (class (x), cname); + class_of_x = class (x); + retval = strcmp (class_of_x, cname); + if (! retval && isobject (x)) + parent_classes_of_x = __parent_classes__ (x); + if (! isempty (parent_classes_of_x)) + retval = any (strcmp (parent_classes_of_x, cname)); + endif + endif endif endfunction
--- a/scripts/gethelp.cc +++ b/scripts/gethelp.cc @@ -135,20 +135,26 @@ main (int argc, char **argv) { std::string name; + std::string file_name; - if (argc != 2) + if (argc != 3) { - std::cerr << "usage: gethelp name\n"; + std::cerr << "usage: gethelp name file-name\n"; return 1; } else - name = argv[1]; + { + name = argv[1]; + file_name = argv[2]; + } std::string help_text = extract_help_text (); if (! help_text.empty ()) { - std::cout << "" << name << "\n" << help_text; + std::cout << "" << name << "\n" + << "@c " << file_name << "\n" + << help_text; if (help_text[help_text.length () - 1] != '\n') std::cout << "\n";
--- a/scripts/mkdoc +++ b/scripts/mkdoc @@ -42,15 +42,17 @@ EOF $FIND $d -name '*.m' | \ - $PERL -ne 'm{(.*)/(.*)\.m}; - for (qx{./gethelp $2 < $_}) { + $PERL -n -e 'chop; + $f = "$_"; + m{(.*)/(.*)\.m}; + for (qx{./gethelp $2 "$f" < "$f"}) { s/^\s+\@/\@/ unless $i_am_in_example; s/^\s+\@group/\@group/; s/^\s+\@end\s+group/\@end\s+group/; $i_am_in_example = 1 if /\s*\@example/; $i_am_in_example = 0 if /\s*\@end\s+example/; print; - }' + }' else echo "gethelp program seems to be missing!" 1>&2 exit 1
--- a/scripts/optimization/__dogleg__.m +++ b/scripts/optimization/__dogleg__.m @@ -53,6 +53,7 @@ alpha = 0; endif else + alpha = delta / xn; snm = 0; endif ## Form the appropriate convex combination.
--- a/scripts/optimization/__fdjac__.m +++ b/scripts/optimization/__fdjac__.m @@ -17,21 +17,19 @@ ## <http://www.gnu.org/licenses/>. ## -*- texinfo -*- -## @deftypefn{Function File} {[@var{fval}, @var{fjac}]} = __fdjac__ (@var{fcn}, @var{x}, @var{err}) +## @deftypefn{Function File} {} __fdjac__ (@var{fcn}, @var{x}, @var{fvec}, @var{err}) ## Undocumented internal function. ## @end deftypefn -function [fvec, fjac] = __fdjac__ (fcn, x, err = 0) +function fjac = __fdjac__ (fcn, x, fvec, err = 0) err = sqrt (max (eps, err)); - fvec = fcn (x); - fv = fvec(:); h = abs (x(:))*err; h(h == 0) = err; - fjac = zeros (length (fv), numel (x)); + fjac = zeros (length (fvec), numel (x)); for i = 1:numel (x) x1 = x; x1(i) += h(i); - fjac(:,i) = (fcn (x1)(:) - fv) / h(i); + fjac(:,i) = (fcn (x1)(:) - fvec) / h(i); endfor endfunction
--- a/scripts/optimization/fsolve.m +++ b/scripts/optimization/fsolve.m @@ -72,10 +72,49 @@ ## @item -3 ## The trust region radius became excessively small. ## @end table -## +## ## Note: If you only have a single nonlinear equation of one variable, using ## @code{fzero} is usually a much better idea. ## @seealso{fzero, optimset} +## +## Note about user-supplied jacobians: +## As an inherent property of the algorithm, jacobian is always requested for a +## solution vector whose residual vector is already known, and it is the last +## accepted successful step. Often this will be one of the last two calls, but +## not always. If the savings by reusing intermediate results from residual +## calculation in jacobian calculation are significant, the best strategy is to +## employ OutputFcn: After a vector is evaluated for residuals, if OutputFcn is +## called with that vector, then the intermediate results should be saved for +## future jacobian evaluation, and should be kept until a jacobian evaluation +## is requested or until outputfcn is called with a different vector, in which +## case they should be dropped in favor of this most recent vector. A short +## example how this can be achieved follows: +## +## @example +## function [fvec, fjac] = my_optim_func (x, optimvalues, state) +## persistent sav = [], sav0 = []; +## if (nargin == 1) +## ## evaluation call +## if (nargout == 1) +## sav0.x = x; # mark saved vector +## ## calculate fvec, save results to sav0. +## elseif (nargout == 2) +## ## calculate fjac using sav. +## endif +## else +## ## outputfcn call. +## if (all (x == sav0.x)) +## sav = sav0; +## endif +## ## maybe output iteration status etc. +## endif +## endfunction +## +## ## .... +## +## fsolve (@@my_optim_func, x0, optimset ("OutputFcn", @@my_optim_func, @dots{})) +## @end example +### ## @end deftypefn ## PKG_ADD: __all_opts__ ("fsolve"); @@ -125,36 +164,53 @@ tolf = optimget (options, "TolFun", sqrt (macheps)); factor = 100; - ## FIXME: TypicalX corresponds to user scaling (???) - autodg = true; niter = 1; - nfev = 0; + nfev = 1; x = x0(:); info = 0; + ## Initial evaluation. + ## Handle arbitrary shapes of x and f and remember them. + fvec = fcn (reshape (x, xsiz)); + fsiz = size (fvec); + fvec = fvec(:); + fn = norm (fvec); + m = length (fvec); + n = length (x); + + if (! isempty (outfcn)) + optimvalues.iter = niter; + optimvalues.funccount = nfev; + optimvalues.fval = fn; + optimvalues.searchdirection = zeros (n, 1); + state = 'init'; + stop = outfcn (x, optimvalues, state); + if (stop) + info = -1; + break; + endif + endif + + nsuciter = 0; + ## Outer loop. while (niter < maxiter && nfev < maxfev && ! info) ## Calculate function value and Jacobian (possibly via FD). - ## Handle arbitrary shapes of x and f and remember them. if (has_jac) [fvec, fjac] = fcn (reshape (x, xsiz)); ## If the jacobian is sparse, disable Broyden updating. if (issparse (fjac)) updating = false; endif + fvec = fvec(:); nfev ++; else - [fvec, fjac] = __fdjac__ (fcn, reshape (x, xsiz)); - nfev += 1 + length (x); + fjac = __fdjac__ (fcn, reshape (x, xsiz), fvec); + nfev += length (x); endif - fsiz = size (fvec); - fvec = fvec(:); - fn = norm (fvec); - m = length (fvec); - n = length (x); ## For square and overdetermined systems, we update a QR ## factorization of the jacobian to avoid solving a full system in each @@ -173,45 +229,34 @@ ## Get column norms, use them as scaling factors. jcn = norm (fjac, 'columns').'; if (niter == 1) - if (autodg) - dg = jcn; - dg(dg == 0) = 1; - endif + dg = jcn; + dg(dg == 0) = 1; xn = norm (dg .* x); ## FIXME: something better? delta = max (factor * xn, 1); endif - ## Rescale if necessary. - if (autodg) - ## FIXME: the original minpack used the following rescaling strategy: - ## dg = max (dg, jcn); - ## but it seems not good if we start with a bad guess yielding jacobian - ## columns with large norms that later decrease, because the corresponding - ## variable will still be overscaled. So instead, we only give the old - ## scaling a small momentum, but do not honor it. + ## Rescale adaptively. + ## FIXME: the original minpack used the following rescaling strategy: + ## dg = max (dg, jcn); + ## but it seems not good if we start with a bad guess yielding jacobian + ## columns with large norms that later decrease, because the corresponding + ## variable will still be overscaled. So instead, we only give the old + ## scaling a small momentum, but do not honor it. - dg = max (0.1*dg, jcn); + dg = max (0.1*dg, jcn); - ## It also seems that in the case of fast (and inhomogeneously) changing - ## jacobian, the Broyden updates are of little use, so maybe we could - ## skip them if a big disproportional change is expected. The question is, - ## of course, how to define the above terms :) - endif + ## It also seems that in the case of fast (and inhomogeneously) changing + ## jacobian, the Broyden updates are of little use, so maybe we could + ## skip them if a big disproportional change is expected. The question is, + ## of course, how to define the above terms :) + lastratio = 0; nfail = 0; nsuc = 0; ## Inner loop. while (niter <= maxiter && nfev < maxfev && ! info) - if (useqr) - tr_mat = r; - tr_vec = q'*fvec; - else - tr_mat = fjac; - tr_vec = fvec; - endif - ## Get trust-region model (dogleg) minimizer. if (useqr) qtf = q'*fvec; @@ -249,7 +294,7 @@ endif ## Update delta. - if (ratio < 0.1) + if (ratio < min(max(0.1, lastratio), 0.9)) nsuc = 0; nfail ++; delta *= 0.5; @@ -259,6 +304,7 @@ break; endif else + lastratio = ratio; nfail = 0; nsuc ++; if (abs (1-ratio) <= 0.1) @@ -274,6 +320,7 @@ xn = norm (dg .* x); fvec = fvec1; fn = fn1; + nsuciter ++; endif niter ++; @@ -321,7 +368,7 @@ endif ## Criterion for recalculating jacobian. - if (! updating || nfail == 2) + if (! updating || nfail == 2 || nsuciter < 2) break; endif @@ -347,6 +394,7 @@ fvec = reshape (fvec, fsiz); output.iterations = niter; + output.successful = nsuciter; output.funcCount = nfev; endfunction @@ -358,12 +406,12 @@ [fx, jx] = fun (x); else fx = fun (x); - jx = [] + jx = []; endif - if (! complexeqn && ! (all (isreal (fx(:))) && all (isreal (jx(:))))) + if (! complexeqn && ! (isreal (fx) && isreal (jx))) error ("fsolve:notreal", "fsolve: non-real value encountered"); - elseif (complexeqn && ! (all (isnumeric (fx(:))) && all (isnumeric(jx(:))))) + elseif (complexeqn && ! (isnumeric (fx) && isnumeric(jx))) error ("fsolve:notnum", "fsolve: non-numeric value encountered"); elseif (any (isnan (fx(:)))) error ("fsolve:isnan", "fsolve: NaN value encountered");
--- a/scripts/plot/Makefile.in +++ b/scripts/plot/Makefile.in @@ -47,6 +47,7 @@ __errcomm__.m \ __errplot__.m \ __ezplot__.m \ + __gnuplot_get_var__.m \ __gnuplot_has_feature__.m \ __go_close_all__.m \ __go_draw_axes__.m \
new file mode 100644 --- /dev/null +++ b/scripts/plot/__gnuplot_get_var__.m @@ -0,0 +1,151 @@ +## Copyright (C) 2009 Ben Abbott +## +## 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/>. + +## -*- texinfo -*- +## @deftypefn {Function File} {@var{value} =} __gnuplot_get_var__ (@var{h}, @var{name}, @var{fmt}) +## Undocumented internal function. +## @end deftypefn + +## Author: Ben Abbott <bpabbott@mac.com> +## Created: 2009-02-07 + +function gp_var_value = __gnuplot_get_var__ (h, gp_var_name, fmt) + + if (nargin == 0) + h = gcf (); + endif + if (nargin < 2) + print_usage (); + endif + if (nargin < 3) + fmt = ''; + endif + + if (numel (h) == 1 && isfigure (h)) + ostream = get (h, "__plot_stream__"); + else + ostream = h; + endif + if (numel (ostream) < 1) + error ("__gnuplot_get_var__: stream to gnuplot not open"); + elseif (ispc ()) # || true + if (numel (ostream) == 1) + error ("__gnuplot_get_var__: Need mkfifo that is not implemented under Windows"); + endif + use_mkfifo = false; + istream = ostream(2); + ostream = ostream(1); + else + use_mkfifo = true; + ostream = ostream(1); + endif + + if (use_mkfifo) + gpin_name = tmpnam (); + + ## Mode: 6*8*8 == 0600 + [err, msg] = mkfifo (gpin_name, 6*8*8); + + if (err != 0) + error ("__gnuplot_get_var__: Can not open fifo (%s)", msg); + endif + endif + + gp_var_name = strtrim (gp_var_name); + n = min (strfind (gp_var_name, " "), strfind (gp_var_name, ",")) - 1; + if (isempty (n)) + n = numel (gp_var_name); + endif + + unwind_protect + + ## Notes: Variables may be undefined if user closes gnuplot by "q" + ## or Alt-F4. Further, this abrupt close also requires the leading + ## "\n" on the next line. + if (use_mkfifo) + fprintf (ostream, "\nset print \"%s\";\n", gpin_name); + fflush (ostream); + [gpin, err] = fopen (gpin_name, "r"); + if (err != 0) + error ("__gnuplot_get_var__: Can not open fifo (%s)", msg); + endif + gp_cmd = sprintf ("\nif (exists(\"%s\")) print %s; else print NaN\n", + gp_var_name(1:n), gp_var_name); + fputs (ostream, gp_cmd); + + ## Close output file, to force it to be flushed + fputs (ostream, "set print;\n"); + fflush (ostream); + + ## Now read from fifo. + reading = true; + str = {}; + while (reading) + str{end+1} = fgets (gpin); + if (isnumeric (str{end}) && (str{end} == -1)) + reading = false; + str = str(1:(end-1)); + endif + endwhile + str = strcat (str{:}); + fclose (gpin); + else + ## Direct gnuplot to print to <STDOUT> + fprintf (ostream, "set print \"-\";\n"); + fflush (ostream); + gp_cmd = sprintf ("\nif (exists(\"%s\")) print \"OCTAVE: \", %s; else print NaN\n", + gp_var_name(1:n), gp_var_name); + fputs (ostream, gp_cmd); + fflush (ostream); + ## Direct gnuplot to print to <STDERR> + fputs (ostream, "set print;\n"); + fflush (ostream); + + str = {}; + while (isempty (str)) + str = char (fread (istream)'); + if (! isempty (str)) + str = regexp (str, "OCTAVE:.*", "match") + str = str{end}(8:end); + endif + fclear (istream); + endwhile + endif + + ## Strip out EOLs and the continuation character "|" + str(str=="\n") = ""; + str(str=="\r") = ""; + n_continue = strfind (str, " \\ "); + if (! isempty (n_continue)) + str(n_continue+1) = ""; + endif + + if (isempty (fmt)) + gp_var_value = strtrim (str); + else + gp_var_value = sscanf (str, fmt); + endif + + unwind_protect_cleanup + if (use_mkfifo) + unlink (gpin_name); + endif + end_unwind_protect + +endfunction +
--- a/scripts/plot/__gnuplot_has_feature__.m +++ b/scripts/plot/__gnuplot_has_feature__.m @@ -29,12 +29,13 @@ "transparent_surface", "epslatex_implies_eps_filesuffix", "epslatexstandalone_terminal", - "screen_coordinates_for_{lrtb}margin",}; + "screen_coordinates_for_{lrtb}margin", + "variable_GPVAL_TERMINALS"}; if (isempty (has_features)) gnuplot_version = __gnuplot_version__ (); - versions = {"4.2.4", "4.3", "4.3", "4.2", "4.2", "4.3"}; - operators = {">", ">=", ">=", ">=", ">=", ">="}; + versions = {"4.2.4", "4.3", "4.3", "4.2", "4.2", "4.3", "4.3"}; + operators = {">", ">=", ">=", ">=", ">=", ">=", ">="}; have_features = logical (zeros (size (features))); for n = 1 : numel (have_features) has_features(n) = compare_versions (gnuplot_version, versions{n}, operators{n});
--- a/scripts/plot/gnuplot_drawnow.m +++ b/scripts/plot/gnuplot_drawnow.m @@ -40,17 +40,34 @@ fid = []; printing = ! output_to_screen (gnuplot_trim_term (term)); unwind_protect - plot_stream = open_gnuplot_stream (1, []); - [enhanced, implicit_margin] = gnuplot_set_term (plot_stream (1), true, h, term, file); - __go_draw_figure__ (h, plot_stream, enhanced, mono, printing, implicit_margin); - if (nargin == 5) - fid = fopen (debug_file, "wb"); - enhanced = gnuplot_set_term (fid, true, h, term, file); - __go_draw_figure__ (h, fid, enhanced, mono, printing, implicit_margin); + plot_stream = open_gnuplot_stream (2, []); + if (__gnuplot_has_feature__ ("variable_GPVAL_TERMINALS")) + available_terminals = __gnuplot_get_var__ (plot_stream, "GPVAL_TERMINALS"); + available_terminals = regexp (available_terminals, "\\b\\w+\\b", "match"); + gnuplot_supports_term = any (strcmpi (available_terminals, + gnuplot_trim_term (term))); + else + gnuplot_supports_term = true; + endif + if (gnuplot_supports_term) + [enhanced, implicit_margin] = gnuplot_set_term (plot_stream (1), true, + h, term, file); + __go_draw_figure__ (h, plot_stream(1), enhanced, mono, printing, implicit_margin); + if (nargin == 5) + fid = fopen (debug_file, "wb"); + [enhanced, implicit_margin] = gnuplot_set_term (fid, true, h, term, file); + __go_draw_figure__ (h, fid, enhanced, mono, printing, implicit_margin); + endif + else + error ("gnuplot_drawnow: the gnuplot terminal, \"%s\", is not available.", + gnuplot_trim_term (term)) endif unwind_protect_cleanup if (! isempty (plot_stream)) - pclose (plot_stream); + pclose (plot_stream(1)); + if (numel (plot_stream) == 2) + pclose (plot_stream(2)); + endif endif if (! isempty (fid)) fclose (fid);
--- a/scripts/plot/print.m +++ b/scripts/plot/print.m @@ -427,11 +427,36 @@ endif set (gcf, "__pixels_per_inch__", resolution) - if (debug) - drawnow (new_terminal, name, mono, debug_file); - else - drawnow (new_terminal, name, mono); - endif + unwind_protect + if (! isempty (size)) + size_in_pixels = sscanf (size ,"%d, %d"); + size_in_pixels = reshape (size_in_pixels, [1, numel(size_in_pixels)]); + size_in_inches = size_in_pixels ./ resolution; + p.paperunits = get (gcf, "paperunits"); + p.papertype = get (gcf, "papertype"); + p.papersize = get (gcf, "papersize"); + p.paperposition = get (gcf, "paperposition"); + p.paperpositionmode = get (gcf, "paperpositionmode"); + set (gcf, "paperunits", "inches"); + set (gcf, "papertype", "<custom>"); + set (gcf, "papersize", size_in_inches); + set (gcf, "paperposition", [0, 0, size_in_inches]); + set (gcf, "paperpositionmode", "manual"); + endif + if (debug) + drawnow (new_terminal, name, mono, debug_file); + else + drawnow (new_terminal, name, mono); + endif + unwind_protect_cleanup + ## FIXME - it would be preferred to delete the added properties here. + if (! isempty (size)) + props = fieldnames (p); + for n = 1:numel(props) + set (gcf, props{n}, p.(props{n})) + endfor + endif + end_unwind_protect if (! isempty (convertname)) command = sprintf ("convert '%s' '%s'", name, convertname);
--- a/scripts/testfun/assert.m +++ b/scripts/testfun/assert.m @@ -141,7 +141,7 @@ iserror = 1; coda = "Dimensions don't match"; - elseif (tol == 0 && ! strcmp (typeinfo (cond), typeinfo (expected))) + elseif (nargin < 3 && ! strcmp (typeinfo (cond), typeinfo (expected))) iserror = 1; coda = cstrcat ("Type ", typeinfo (cond), " != ", typeinfo (expected));
--- a/src/ChangeLog +++ b/src/ChangeLog @@ -1,3 +1,133 @@ +2009-03-25 John W. Eaton <jwe@octave.org> + + * DLD-FUNCTIONS/find.cc (find_nonzero_elem_idx): Also return + [](0x0) if the array has 0 rows and it is not a column vector. + + * oct-stream.cc (octave_stream::write (const Array<T>&, + octave_idx_type, oct_data_conv::data_type, octave_idx_type, + oct_mach_info::float_format)): Seek to skip if still inside bounds + of existing file. Otherwise, write NUL to skip. + + * Makefile.in (%.df : %.cc): Write source file name to output, + wrapped in XDEFUN_FILE_NAME macro. + * mkbuiltins: Provide definition for XDEFUN_FILE_NAME. + * mkgendoc: Likeiwse. + (XDEFUN_DLD_INTERNAL, XDEFUNX_DLD_INTERNAL, XDEFUN_INTERNAL, + XDEFCONSTFUN_INTERNAL, XDEFUNX_INTERNAL, XDEFVAR_INTERNAL, + XDEFCONST_INTERNAL): Pass file_name to print_doc_string. + (print_doc_string): New arg, FILE_NAME. Print file name as + comment. + +2009-03-24 John W. Eaton <jwe@octave.org> + + * ov-class.cc (F__parent_classes__): New function. + +2009-03-24 Robert T. Short <octave@phaselockedsystems.com> + + * ov-class.h, ov-class.cc (octave_class::octave_class (const + Octave_map&, const std::string&, const octave_value_list&)): + New constructor. + (octave_class::find_parent_class, octave_class::parent_classes): + New functions. + (octave_class::dotref): Also look in parent class. + (Fclass): Handle parent class arguments. + + * ov-base.h (octave_base_value::get_parent_list, + octave_base_value::parent_classes): New virtual functions. + + * load-path.h, load-path.cc (load_path::parent_map): New data member. + (load_path::add_to_parent_map): New static function. + (load_path::do_add_to_parent_map): New member function. + (load_path::do_find_method): Also look in parent classes for methods. + (load_path::parent_map_type, load_path::const_parent_map_iterator, + load_path::parent_map_iterator): New typedefs. + +2009-03-23 John W. Eaton <jwe@octave.org> + + * symtab.h + (symbol_table::symbol_record::symobl_recoord_rep::is_variable): + Also return true if symbol is tagged as a variable. + (symbol_table::symbol_record::symobl_recoord_rep::force_variable): + Don't set variable value. + (symbol_table::symbol_record::symobl_recoord_rep::clear_forced, + symbol_table::symbol_record::clear_forced): Delete. + (symbol_table::unmark_forced_variables): Rename from + symbol_table::clear_forced_variables. + (symbol_table::do_unmark_forced_variables): Rename from + symbol_table::do_clear_forced_variables. + * parse.y (make_script, finish_function): Call + symbol_table::unmark_forced_variables instead of + symbol_table::clear_forced_variables. + * octave.cc (unmark_forced_vars): New function. + (execute_eval_option_code): Add it to the unwind-protect stack. + +2009-03-22 Jaroslav Hajek <highegg@gmail.com> + + * pt-eval.cc (tree_evaluator::visit_simple_for_command): + Remove struct branch, handle structs by the generic code. + (tree_evaluator::visit_complex_for_command): + Add missing const qualifiers. + +2009-03-21 Jaroslav Hajek <highegg@gmail.com> + + * ov-base-diag.cc: Add missing include. + * sparse-xdiv.cc: Ditto. + * DLD-FUNCTIONS/__glpk__.cc: Ditto. + * DLD-FUNCTIONS/__lin_interpn__.cc: Ditto. + * DLD-FUNCTIONS/__voronoi__.cc: Ditto. + +2009-03-20 Jaroslav Hajek <highegg@gmail.com> + + * ov-re-mat.cc (octave_matrix::load_ascii): Simplify. + * ov-flt-re-mat.cc (octave_float_matrix::load_ascii): Simplify. + * ov-cx-mat.cc (octave_complex_matrix::load_ascii): Simplify. + * ov-flt-cx-mat.cc (octave_float_complex_matrix::load_ascii): Simplify. + +2009-03-20 Jaroslav Hajek <highegg@gmail.com> + + * ov-re-mat.cc (octave_matrix::isnan, octave_matrix::isinf, + octave_matrix::finite): Simplify. + * ov-flt-re-mat.cc (octave_float_matrix::isnan, + octave_float_matrix::isinf, octave_float_matrix::finite): Simplify. + * ov-cx-mat.cc (octave_complex_matrix::isnan, + octave_complex_matrix::isinf, octave_complex_matrix::finite): + Simplify. + * ov-flt-cx-mat.cc (octave_float_complex_matrix::isnan, + octave_float_complex_matrix::isinf, + octave_float_complex_matrix::finite): Simplify. + +2009-03-19 Benjamin Lindner <lindnerb@users.sourceforge.net> + + * ls-oct-ascii.cc (extract_keyword): Replace loop with call to + read_until_newline to avoid leaving stray '\r' in stream when + reading files with CRLF line endings. + +2009-03-18 Jaroslav Hajek <highegg@gmail.com> + + * DLD-FUNCTIONS/cellfun.cc + (scalar_col_helper, scalar_col_helper_def, scalar_col_helper_nda, + scalar_query_helper): New classes. + (make_col_helper): New function. + (Fcellfun): Use col_helper classes for UniformOutput. Avoid copying by + using const variables. + +2009-03-17 Jaroslav Hajek <highegg@gmail.com> + + * OPERATORS/op-b-bm.cc: Add missing INSTALL_BINOPs. + * OPERATORS/op-bm-b.cc: Ditto. + * OPERATORS/op-fm-fm.cc: Ditto. + * OPERATORS/op-m-m.cc: Ditto. + +2009-03-17 Jaroslav Hajek <highegg@gmail.com> + + * ov.cc (octave_value::octave_value): Move to ov.h + * ov.h (octave_value::octave_value): Implement fast inline constructor + using a static rep. + +2009-03-15 Jaroslav Hajek <highegg@gmail.com> + + * DLD-FUNCTIONS/cellfun.cc (Fcellslices): New DLD function. + 2009-03-13 Jaroslav Hajek <highegg@gmail.com> * ov.h (octave_value::compound_binary_op): Support bool compound ops.
--- a/src/DLD-FUNCTIONS/__glpk__.cc +++ b/src/DLD-FUNCTIONS/__glpk__.cc @@ -28,6 +28,8 @@ #include <csetjmp> #include <ctime> +#include "lo-ieee.h" + #include "defun-dld.h" #include "error.h" #include "gripes.h"
--- a/src/DLD-FUNCTIONS/__lin_interpn__.cc +++ b/src/DLD-FUNCTIONS/__lin_interpn__.cc @@ -24,6 +24,7 @@ #include <config.h> #endif +#include "lo-ieee.h" #include "dNDArray.h" #include "oct-locbuf.h"
--- a/src/DLD-FUNCTIONS/__voronoi__.cc +++ b/src/DLD-FUNCTIONS/__voronoi__.cc @@ -34,6 +34,7 @@ #ifdef HAVE_CONFIG_H #include <config.h> #endif +#include "lo-ieee.h" #include "oct.h" #include "Cell.h"
--- a/src/DLD-FUNCTIONS/cellfun.cc +++ b/src/DLD-FUNCTIONS/cellfun.cc @@ -1,6 +1,7 @@ /* Copyright (C) 2005, 2006, 2007, 2008, 2009 Mohamed Kamoun +Copyright (C) 2009 Jaroslav Hajek This file is part of Octave. @@ -27,6 +28,7 @@ #include <string> #include <vector> #include <list> +#include <memory> #include "lo-mappers.h" #include "oct-locbuf.h" @@ -39,6 +41,142 @@ #include "ov-colon.h" #include "unwind-prot.h" +// Rationale: +// The octave_base_value::subsasgn method carries too much overhead for +// per-element assignment strategy. +// This class will optimize the most optimistic and most likely case +// when the output really is scalar by defining a hierarchy of virtual +// collectors specialized for some scalar types. + +class scalar_col_helper +{ +public: + virtual bool collect (octave_idx_type i, const octave_value& val) = 0; + virtual octave_value result (void) = 0; + virtual ~scalar_col_helper (void) { } +}; + +// The default collector represents what was previously done in the main loop. +// This reuses the existing assignment machinery via octave_value::subsasgn, +// which can perform all sorts of conversions, but is relatively slow. + +class scalar_col_helper_def : public scalar_col_helper +{ + std::list<octave_value_list> idx_list; + octave_value resval; +public: + scalar_col_helper_def (const octave_value& val, const dim_vector& dims) + : idx_list (1), resval (val) + { + idx_list.front ().resize (1); + if (resval.dims () != dims) + resval.resize (dims); + } + ~scalar_col_helper_def (void) { } + + bool collect (octave_idx_type i, const octave_value& val) + { + if (val.numel () == 1) + { + idx_list.front ()(0) = static_cast<double> (i + 1); + resval = resval.subsasgn ("(", idx_list, val); + } + else + error ("cellfun: expecting all values to be scalars for UniformOutput = true"); + + return true; + } + octave_value result (void) + { + return resval; + } +}; + +template <class T> +struct scalar_query_helper { }; + +#define DEF_QUERY_HELPER(T, TEST, QUERY) \ +template <> \ +struct scalar_query_helper<T> \ +{ \ + static bool has_value (const octave_value& val) \ + { return TEST; } \ + static T get_value (const octave_value& val) \ + { return QUERY; } \ +} + +DEF_QUERY_HELPER (double, val.is_real_scalar (), val.scalar_value ()); +DEF_QUERY_HELPER (Complex, val.is_complex_scalar (), val.complex_value ()); +DEF_QUERY_HELPER (float, val.is_single_type () && val.is_real_scalar (), + val.float_scalar_value ()); +DEF_QUERY_HELPER (FloatComplex, val.is_single_type () && val.is_complex_scalar (), + val.float_complex_value ()); +DEF_QUERY_HELPER (bool, val.is_bool_scalar (), val.bool_value ()); +// FIXME: More? + +// This specializes for collecting elements of a single type, by accessing +// an array directly. If the scalar is not valid, it returns false. + +template <class NDA> +class scalar_col_helper_nda : public scalar_col_helper +{ + NDA arrayval; + typedef typename NDA::element_type T; +public: + scalar_col_helper_nda (const octave_value& val, const dim_vector& dims) + : arrayval (dims) + { + arrayval(0) = scalar_query_helper<T>::get_value (val); + } + ~scalar_col_helper_nda (void) { } + + bool collect (octave_idx_type i, const octave_value& val) + { + bool retval = scalar_query_helper<T>::has_value (val); + if (retval) + arrayval(i) = scalar_query_helper<T>::get_value (val); + return retval; + } + octave_value result (void) + { + return arrayval; + } +}; + +template class scalar_col_helper_nda<NDArray>; +template class scalar_col_helper_nda<FloatNDArray>; +template class scalar_col_helper_nda<ComplexNDArray>; +template class scalar_col_helper_nda<FloatComplexNDArray>; +template class scalar_col_helper_nda<boolNDArray>; + +// the virtual constructor. +scalar_col_helper * +make_col_helper (const octave_value& val, const dim_vector& dims) +{ + scalar_col_helper *retval; + + if (val.is_bool_scalar ()) + retval = new scalar_col_helper_nda<boolNDArray> (val, dims); + else if (val.is_complex_scalar ()) + { + if (val.is_single_type ()) + retval = new scalar_col_helper_nda<FloatComplexNDArray> (val, dims); + else + retval = new scalar_col_helper_nda<ComplexNDArray> (val, dims); + } + else if (val.is_real_scalar ()) + { + if (val.is_single_type ()) + retval = new scalar_col_helper_nda<FloatNDArray> (val, dims); + else + retval = new scalar_col_helper_nda<NDArray> (val, dims); + } + else + retval = new scalar_col_helper_def (val, dims); + + return retval; +} + DEFUN_DLD (cellfun, args, nargout, "-*- texinfo -*-\n\ @deftypefn {Loadable Function} {} cellfun (@var{name}, @var{c})\n\ @@ -163,7 +301,7 @@ return retval; } - Cell f_args = args(1).cell_value (); + const Cell f_args = args(1).cell_value (); octave_idx_type k = f_args.numel (); @@ -202,7 +340,7 @@ result(count) = static_cast<double> (f_args.elem(count).ndims ()); retval(0) = result; } - else if (name == "prodofsize") + else if (name == "prodofsize" || name == "numel") { NDArray result (f_args.dims ()); for (octave_idx_type count = 0; count < k ; count++) @@ -270,7 +408,6 @@ error ("unknown function"); else { - octave_value_list idx; octave_value_list inputlist; bool uniform_output = true; bool have_error_handler = false; @@ -279,6 +416,8 @@ int offset = 1; int i = 1; OCTAVE_LOCAL_BUFFER (Cell, inputs, nargin); + // This is to prevent copy-on-write. + const Cell *cinputs = inputs; while (i < nargin) { @@ -344,19 +483,20 @@ } } - inputlist.resize(nargin-offset); + nargin -= offset; + inputlist.resize(nargin); if (have_error_handler) buffer_error_messages++; if (uniform_output) { - retval.resize(nargout); + OCTAVE_LOCAL_BUFFER (std::auto_ptr<scalar_col_helper>, retptr, nargout); for (octave_idx_type count = 0; count < k ; count++) { - for (int j = 0; j < nargin-offset; j++) - inputlist(j) = inputs[j](count); + for (int j = 0; j < nargin; j++) + inputlist(j) = cinputs[j](count); octave_value_list tmp = feval (func, inputlist, nargout); @@ -390,38 +530,46 @@ { for (int j = 0; j < nargout; j++) { - octave_value val; - val = tmp(j); + octave_value val = tmp(j); - if (error_state) - goto cellfun_err; - - retval(j) = val.resize(f_args.dims()); + if (val.numel () == 1) + retptr[j].reset (make_col_helper (val, f_args.dims ())); + else + { + error ("cellfun: expecting all values to be scalars for UniformOutput = true"); + break; + } } } else { - idx(0) = octave_value (static_cast<double>(count+1)); for (int j = 0; j < nargout; j++) { - // FIXME -- need an easier way to express - // this test. octave_value val = tmp(j); - if (val.ndims () == 2 - && val.rows () == 1 && val.columns () == 1) - retval(j) = - retval(j).subsasgn ("(", - std::list<octave_value_list> - (1, idx(0)), val); - else - error ("cellfun: expecting all values to be scalars for UniformOutput = true"); - } + if (! retptr[j]->collect (count, val)) + { + // FIXME: A more elaborate structure would allow again a virtual + // constructor here. + retptr[j].reset (new scalar_col_helper_def (retptr[j]->result (), + f_args.dims ())); + retptr[j]->collect (count, val); + } + } } if (error_state) break; } + + retval.resize (nargout); + for (int j = 0; j < nargout; j++) + { + if (retptr[j].get ()) + retval(j) = retptr[j]->result (); + else + retval(j) = Matrix (); + } } else { @@ -431,8 +579,8 @@ for (octave_idx_type count = 0; count < k ; count++) { - for (int j = 0; j < nargin-offset; j++) - inputlist(j) = inputs[j](count); + for (int j = 0; j < nargin; j++) + inputlist(j) = cinputs[j](count); octave_value_list tmp = feval (func, inputlist, nargout); @@ -1069,6 +1217,97 @@ %! assert(c,{empty1by0str,'abcd','ef',empty1by0str,'ghij',empty1by0str}) */ + +template <class NDA> +Cell +do_cellslices_nda (const NDA& array, const idx_vector& lb, const idx_vector& ub) +{ + octave_idx_type n = lb.length (0); + Cell retval (1, n); + for (octave_idx_type i = 0; i < n && ! error_state; i++) + retval(i) = array.index (idx_vector (lb(i), ub(i) + 1)); + return retval; +} + +DEFUN_DLD (cellslices, args, , + "-*- texinfo -*-\n\ +@deftypefn {Loadable Function} {@var{sl} =} cellslices (@var{x}, @var{lb}, @var{ub})\n\ +Given a vector @var{x}, this function produces a cell array of slices from the vector\n\ +determined by the index vectors @var{lb}, @var{ub}, for lower and upper bounds, respectively.\n\ +In other words, it is equivalent to the following code:\n\ +\n\ +@example\n\ +n = length (lb);\n\ +sl = cell (1, n);\n\ +for i = 1:length (lb)\n\ + sl@{i@} = x(lb(i):ub(i));\n\ +endfor\n\ +@end example\n\ +\n\ +If @var{X} is a matrix, linear indexing will be used.\n\ +@seealso{mat2cell}\n\ +@end deftypefn") +{ + octave_value retval; + int nargin = args.length (); + if (nargin == 3) + { + octave_value x = args(0); + idx_vector lb = args(1).index_vector (), ub = args(2).index_vector (); + if (! error_state) + { + if (lb.is_colon () || ub.is_colon ()) + error ("cellslices: invalid use of colon"); + else if (lb.length (0) != ub.length (0)) + error ("cellslices: the lengths of lb and ub must match"); + else + { + Cell retcell; + if (! x.is_sparse_type () && x.is_matrix_type ()) + { + // specialize for some dense arrays. + if (x.is_bool_type ()) + retcell = do_cellslices_nda (x.bool_array_value (), lb, ub); + else if (x.is_char_matrix ()) + retcell = do_cellslices_nda (x.char_array_value (), lb, ub); + else if (x.is_complex_type ()) + { + if (x.is_single_type ()) + retcell = do_cellslices_nda (x.float_complex_array_value (), lb, ub); + else + retcell = do_cellslices_nda (x.complex_array_value (), lb, ub); + } + else + { + if (x.is_single_type ()) + retcell = do_cellslices_nda (x.float_array_value (), lb, ub); + else + retcell = do_cellslices_nda (x.array_value (), lb, ub); + } + } + else + { + // generic code. + octave_idx_type n = lb.length (0); + retcell = Cell (1, n); + octave_value_list idx (1, octave_value ()); + for (octave_idx_type i = 0; i < n && ! error_state; i++) + { + idx(0) = Range (static_cast<double> (lb(i)) + 1, + static_cast<double> (ub(i)) + 1); + retcell(i) = x.do_index_op (idx); + } + } + if (! error_state) + retval = retcell; + } + } + } + else + print_usage (); + + return retval; +} /* ;;; Local Variables: ***
--- a/src/DLD-FUNCTIONS/find.cc +++ b/src/DLD-FUNCTIONS/find.cc @@ -89,14 +89,23 @@ octave_idx_type result_nr = count; octave_idx_type result_nc = 1; + bool column_vector_arg = false; bool scalar_arg = false; - if (nda.ndims () == 2 && nda.rows () == 1) + if (nda.ndims () == 2) { - result_nr = 1; - result_nc = count; + octave_idx_type nr = nda.rows (); + octave_idx_type nc = nda.columns (); - scalar_arg = (nda.columns () == 1); + if (nr == 1) + { + result_nr = 1; + result_nc = count; + + scalar_arg = (nc == 1); + } + else if (nc == 1) + column_vector_arg = true; } Matrix idx (result_nr, result_nc); @@ -141,7 +150,7 @@ i++; } } - else if (scalar_arg) + else if (scalar_arg || (nda.rows () == 0 && ! column_vector_arg)) { idx.resize (0, 0);
--- a/src/Makefile.in +++ b/src/Makefile.in @@ -40,6 +40,7 @@ %.df : %.cc @echo making $@ from $< @(echo "// DO NOT EDIT! Generated automatically by mkdefs." ; \ + echo " XDEFUN_FILE_NAME (\"$<\")" ; \ egrep '^(///*|/\*) *PKG_ADD:' $< ; \ $(CXXCPP) $(CPPFLAGS) $(CXXFLAGS_NO_PT_FLAGS) -DMAKE_BUILTINS $< \ | $(srcdir)/mkdefs) > $@-t
--- a/src/OPERATORS/op-b-bm.cc +++ b/src/OPERATORS/op-b-bm.cc @@ -65,6 +65,8 @@ { INSTALL_BINOP (op_el_and, octave_bool, octave_bool_matrix, el_and); INSTALL_BINOP (op_el_or, octave_bool, octave_bool_matrix, el_or); + INSTALL_BINOP (op_el_and_not, octave_bool, octave_bool_matrix, el_and_not); + INSTALL_BINOP (op_el_or_not, octave_bool, octave_bool_matrix, el_or_not); INSTALL_CATOP (octave_bool, octave_bool_matrix, b_bm); INSTALL_CATOP (octave_bool, octave_matrix, b_m);
--- a/src/OPERATORS/op-bm-b.cc +++ b/src/OPERATORS/op-bm-b.cc @@ -86,6 +86,8 @@ { INSTALL_BINOP (op_el_and, octave_bool_matrix, octave_bool, el_and); INSTALL_BINOP (op_el_or, octave_bool_matrix, octave_bool, el_or); + INSTALL_BINOP (op_el_not_and, octave_bool_matrix, octave_bool, el_not_and); + INSTALL_BINOP (op_el_not_or, octave_bool_matrix, octave_bool, el_not_or); INSTALL_CATOP (octave_bool_matrix, octave_bool, bm_b); INSTALL_CATOP (octave_bool_matrix, octave_scalar, bm_s);
--- a/src/OPERATORS/op-fm-fm.cc +++ b/src/OPERATORS/op-fm-fm.cc @@ -202,6 +202,10 @@ INSTALL_BINOP (op_el_ldiv, octave_float_matrix, octave_float_matrix, el_ldiv); INSTALL_BINOP (op_el_and, octave_float_matrix, octave_float_matrix, el_and); INSTALL_BINOP (op_el_or, octave_float_matrix, octave_float_matrix, el_or); + INSTALL_BINOP (op_el_and_not, octave_float_matrix, octave_float_matrix, el_and_not); + INSTALL_BINOP (op_el_or_not, octave_float_matrix, octave_float_matrix, el_or_not); + INSTALL_BINOP (op_el_not_and, octave_float_matrix, octave_float_matrix, el_not_and); + INSTALL_BINOP (op_el_not_or, octave_float_matrix, octave_float_matrix, el_not_or); INSTALL_BINOP (op_trans_mul, octave_float_matrix, octave_float_matrix, trans_mul); INSTALL_BINOP (op_mul_trans, octave_float_matrix, octave_float_matrix, mul_trans); INSTALL_BINOP (op_herm_mul, octave_float_matrix, octave_float_matrix, trans_mul);
--- a/src/OPERATORS/op-m-m.cc +++ b/src/OPERATORS/op-m-m.cc @@ -175,6 +175,10 @@ INSTALL_BINOP (op_el_ldiv, octave_matrix, octave_matrix, el_ldiv); INSTALL_BINOP (op_el_and, octave_matrix, octave_matrix, el_and); INSTALL_BINOP (op_el_or, octave_matrix, octave_matrix, el_or); + INSTALL_BINOP (op_el_and_not, octave_matrix, octave_matrix, el_and_not); + INSTALL_BINOP (op_el_or_not, octave_matrix, octave_matrix, el_or_not); + INSTALL_BINOP (op_el_not_and, octave_matrix, octave_matrix, el_not_and); + INSTALL_BINOP (op_el_not_or, octave_matrix, octave_matrix, el_not_or); INSTALL_BINOP (op_trans_mul, octave_matrix, octave_matrix, trans_mul); INSTALL_BINOP (op_mul_trans, octave_matrix, octave_matrix, mul_trans); INSTALL_BINOP (op_herm_mul, octave_matrix, octave_matrix, trans_mul);
--- a/src/graphics.h.in +++ b/src/graphics.h.in @@ -35,6 +35,8 @@ #include <set> #include <string> +#include "lo-ieee.h" + #include "gripes.h" #include "oct-map.h" #include "oct-mutex.h"
--- a/src/load-path.cc +++ b/src/load-path.cc @@ -512,6 +512,7 @@ fcn_map.clear (); private_fcn_map.clear (); method_map.clear (); + parent_map.clear (); do_append (".", false); } @@ -1016,7 +1017,29 @@ const_fcn_map_iterator p = m.find (meth); - if (p != m.end ()) + if (p == m.end ()) + { + // Look in parent classes. + + const_parent_map_iterator r = parent_map.find (class_name); + + if (r != parent_map.end ()) + { + const std::list<std::string>& plist = r->second; + std::list<std::string>::const_iterator it = plist.begin (); + + while (it != plist.end ()) + { + retval = do_find_method (*it, meth, dir_name, type); + + if (retval != "") + break; + + it++; + } + } + } + else { const file_info_list_type& file_info_list = p->second; @@ -1538,6 +1561,13 @@ } void +load_path::do_add_to_parent_map (const std::string& classname, + const std::list<std::string>& parent_list) const +{ + parent_map[classname] = parent_list; +} + +void load_path::add_to_fcn_map (const dir_info& di, bool at_end) const { std::string dir_name = di.dir_name;
--- a/src/load-path.h +++ b/src/load-path.h @@ -37,7 +37,8 @@ { protected: - load_path (void) : dir_info_list (), fcn_map (), method_map () { } + load_path (void) + : dir_info_list (), fcn_map (), method_map (), parent_map () { } public: @@ -228,6 +229,13 @@ return instance_ok () ? instance->do_system_path () : std::string (); } + static void add_to_parent_map (const std::string& classname, + const std::list<std::string>& parent_list) + { + if (instance_ok ()) + instance->do_add_to_parent_map (classname, parent_list); + } + private: static const int M_FILE = 1; @@ -386,6 +394,12 @@ typedef method_map_type::const_iterator const_method_map_iterator; typedef method_map_type::iterator method_map_iterator; + + // <CLASS_NAME, PARENT_LIST>> + typedef std::map<std::string, std::list<std::string> > parent_map_type; + + typedef parent_map_type::const_iterator const_parent_map_iterator; + typedef parent_map_type::iterator parent_map_iterator; mutable dir_info_list_type dir_info_list; @@ -395,6 +409,8 @@ mutable method_map_type method_map; + mutable parent_map_type parent_map; + static load_path *instance; static hook_fcn_ptr add_hook; @@ -493,6 +509,9 @@ std::string do_get_command_line_path (void) const { return command_line_path; } + void do_add_to_parent_map (const std::string& classname, + const std::list<std::string>& parent_list) const; + void add_to_fcn_map (const dir_info& di, bool at_end) const; void add_to_private_fcn_map (const dir_info& di) const;
--- a/src/ls-oct-ascii.cc +++ b/src/ls-oct-ascii.cc @@ -110,14 +110,8 @@ while (is.get (c) && (c == ' ' || c == '\t' || c == ':')) ; // Skip whitespace and the colon. - if (c != '\n' && c != '\r') - { - value << c; - while (is.get (c) && c != '\n' && c != '\r') - value << c; - } - - retval = value.str (); + is.putback(c); + retval = read_until_newline (is, false); break; } else if (next_only)
--- a/src/mkbuiltins +++ b/src/mkbuiltins @@ -65,6 +65,8 @@ #endif +#define XDEFUN_FILE_NAME(name) + #define XDEFUN_INTERNAL(name, args_name, nargout_name, doc) \ extern DECLARE_FUN (name, args_name, nargout_name); \ install_builtin_function (F ## name, #name, doc); \
--- a/src/mkgendoc +++ b/src/mkgendoc @@ -43,42 +43,47 @@ #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, doc); + print_doc_string (#name, file_name, doc); #define XDEFUNX_DLD_INTERNAL(name, fname, args_name, nargout_name, doc) \ - print_doc_string (name, doc); + print_doc_string (name, file_name, doc); #define XDEFUN_INTERNAL(name, args_name, nargout_name, doc) \ - print_doc_string (#name, doc); + print_doc_string (#name, file_name, doc); #define XDEFCONSTFUN_INTERNAL(name, args_name, nargout_name, doc) \ - print_doc_string (#name, doc); + print_doc_string (#name, file_name, doc); #define XDEFUNX_INTERNAL(name, fname, args_name, nargout_name, doc) \ - print_doc_string (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, doc); + print_doc_string (#name, file_name, doc); #define XDEFCONST_INTERNAL(name, defn, doc) \ - print_doc_string (#name, doc); + print_doc_string (#name, file_name, doc); static void -print_doc_string (const std::string& name, const std::string& doc) +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); + std::cout << name.substr (1, len-2) << "\n"; else - std::cout << name; + std::cout << name << "\n"; - std::cout << "\n" << doc << "\n"; + std::cout << "@c " << file_name << "\n" + << doc << "\n"; } EOF
--- a/src/oct-stream.cc +++ b/src/oct-stream.cc @@ -3595,17 +3595,36 @@ { std::ostream& os = *osp; - // It seems that Matlab writes zeros instead of actually - // seeking. Hmm... - if (skip != 0 && (i % block_size) == 0) { - // FIXME -- probably should try to write larger - // blocks... - - unsigned char zero = 0; - for (int j = 0; j < skip; j++) - os.write (reinterpret_cast<const char *> (&zero), 1); + // Seek to skip when inside bounds of existing file. + // Otherwise, write NUL to skip. + + long orig_pos = tell (); + + seek (0, SEEK_END); + + long eof_pos = tell (); + + // Is it possible for this to fail to return us to the + // original position? + seek (orig_pos, SEEK_SET); + + long remaining = eof_pos - orig_pos; + + if (remaining < skip) + { + seek (0, SEEK_END); + + // FIXME -- probably should try to write larger + // blocks... + + unsigned char zero = 0; + for (octave_idx_type j = 0; j < skip - remaining; j++) + os.write (reinterpret_cast<const char *> (&zero), 1); + } + else + seek (skip, SEEK_CUR); } if (os)
--- a/src/octave.cc +++ b/src/octave.cc @@ -367,6 +367,18 @@ unwind_protect::run_frame ("execute_startup_files"); } +static void +unmark_forced_vars (void *arg) +{ + // Unmark any symbols that may have been tagged as local variables + // while parsing (for example, by force_local_variable in lex.l). + + symbol_table::scope_id *pscope = static_cast <symbol_table::scope_id *> (arg); + + if (pscope) + symbol_table::unmark_forced_variables (*pscope); +} + static int execute_eval_option_code (const std::string& code) { @@ -386,6 +398,11 @@ unwind_protect_bool (interactive); + // Do this with an unwind-protect cleanup function so that the + // forced variables will be unmarked in the event of an interrupt. + symbol_table::scope_id scope = symbol_table::top_scope (); + unwind_protect::add (unmark_forced_vars, &scope); + interactive = false; int parse_status = 0;
--- a/src/ov-base-diag.cc +++ b/src/ov-base-diag.cc @@ -27,6 +27,7 @@ #include <iostream> #include "mach-info.h" +#include "lo-ieee.h" #include "ov-base.h" #include "ov-base-mat.h"
--- a/src/ov-base-mat.h +++ b/src/ov-base-mat.h @@ -98,6 +98,8 @@ dim_vector dims (void) const { return matrix.dims (); } + octave_idx_type numel (void) const { return matrix.numel (); } + octave_idx_type nnz (void) const { return matrix.nnz (); } octave_value reshape (const dim_vector& new_dims) const
--- a/src/ov-base.h +++ b/src/ov-base.h @@ -450,6 +450,13 @@ virtual string_vector map_keys (void) const; + virtual string_vector parent_class_names (void) const + { return string_vector (); } + + // FIXME -- should this warn if called for a non-class type? + virtual octave_base_value *find_parent_class (const std::string&) + { return 0; } + virtual octave_function *function_value (bool silent = false); virtual const octave_function *function_value (bool silent = false) const;
--- a/src/ov-class.cc +++ b/src/ov-class.cc @@ -62,6 +62,85 @@ (octave_class::t_name, "<unknown>", octave_value (new octave_class ())); } +octave_class::octave_class (const Octave_map& m, const std::string& id, + const octave_value_list& parents) + : octave_base_value (), map (m), c_name (id) +{ + octave_idx_type n = parents.length (); + + for (octave_idx_type idx = 0; idx < n; idx++) + { + octave_value parent = parents(idx); + + if (! parent.is_object ()) + error ("parents must be objects"); + else + { + std::string cnm = parent.class_name (); + + parent_list.push_back (cnm); + + map.assign (cnm, parent); + } + } + + load_path::add_to_parent_map (id, parent_list); +} + +octave_base_value * +octave_class::find_parent_class (const std::string& parent_class_name) +{ + octave_base_value* retval = 0; + std::string dbg_clsnm = class_name (); + + if (parent_class_name == class_name ()) + retval = this; + else + { + // Search in the list of immediate parents first, then in the + // ancestor tree. + + std::list<std::string>::iterator + p = find (parent_list.begin (), parent_list.end (), parent_class_name); + + if (p != parent_list.end ()) + { + Octave_map::const_iterator pmap = map.seek (parent_class_name); + + if (pmap != map.end ()) + { + const Cell& tmp = pmap->second; + + octave_value vtmp = tmp(0); + + retval = vtmp.internal_rep (); + } + } + else + { + for (std::list<std::string>::iterator pit = parent_list.begin (); + pit != parent_list.end (); + pit++) + { + Octave_map::const_iterator smap = map.seek (*pit); + + const Cell& tmp = smap->second; + + octave_value vtmp = tmp(0); + + octave_base_value *obvp = vtmp.internal_rep (); + + retval = obvp->find_parent_class (parent_class_name); + + if (retval) + break; + } + } + } + + return retval; +} + Cell octave_class::dotref (const octave_value_list& idx) { @@ -69,12 +148,27 @@ assert (idx.length () == 1); + // FIXME -- Is there a "proper" way to do this? + octave_function* fcn = octave_call_stack::current (); + std::string my_dir = fcn->dir_name (); + int ipos = my_dir.find_last_of ("@"); + std::string method_class = my_dir.substr (ipos+1); + + // Find the class in which this method resides before attempting to access + // the requested field. + + octave_base_value *obvp = find_parent_class (method_class); + + Octave_map my_map; + + my_map = obvp ? obvp->map_value () : map; + std::string nm = idx(0).string_value (); - Octave_map::const_iterator p = map.seek (nm); + Octave_map::const_iterator p = my_map.seek (nm); - if (p != map.end ()) - retval = map.contents (p); + if (p != my_map.end ()) + retval = my_map.contents (p); else error ("class has no member `%s'", nm.c_str ()); @@ -1168,18 +1262,22 @@ "-*- texinfo -*-\n\ @deftypefn {Built-in Function} {} class (@var{expr})\n\ @deftypefnx {Built-in Function} {} class (@var{s}, @var{id})\n\ -\n\ -Return the class of the expression @var{expr}, as a string or\n\ -create a class object from the structure @var{s} with name @var{id}.\n\ +@deftypefnx {Built-in Function} {} class (@var{s}, @var{id}, @var{p}, @dots{})\n\ +Return the class of the expression @var{expr} or create a class with\n\ +fields from structure @var{s} and name (string) @var{id}. Additional\n\ +arguments name a list of parent classes from which the new class is\n\ +derived.\n\ @end deftypefn") { octave_value retval; int nargin = args.length (); - if (nargin == 1) + if (nargin == 0) + print_usage (); + else if (nargin == 1) retval = args(0).class_name (); - else if (nargin == 2) + else { Octave_map m = args(0).map_value (); @@ -1192,7 +1290,16 @@ octave_function *fcn = octave_call_stack::caller (); if (fcn && fcn->is_class_constructor ()) - retval = octave_value (new octave_class (m, id)); + { + if (nargin == 2) + retval = octave_value (new octave_class (m, id)); + else + { + octave_value_list parents = args.slice (2, nargin-2); + + retval = octave_value (new octave_class (m, id, parents)); + } + } else error ("class: invalid call from outside class constructor"); } @@ -1202,13 +1309,31 @@ else error ("class: expecting structure as first argument"); } + + return retval; +} + +DEFUN (__parent_classes__, args, , + "-*- texinfo -*-\n\ +@deftypefn {Built-in Function} {} __parent_classes__ (@var{x})\n\ +Undocumented internal function.\n\ +@end deftypefn") +{ + octave_value retval = Cell (); + + if (args.length () == 1) + { + octave_value arg = args(0); + + if (arg.is_object ()) + retval = Cell (arg.parent_class_names ()); + } else print_usage (); return retval; } - DEFUN (isobject, args, , "-*- texinfo -*-\n\ @deftypefn {Built-in Function} {} isobject (@var{x})\n\
--- a/src/ov-class.h +++ b/src/ov-class.h @@ -55,7 +55,11 @@ : octave_base_value (), map (m), c_name (id) { } octave_class (const octave_class& s) - : octave_base_value (s), map (s.map), c_name (s.c_name) { } + : octave_base_value (s), map (s.map), c_name (s.c_name), + parent_list (s.parent_list) { } + + octave_class (const Octave_map& m, const std::string& id, + const octave_value_list& parents); ~octave_class (void) { } @@ -118,6 +122,11 @@ string_vector map_keys (void) const; + string_vector parent_class_names (void) const + { return string_vector (parent_list); } + + octave_base_value *find_parent_class (const std::string&); + void print (std::ostream& os, bool pr_as_read_syntax = false) const; void print_raw (std::ostream& os, bool pr_as_read_syntax = false) const; @@ -165,6 +174,7 @@ static const std::string t_name; std::string c_name; + std::list<std::string> parent_list; bool in_class_method (void) const; };
--- a/src/ov-cx-mat.cc +++ b/src/ov-cx-mat.cc @@ -353,20 +353,15 @@ { ComplexNDArray tmp(dv); - if (tmp.is_empty ()) - matrix = tmp; - else - { - is >> tmp; + is >> tmp; - if (is) - matrix = tmp; - else - { - error ("load: failed to load matrix constant"); - success = false; - } - } + if (is) + matrix = tmp; + else + { + error ("load: failed to load matrix constant"); + success = false; + } } else { @@ -856,6 +851,24 @@ return ::imag (matrix); } +octave_value +octave_complex_matrix::isnan (void) const +{ + return matrix.isnan (); +} + +octave_value +octave_complex_matrix::isinf (void) const +{ + return matrix.isinf (); +} + +octave_value +octave_complex_matrix::finite (void) const +{ + return matrix.isfinite (); +} + DARRAY_MAPPER (erf, NDArray::dmapper, ::erf) DARRAY_MAPPER (erfc, NDArray::dmapper, ::erfc) DARRAY_MAPPER (gamma, NDArray::dmapper, xgamma) @@ -888,10 +901,7 @@ ARRAY_MAPPER (sqrt, ComplexNDArray::cmapper, std::sqrt) ARRAY_MAPPER (tan, ComplexNDArray::cmapper, std::tan) ARRAY_MAPPER (tanh, ComplexNDArray::cmapper, std::tanh) -ARRAY_MAPPER (finite, ComplexNDArray::bmapper, xfinite) -ARRAY_MAPPER (isinf, ComplexNDArray::bmapper, xisinf) ARRAY_MAPPER (isna, ComplexNDArray::bmapper, octave_is_NA) -ARRAY_MAPPER (isnan, ComplexNDArray::bmapper, xisnan) /* ;;; Local Variables: ***
--- a/src/ov-flt-cx-mat.cc +++ b/src/ov-flt-cx-mat.cc @@ -342,20 +342,15 @@ { FloatComplexNDArray tmp(dv); - if (tmp.is_empty ()) - matrix = tmp; - else - { - is >> tmp; + is >> tmp; - if (is) - matrix = tmp; - else - { - error ("load: failed to load matrix constant"); - success = false; - } - } + if (is) + matrix = tmp; + else + { + error ("load: failed to load matrix constant"); + success = false; + } } else { @@ -823,6 +818,24 @@ return ::imag (matrix); } +octave_value +octave_float_complex_matrix::isnan (void) const +{ + return matrix.isnan (); +} + +octave_value +octave_float_complex_matrix::isinf (void) const +{ + return matrix.isinf (); +} + +octave_value +octave_float_complex_matrix::finite (void) const +{ + return matrix.isfinite (); +} + DARRAY_MAPPER (erf, FloatNDArray::dmapper, ::erff) DARRAY_MAPPER (erfc, FloatNDArray::dmapper, ::erfcf) DARRAY_MAPPER (gamma, FloatNDArray::dmapper, xgamma) @@ -855,10 +868,7 @@ ARRAY_MAPPER (sqrt, FloatComplexNDArray::cmapper, std::sqrt) ARRAY_MAPPER (tan, FloatComplexNDArray::cmapper, std::tan) ARRAY_MAPPER (tanh, FloatComplexNDArray::cmapper, std::tanh) -ARRAY_MAPPER (finite, FloatComplexNDArray::bmapper, xfinite) -ARRAY_MAPPER (isinf, FloatComplexNDArray::bmapper, xisinf) ARRAY_MAPPER (isna, FloatComplexNDArray::bmapper, octave_is_NA) -ARRAY_MAPPER (isnan, FloatComplexNDArray::bmapper, xisnan) /* ;;; Local Variables: ***
--- a/src/ov-flt-re-mat.cc +++ b/src/ov-flt-re-mat.cc @@ -372,20 +372,15 @@ { FloatNDArray tmp(dv); - if (tmp.is_empty ()) - matrix = tmp; - else - { - is >> tmp; + is >> tmp; - if (is) - matrix = tmp; - else - { - error ("load: failed to load matrix constant"); - success = false; - } - } + if (is) + matrix = tmp; + else + { + error ("load: failed to load matrix constant"); + success = false; + } } else { @@ -766,6 +761,24 @@ return FloatNDArray (matrix.dims (), 0.0); } +octave_value +octave_float_matrix::isnan (void) const +{ + return matrix.isnan (); +} + +octave_value +octave_float_matrix::isinf (void) const +{ + return matrix.isinf (); +} + +octave_value +octave_float_matrix::finite (void) const +{ + return matrix.isfinite (); +} + ARRAY_MAPPER (erf, FloatNDArray::dmapper, ::erff) ARRAY_MAPPER (erfc, FloatNDArray::dmapper, ::erfcf) ARRAY_MAPPER (gamma, FloatNDArray::dmapper, xgamma) @@ -797,10 +810,7 @@ CD_ARRAY_MAPPER (sqrt, ::sqrtf, std::sqrt, 0.0, octave_Float_Inf) ARRAY_MAPPER (tan, FloatNDArray::dmapper, ::tanf) ARRAY_MAPPER (tanh, FloatNDArray::dmapper, ::tanhf) -ARRAY_MAPPER (finite, FloatNDArray::bmapper, xfinite) -ARRAY_MAPPER (isinf, FloatNDArray::bmapper, xisinf) ARRAY_MAPPER (isna, FloatNDArray::bmapper, octave_is_NA) -ARRAY_MAPPER (isnan, FloatNDArray::bmapper, xisnan) DEFUN (single, args, , "-*- texinfo -*-\n\
--- a/src/ov-re-mat.cc +++ b/src/ov-re-mat.cc @@ -380,20 +380,15 @@ { NDArray tmp(dv); - if (tmp.is_empty ()) - matrix = tmp; - else - { - is >> tmp; + is >> tmp; - if (is) - matrix = tmp; - else - { - error ("load: failed to load matrix constant"); - success = false; - } - } + if (is) + matrix = tmp; + else + { + error ("load: failed to load matrix constant"); + success = false; + } } else { @@ -794,6 +789,24 @@ return NDArray (matrix.dims (), 0.0); } +octave_value +octave_matrix::isnan (void) const +{ + return matrix.isnan (); +} + +octave_value +octave_matrix::isinf (void) const +{ + return matrix.isinf (); +} + +octave_value +octave_matrix::finite (void) const +{ + return matrix.isfinite (); +} + ARRAY_MAPPER (erf, NDArray::dmapper, ::erf) ARRAY_MAPPER (erfc, NDArray::dmapper, ::erfc) ARRAY_MAPPER (gamma, NDArray::dmapper, xgamma) @@ -825,10 +838,7 @@ CD_ARRAY_MAPPER (sqrt, ::sqrt, std::sqrt, 0.0, octave_Inf) ARRAY_MAPPER (tan, NDArray::dmapper, ::tan) ARRAY_MAPPER (tanh, NDArray::dmapper, ::tanh) -ARRAY_MAPPER (finite, NDArray::bmapper, xfinite) -ARRAY_MAPPER (isinf, NDArray::bmapper, xisinf) ARRAY_MAPPER (isna, NDArray::bmapper, octave_is_NA) -ARRAY_MAPPER (isnan, NDArray::bmapper, xisnan) DEFUN (double, args, , "-*- texinfo -*-\n\
--- a/src/ov.cc +++ b/src/ov.cc @@ -474,11 +474,6 @@ return retval; } -octave_value::octave_value (void) - : rep (new octave_base_value ()) -{ -} - octave_value::octave_value (short int i) : rep (new octave_scalar (i)) {
--- a/src/ov.h +++ b/src/ov.h @@ -156,7 +156,13 @@ enum magic_colon { magic_colon_t }; - octave_value (void); + octave_value (void) + { + static octave_base_value nil_rep; + rep = &nil_rep; + rep->count++; + } + octave_value (short int i); octave_value (unsigned short int i); octave_value (int i); @@ -824,6 +830,13 @@ string_vector map_keys (void) const { return rep->map_keys (); } + string_vector parent_class_names (void) const + { return rep->parent_class_names (); } + + octave_base_value * + find_parent_class (const std::string& parent_class_name) + { return rep->find_parent_class (parent_class_name); } + octave_function *function_value (bool silent = false); const octave_function *function_value (bool silent = false) const;
--- a/src/parse.y +++ b/src/parse.y @@ -2528,7 +2528,10 @@ curr_fcn_ptr = script; - symbol_table::clear_forced_variables (); + // Unmark any symbols that may have been tagged as local variables + // while parsing (for example, by force_local_variable in lex.l). + + symbol_table::unmark_forced_variables (); } // Begin defining a function. @@ -2708,10 +2711,11 @@ retval = new tree_function_def (fcn); } - // Clear any local variables that may have been added while - // parsing (for example, by force_local_variable in lex.l). - - symbol_table::clear_forced_variables (fcn->scope ()); + // Unmark any symbols that may have been tagged as local + // variables while parsing (for example, by force_local_variable + // in lex.l). + + symbol_table::unmark_forced_variables (fcn->scope ()); } return retval;
--- a/src/pt-eval.cc +++ b/src/pt-eval.cc @@ -339,8 +339,8 @@ DO_SIMPLE_FOR_LOOP_ONCE (rhs); } - else if (rhs.is_matrix_type () - || rhs.is_cell () || rhs.is_string ()) + else if (rhs.is_matrix_type () || rhs.is_cell () || rhs.is_string () + || rhs.is_map ()) { // A matrix or cell is reshaped to 2 dimensions and iterated by // columns. @@ -384,27 +384,6 @@ } } } - else if (rhs.is_map ()) - { - Octave_map tmp_val (rhs.map_value ()); - - bool quit = false; - - for (Octave_map::iterator p = tmp_val.begin (); - p != tmp_val.end (); - p++) - { - Cell val_lst = tmp_val.contents (p); - - octave_value val - = (val_lst.length () == 1) ? val_lst(0) : octave_value (val_lst); - - DO_SIMPLE_FOR_LOOP_ONCE (val); - - if (quit) - break; - } - } else { ::error ("invalid type in for loop expression near line %d, column %d", @@ -456,15 +435,15 @@ octave_lvalue key_ref = elt->lvalue (); - Octave_map tmp_val (rhs.map_value ()); + const Octave_map tmp_val (rhs.map_value ()); tree_statement_list *loop_body = cmd.body (); - for (Octave_map::iterator q = tmp_val.begin (); q != tmp_val.end (); q++) + for (Octave_map::const_iterator q = tmp_val.begin (); q != tmp_val.end (); q++) { octave_value key = tmp_val.key (q); - Cell val_lst = tmp_val.contents (q); + const Cell val_lst = tmp_val.contents (q); octave_idx_type n = tmp_val.numel ();
--- a/src/sparse-xdiv.cc +++ b/src/sparse-xdiv.cc @@ -31,6 +31,7 @@ #include "oct-cmplx.h" #include "quit.h" #include "error.h" +#include "lo-ieee.h" #include "dSparse.h" #include "dDiagMatrix.h"
--- a/src/symtab.h +++ b/src/symtab.h @@ -204,10 +204,7 @@ octave_value& val = varref (context); if (! val.is_defined ()) - { - val = Matrix (); - mark_forced (); - } + mark_forced (); } octave_value& varref (context_id context) @@ -293,15 +290,6 @@ } } - void clear_forced (void) - { - if (is_forced ()) - { - varref (xcurrent_context) = octave_value (); - unmark_forced (); - } - } - bool is_defined (context_id context) const { return varval (context).is_defined (); @@ -309,7 +297,7 @@ bool is_variable (context_id context) const { - return (storage_class != local || is_defined (context)); + return (storage_class != local || is_defined (context) || is_forced ()); } bool is_local (void) const { return storage_class & local; } @@ -457,8 +445,6 @@ void clear (void) { rep->clear (); } - void clear_forced (void) { rep->clear_forced (); } - bool is_defined (context_id context = xcurrent_context) const { return rep->is_defined (context); @@ -1353,12 +1339,12 @@ inst->do_clear_variables (); } - static void clear_forced_variables (scope_id scope = xcurrent_scope) + static void unmark_forced_variables (scope_id scope = xcurrent_scope) { symbol_table *inst = get_instance (scope); if (inst) - inst->do_clear_forced_variables (); + inst->do_unmark_forced_variables (); } // For unwind_protect. @@ -2133,10 +2119,10 @@ p->second.clear (); } - void do_clear_forced_variables (void) + void do_unmark_forced_variables (void) { for (table_iterator p = table.begin (); p != table.end (); p++) - p->second.clear_forced (); + p->second.unmark_forced (); } void do_clear_global (const std::string& name)