Mercurial > hg > octave-max

changeset 13477:ca767e4055c5

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Trivial merge with Savannah
author Jordi Gutiérrez Hermoso <jordigh@gmail.com>
date Wed, 20 Apr 2011 09:13:54 -0500
parents bf51c1bb7033 (current diff) 57e5c31c28d5 (diff)
children a51ad78b51c5
files ChangeLog ChangeLog.1 doc/ChangeLog libcruft/ChangeLog liboctave/ChangeLog scripts/ChangeLog src/ChangeLog test/ChangeLog
diffstat 38 files changed, 753 insertions(+), 394 deletions(-) [+]
line wrap: on
line diff
--- a/Makefile.am
+++ b/Makefile.am
@@ -26,20 +26,27 @@
 
 ACLOCAL_AMFLAGS = -I m4
 
-BUILT_DISTFILES = AUTHORS BUGS INSTALL.OCTAVE
+BUILT_DISTFILES = AUTHORS BUGS ChangeLog INSTALL.OCTAVE
 
 EXTRA_DIST = \
   AUTHORS \
   BUGS \
   COPYING \
   ChangeLog \
-  ChangeLog.1 \
   INSTALL \
   INSTALL.OCTAVE \
   NEWS \
   NEWS.1 \
   NEWS.2 \
   NEWS.3 \
+  OLD-ChangeLogs/ChangeLog \
+  OLD-ChangeLogs/ChangeLog.1 \
+  OLD-ChangeLogs/doc-ChangeLog \
+  OLD-ChangeLogs/libcruft-ChangeLog \
+  OLD-ChangeLogs/liboctave-ChangeLog \
+  OLD-ChangeLogs/scripts-ChangeLog \
+  OLD-ChangeLogs/src-ChangeLog \
+  OLD-ChangeLogs/test-ChangeLog \
   PROJECTS \
   README \
   README.Cygwin \
@@ -156,6 +163,11 @@
 	$(MAKE) -C doc/interpreter ../../$@
 .PHONY: AUTHORS BUGS INSTALL.OCTAVE
 
+ChangeLog:
+	(cd $(srcdir); hg log --style=changelog.tmpl --prune=b0e60ad4ae26ec2ca3567a17b29a53e0cd2051d8 --branch=`hg branch`; echo ""; echo "See the files in the OLD-ChangeLogs directory for older changes") > $@.t
+	mv $@.t $@
+.PHONY: ChangeLog
+
 octetc_DATA = NEWS
 
 DIRS_TO_MAKE = \
--- a/NEWS
+++ b/NEWS
@@ -352,16 +352,15 @@
     exponent where the exponent is a multiple of 3.
 
  ** The following functions are new in Octave 3.4:
-
-      accumdim    erfcx        nfields      pqpnonneg  uigetdir 
-      bitpack     fileread     nth_element  quadcc     uigetfile  
-      bitunpack   fminbnd      onCleanup    randi      uiputfile    
-      blkmm       fskipl       pbaspect     repelems   uimenu  
-      cbrt        ifelse       pie3         reset      whitebg
-      curl        ishermitian  powerset     rsf2csf    
-      chop        isindex      ppder        saveas          
-      daspect     luupdate     ppint        strread          
-      divergence  merge        ppjumps      textread 
+      accumdim   divergence   merge        ppjumps    textread
+      bitpack    erfcx        nfields      pqpnonneg  uigetdir
+      bitunpack  fileread     nth_element  quadcc     uigetfile
+      blkmm      fminbnd      onCleanup    randi      uiputfile
+      cbrt       fskipl       pbaspect     repelems   uimenu
+      curl       ifelse       pie3         reset      whitebg
+      chop       ishermitian  powerset     rsf2csf 
+      colstyle   isindex      ppder        saveas  
+      daspect    luupdate     ppint        strread 
 
  ** Using the image function to view images with external programs such
     as display, xv, and xloadimage is no longer supported.  The
rename from ChangeLog
rename to OLD-ChangeLogs/ChangeLog
--- a/ChangeLog
+++ b/OLD-ChangeLogs/ChangeLog
@@ -1,3 +1,7 @@
+2011-04-14  Rik  <octave@nomad.inbox5.com>
+
+	* NEWS: Add colstyle to list of new functions for 3.4
+
 2011-04-08  Rik  <octave@nomad.inbox5.com>
 
 	* NEWS: Deprecate studentize(), add new function zscore().
rename from ChangeLog.1
rename to OLD-ChangeLogs/ChangeLog.1
rename from doc/ChangeLog
rename to OLD-ChangeLogs/doc-ChangeLog
--- a/doc/ChangeLog
+++ b/OLD-ChangeLogs/doc-ChangeLog
@@ -2,6 +2,13 @@
 
 	* interpreter/plot.txi: Clarify that inheritance of default property
 	values only applies to the named object type.
+2011-04-14  Rik  <octave@nomad.inbox5.com>
+
+	* interpreter/plot.txi: Add colstyle function to documentation.
+
+2011-04-12  Rik  <octave@nomad.inbox5.com>
+
+	* interpreter/expr.txi: Correct use of it's -> its in documentation.
 
 2011-04-12  Ben Abbott  <bpabbott@mac.com>
 
@@ -44,7 +51,7 @@
 
 2011-04-04  Rik  <octave@nomad.inbox5.com>
 
-	* interpreter/doccheck/aspell-octave.en.pws, interpreter/nonlin.txi, 
+	* interpreter/doccheck/aspell-octave.en.pws, interpreter/nonlin.txi,
 	interpreter/tips.txi: Spellcheck documentation for 3.4.1 release.
 
 2011-04-04  Rik  <octave@nomad.inbox5.com>
rename from libcruft/ChangeLog
rename to OLD-ChangeLogs/libcruft-ChangeLog
rename from liboctave/ChangeLog
rename to OLD-ChangeLogs/liboctave-ChangeLog
rename from scripts/ChangeLog
rename to OLD-ChangeLogs/scripts-ChangeLog
--- a/scripts/ChangeLog
+++ b/OLD-ChangeLogs/scripts-ChangeLog
@@ -1,3 +1,7 @@
+2011-04-18  Paul Boven  <p.boven@xs4all.nl>
+
+	* image/image.m: Fixed naming of variables in texinfo
+
 2011-04-17  Patrick Häcker  <magicmuscleman>
 
 	* strings/mat2str.m: Limit the number of digits to one less than
@@ -13,6 +17,11 @@
 	Handle nD-arguments correctly. Tests added.
 	(bugs #32040, #32045)
 
+2011-04-13  David Bateman  <dbateman@free.fr>
+
+	* plot/colstyle.m : New function.
+	* plot/module.mk plot_FCN_FILES) : Add it here.
+
 2011-04-13  Rik  <octave@nomad.inbox5.com>
 
 	* help/__makeinfo__.m: Simplify function by using regular expressions.
rename from src/ChangeLog
rename to OLD-ChangeLogs/src-ChangeLog
--- a/src/ChangeLog
+++ b/OLD-ChangeLogs/src-ChangeLog
@@ -1,3 +1,17 @@
+2011-04-19  Kai Habel  <kai.habel@gmx.de>
+
+	* src/DLD-FUNCTIONS/__init_fltk__.cc(plot_window::plot_window):
+	Instantiate canvas before uimenu.
+
+2011-04-13  Rik  <octave@nomad.inbox5.com>
+
+	* help.cc: Add spaces after commas in @seealso blocks.
+
+2011-04-12  Rik  <octave@nomad.inbox5.com>
+
+	* load-path.cc (restoredefaultpath): Correct use of it's -> its in
+	documentation.
+
 2011-04-10  John Eaton  <jwe@octave.org>
 
 	* graphics.cc (Fishandle) Accept vector of handles (bug #33025).
@@ -81,7 +95,7 @@
 
 	* DLD-FUNCTIONS/inv.cc (inv, inverse), DLD-FUNCTIONS/tril.cc (tril),
 	data.cc (cumsum, szie), file-io.cc (fgets), ov-typeinfo.cc (typeinfo),
-	ov-usr-fcn.cc (nargout), utils.cc (make_absolute_filename), 
+	ov-usr-fcn.cc (nargout), utils.cc (make_absolute_filename),
 	variables.cc (who): Improve docstrings
 
 2011-03-25  John W. Eaton  <jwe@octave.org>
rename from test/ChangeLog
rename to OLD-ChangeLogs/test-ChangeLog
--- a/autogen.sh
+++ b/autogen.sh
@@ -15,7 +15,7 @@
 ## building the rest of Octave, and INSTALL, which is linked from
 ## gnulib/doc/INSTALL by the bootstrap script.
 
-for f in NEWS README ChangeLog COPYING; do
+for f in NEWS README COPYING; do
   if ! test -f $f; then
     echo "required file $f is missing" 2>&1
     exit 1
new file mode 100644
--- /dev/null
+++ b/changelog.tmpl
@@ -0,0 +1,2 @@
+header = '{date|shortdate}  {author|person}  <{author|email}>\n\n'
+changeset = '\t{desc|tabindent|strip}\n\n\tFiles: {files|stringify|fill68|tabindent|strip}\n\n'
--- a/doc/interpreter/doccheck/aspell-octave.en.pws
+++ b/doc/interpreter/doccheck/aspell-octave.en.pws
@@ -179,6 +179,7 @@
 delaunayn
 DeleteFcn
 delim
+deltaX
 demi
 Demmel
 DeskJet
@@ -575,6 +576,7 @@
 nocompute
 nolabel
 noncommercially
+nonsmooth
 nonzeros
 noperm
 normcdf
@@ -674,7 +676,11 @@
 QQ
 QRUPDATE
 qrupdate
+quadcc
+quadgk
+quadl
 quadpack
+quadv
 Quantile
 quantile
 quantiles
@@ -806,6 +812,8 @@
 Subfunctions
 subfunctions
 subinterval
+Subintervals
+subintervals
 sublicenses
 Sublicensing
 submatrices
@@ -866,6 +874,7 @@
 tp
 tpdf
 traceback
+trapz
 treelayout
 treeplot
 tridiagonal
--- a/doc/interpreter/expr.txi
+++ b/doc/interpreter/expr.txi
@@ -912,7 +912,7 @@
 @sc{matlab} has special behavior that allows the operators @samp{&} and
 @samp{|} to short-circuit when used in the truth expression for @code{if} and 
 @code{while} statements.  The Octave parser may be instructed to behave in the
-same manner, but it's use is strongly discouraged.
+same manner, but its use is strongly discouraged.
 
 @DOCSTRING(do_braindead_shortcircuit_evaluation)
 
--- a/doc/interpreter/octave.texi
+++ b/doc/interpreter/octave.texi
@@ -657,8 +657,8 @@
 Numerical Integration
 
 * Functions of One Variable:: 
+* Orthogonal Collocation::      
 * Functions of Multiple Variables:: 
-* Orthogonal Collocation::      
 
 Differential Equations
 
--- a/doc/interpreter/plot.txi
+++ b/doc/interpreter/plot.txi
@@ -2553,6 +2553,11 @@
 of 2 is twice as large as the default, etc.
 @end table
 
+The @code{colstyle} function will parse a @code{plot}-style specification
+and will return the color, line, and marker values that would result.
+
+@DOCSTRING(colstyle)
+
 @node Callbacks
 @subsection Callbacks
 @cindex callbacks
--- a/doc/interpreter/quad.txi
+++ b/doc/interpreter/quad.txi
@@ -20,19 +20,19 @@
 @chapter Numerical Integration
 
 Octave comes with several built-in functions for computing the integral
-of a function numerically.  These functions all solve 1-dimensional
-integration problems.
+of a function numerically (termed quadrature).  These functions all solve
+1-dimensional integration problems.
 
 @menu
 * Functions of One Variable:: 
+* Orthogonal Collocation::      
 * Functions of Multiple Variables:: 
-* Orthogonal Collocation::      
 @end menu
 
 @node Functions of One Variable
 @section Functions of One Variable
 
-Octave supports three different algorithms for computing the integral
+Octave supports five different algorithms for computing the integral
 @tex
 $$
  \int_a^b f(x) d x
@@ -45,30 +45,42 @@
 @item quad
 Numerical integration based on Gaussian quadrature.
 
+@item quadv
+Numerical integration using an adaptive vectorized Simpson's rule.
+
 @item quadl
 Numerical integration using an adaptive Lobatto rule.
 
 @item quadgk
 Numerical integration using an adaptive Gauss-Konrod rule.
 
-@item quadv
-Numerical integration using an adaptive vectorized Simpson's rule.
-
 @item quadcc
 Numerical integration using adaptive Clenshaw-Curtis rules.
 
-@item trapz
-Numerical integration using the trapezoidal method.
+@item trapz, cumtrapz
+Numerical integration of data using the trapezoidal method.
 @end table
 
 @noindent
-Besides these functions Octave also allows you to perform cumulative
-numerical integration using the trapezoidal method through the
-@code{cumtrapz} function.
+The best quadrature algorithm to use depends on the integrand.  If you have
+empirical data, rather than a function, the choice is @code{trapz} or
+@code{cumtrapz}.  If you are uncertain about the characteristics of the
+integrand, @code{quadcc} will be the most robust as it can handle
+discontinuities, singularities, oscillatory functions, and infinite intervals.
+When the integrand is smooth @code{quadgk} may be the fastest of the
+algorithms.
 
-@DOCSTRING(quad)
+@multitable @columnfractions 0.05 0.15 .80
+@headitem @tab Function @tab Characteristics
+@item @tab quad   @tab Low accuracy with nonsmooth integrands
+@item @tab quadv  @tab Medium accuracy with smooth integrands
+@item @tab quadl  @tab Medium accuracy with smooth integrands.  Slower than quadgk.
+@item @tab quadgk @tab Medium accuracy (@math{1e^{-6}}--@math{1e^{-9}}) with smooth integrands.
+@item @tab        @tab Handles oscillatory functions and infinite bounds
+@item @tab quadcc @tab Low to High accuracy with nonsmooth/smooth integrands
+@item @tab        @tab Handles oscillatory functions, singularities, and infinite bounds
+@end multitable
 
-@DOCSTRING(quad_options)
 
 Here is an example of using @code{quad} to integrate the function
 @tex
@@ -95,21 +107,22 @@
 @example
 @group
 function y = f (x)
-  y = x .* sin (1 ./ x) .* sqrt (abs (1 - x));
+  y = x .* sin (1./x) .* sqrt (abs (1 - x));
 endfunction
 @end group
 @end example
 
-Note the use of the `dot' forms of the operators.  This is not necessary
-for the call to @code{quad}, but it makes it much easier to generate a
-set of points for plotting (because it makes it possible to call the
-function with a vector argument to produce a vector result).
+Note the use of the `dot' forms of the operators.  This is not necessary for
+the @code{quad} integrator, but is required by the other integrators.  In any
+case, it makes it much easier to generate a set of points for plotting because
+it is possible to call the function with a vector argument to produce a vector
+result.
 
-Then we simply call quad:
+The second step is to call quad with the limits of integration:
 
 @example
 @group
-[v, ier, nfun, err] = quad ("f", 0, 3)
+[q, ier, nfun, err] = quad ("f", 0, 3)
      @result{} 1.9819
      @result{} 1
      @result{} 5061
@@ -121,13 +134,83 @@
 is reasonably accurate (to see why, examine what happens to the result
 if you move the lower bound to 0.1, then 0.01, then 0.001, etc.).
 
+The function "f" can be the string name of a function, a function handle, or
+an inline function.  These options make it quite easy to do integration
+without having to fully define a function in an m-file.  For example:
+
+@example
+@group
+# Verify integral (x^3) = x^4/4
+f = inline ("x.^3");
+quadgk (f, 0, 1)
+     @result{} 0.25000
+
+# Verify gamma function = (n-1)! for n = 4
+f = @@(x) x.^3 .* exp (-x);
+quadcc (f, 0, Inf)
+     @result{} 6.0000
+@end group
+@end example
+
+@DOCSTRING(quad)
+
+@DOCSTRING(quad_options)
+
+@DOCSTRING(quadv)
+
 @DOCSTRING(quadl)
 
 @DOCSTRING(quadgk)
 
-@DOCSTRING(quadv)
+@DOCSTRING(quadcc)
+
+Sometimes one does not have the function, but only the raw (x, y) points from
+which to perform an integration.  This can occur when collecting data in an
+experiment.  The @code{trapz} function can integrate these values as shown in
+the following example where "data" has been collected on the cosine function
+over the range [0, pi/2). 
+
+@example
+@group
+x = 0:0.1:pi/2;  # Uniformly spaced points
+y = cos (x);
+trapz (x, y)
+     @result{} 0.99666
+@end group
+@end example
+
+The answer is reasonably close to the exact value of 1.  Ordinary quadrature
+is sensitive to the characteristics of the integrand.  Empirical integration 
+depends not just on the integrand, but also on the particular points chosen to
+represent the function.  Repeating the example above with the sine function
+over the range [0, pi/2) produces far inferior results.
 
-@DOCSTRING(quadcc)
+@example
+@group
+x = 0:0.1:pi/2;  # Uniformly spaced points
+y = sin (x);
+trapz (x, y)
+     @result{} 0.92849
+@end group
+@end example
+
+However, a slightly different choice of data points can change the result
+significantly.  The same integration, with the same number of points, but
+spaced differently produces a more accurate answer.
+
+@example
+@group
+x = linspace (0, pi/2, 16);  # Uniformly spaced, but including endpoint
+y = sin (x);
+trapz (x, y)
+     @result{} 0.99909
+@end group
+@end example
+
+In general there may be no way of knowing the best distribution of points ahead
+of time.  Or the points may come from an experiment where there is no freedom to
+select the best distribution.  In any case, one must remain aware of this
+issue when using @code{trapz}.
 
 @DOCSTRING(trapz)
 
@@ -183,9 +266,9 @@
 @section Functions of Multiple Variables
 
 Octave does not have built-in functions for computing the integral of
-functions of multiple variables directly.  It is however possible to
+functions of multiple variables directly.  It is possible, however, to
 compute the integral of a function of multiple variables using the
-functions for one-dimensional integrals.
+existing functions for one-dimensional integrals.
 
 To illustrate how the integration can be performed, we will integrate
 the function
@@ -205,23 +288,23 @@
 
 The first approach creates a function that integrates @math{f} with
 respect to @math{x}, and then integrates that function with respect to
-@math{y}.  Since @code{quad} is written in Fortran it cannot be called
+@math{y}.  Because @code{quad} is written in Fortran it cannot be called
 recursively.  This means that @code{quad} cannot integrate a function
 that calls @code{quad}, and hence cannot be used to perform the double
-integration.  It is however possible with @code{quadl}, which is what
-the following code does.
+integration.  Any of the other integrators, however, can be used which is
+what the following code demonstrates.
 
 @example
 @group
-function I = g(y)
-  I = ones(1, length(y));
-  for i = 1:length(y)
-    f = @@(x) sin(pi.*x.*y(i)).*sqrt(x.*y(i));
-    I(i) = quadl(f, 0, 1);
+function q = g(y)
+  q = ones (size (y));
+  for i = 1:length (y)
+    f = @@(x) sin (pi*x.*y(i)) .* sqrt (x.*y(i));
+    q(i) = quadgk (f, 0, 1);
   endfor
 endfunction
 
-I = quadl("g", 0, 1)
+I = quadgk ("g", 0, 1)
       @result{} 0.30022
 @end group
 @end example
@@ -232,7 +315,7 @@
 
 @example
 @group
-I =  dblquad (@@(x, y) sin(pi.*x.*y).*sqrt(x.*y), 0, 1, 0, 1)
+I = dblquad (@@(x, y) sin (pi*x.*y) .* sqrt (x.*y), 0, 1, 0, 1)
       @result{} 0.30022
 @end group
 @end example
@@ -241,12 +324,12 @@
 
 @DOCSTRING(triplequad)
 
-The above mentioned approach works but is fairly slow, and that problem
-increases exponentially with the dimensionality the problem.  Another
+The above mentioned approach works, but is fairly slow, and that problem
+increases exponentially with the dimensionality of the integral.  Another
 possible solution is to use Orthogonal Collocation as described in the
-previous section.  The integral of a function @math{f(x,y)} for
-@math{x} and @math{y} between 0 and 1 can be approximated using @math{n}
-points by
+previous section (@pxref{Orthogonal Collocation}).  The integral of a function
+@math{f(x,y)} for @math{x} and @math{y} between 0 and 1 can be approximated
+using @math{n} points by
 @tex
 $$
  \int_0^1 \int_0^1 f(x,y) d x d y \approx \sum_{i=1}^n \sum_{j=1}^n q_i q_j f(r_i, r_j),
@@ -257,13 +340,13 @@
 @end ifnottex
 where @math{q} and @math{r} is as returned by @code{colloc(n)}.  The
 generalization to more than two variables is straight forward.  The
-following code computes the studied integral using @math{n=7} points.
+following code computes the studied integral using @math{n=8} points.
 
 @example
 @group
-f = @@(x,y) sin(pi*x*y').*sqrt(x*y');
-n = 7;
-[t, A, B, q] = colloc(n);
+f = @@(x,y) sin (pi*x*y') .* sqrt (x*y');
+n = 8;
+[t, ~, ~, q] = colloc (n);
 I = q'*f(t,t)*q;
       @result{} 0.30022
 @end group
@@ -272,7 +355,5 @@
 @noindent
 It should be noted that the number of points determines the quality
 of the approximation.  If the integration needs to be performed between
-@math{a} and @math{b} instead of 0 and 1, a change of variables is needed.
+@math{a} and @math{b}, instead of 0 and 1, then a change of variables is needed.
 
-
-
--- a/liboctave/Quad-opts.in
+++ b/liboctave/Quad-opts.in
@@ -34,7 +34,7 @@
   DOC_ITEM
 Non-negative relative tolerance.  If the absolute tolerance is zero,
 the relative tolerance must be greater than or equal to 
-@code{max (50*eps, 0.5e-28)}.
+@w{@code{max (50*eps, 0.5e-28)}}.
 
   END_DOC_ITEM
   TYPE = "double"
@@ -59,7 +59,7 @@
   DOC_ITEM
 Non-negative relative tolerance for single precision.  If the absolute
 tolerance is zero, the relative tolerance must be greater than or equal to 
-@code{max (50*eps, 0.5e-28)}.
+@w{@code{max (50*eps, 0.5e-28)}}.
   END_DOC_ITEM
   TYPE = "float"
   INIT_VALUE = "::sqrt (FLT_EPSILON)"
--- a/mk-opts.pl
+++ b/mk-opts.pl
@@ -126,6 +126,7 @@
           die "mk-opts.pl: unknown command: $_\n"
         }
     }
+  $INCLUDE = "" if not defined $INCLUDE;   # Initialize value if required
 }
 
 sub parse_option_block
@@ -230,11 +231,12 @@
 
   if (not defined $DOC_STRING)
     {
-      $DOC_STRING = "When called with two arguments, this function\\n\\
-allows you set options parameters for the function \@code{$FCN_NAME}.\\n\\
-Given one argument, \@code{$OPT_FCN_NAME} returns the value of the\\n\\
-corresponding option.  If no arguments are supplied, the names of all\\n\\
-the available options and their current values are displayed.\\n\\\n";
+      $DOC_STRING = "Query or set options for the function \@code{$FCN_NAME}.\\n\\
+When called with no arguments, the names of all available options and\\n\\
+their current values are displayed.\\n\\
+Given one argument, return the value of the corresponding option.\\n\\
+When called with two arguments, \@code{$OPT_FCN_NAME} set the option\\n\\
+\@var{opt} to value \@var{val}.";
     }
 }
 
@@ -909,8 +911,11 @@
   print <<"_END_EMIT_OPTIONS_FUNCTION_HDR_";
 DEFUN_DLD ($OPT_FCN_NAME, args, ,
   "-*- texinfo -*-\\n\\
-\@deftypefn {Loadable Function} {} $OPT_FCN_NAME (\@var{opt}, \@var{val})\\n\\
+\@deftypefn  {Loadable Function} {} $OPT_FCN_NAME ()\\n\\
+\@deftypefnx {Loadable Function} {val =} $OPT_FCN_NAME (\@var{opt})\\n\\
+\@deftypefnx {Loadable Function} {} $OPT_FCN_NAME (\@var{opt}, \@var{val})\\n\\
 $DOC_STRING\\n\\
+\\n\\
 Options include\\n\\
 \\n\\
 \@table \@code\\n\\
--- a/scripts/general/cumtrapz.m
+++ b/scripts/general/cumtrapz.m
@@ -17,17 +17,25 @@
 ## <http://www.gnu.org/licenses/>.
 
 ## -*- texinfo -*-
-## @deftypefn  {Function File} {@var{z} =} cumtrapz (@var{y})
-## @deftypefnx {Function File} {@var{z} =} cumtrapz (@var{x}, @var{y})
-## @deftypefnx {Function File} {@var{z} =} cumtrapz (@dots{}, @var{dim})
+## @deftypefn  {Function File} {@var{q} =} cumtrapz (@var{y})
+## @deftypefnx {Function File} {@var{q} =} cumtrapz (@var{x}, @var{y})
+## @deftypefnx {Function File} {@var{q} =} cumtrapz (@dots{}, @var{dim})
 ##
-## Cumulative numerical integration using trapezoidal method.
-## @code{cumtrapz (@var{y})} computes the cumulative integral of the
-## @var{y} along the first non-singleton dimension.  If the argument
-## @var{x} is omitted an equally spaced vector is assumed.  @code{cumtrapz
-## (@var{x}, @var{y})} evaluates the cumulative integral with respect
-## to @var{x}.
+## Cumulative numerical integration of points @var{y} using the trapezoidal
+## method.
+## @w{@code{cumtrapz (@var{y})}} computes the cumulative integral of @var{y}
+## along the first non-singleton dimension.  Where @code{trapz} reports
+## only the overall integral sum, @code{cumtrapz} reports the current partial
+## sum value at each point of @var{y}.  When the argument @var{x} is omitted
+## an equally spaced @var{x} vector with unit spacing (1) is assumed. 
+## @code{cumtrapz (@var{x}, @var{y})} evaluates the integral with respect to
+## the spacing in @var{x} and the values in @var{y}.  This is useful if the
+## points in @var{y} have been sampled unevenly.  If the optional @var{dim}
+## argument is given, operate along this dimension.
 ##
+## If @var{x} is not specified then unit spacing will be used.  To scale
+## the integral to the correct value you must multiply by the actual spacing
+## value (deltaX).
 ## @seealso{trapz, cumsum}
 ## @end deftypefn
 
--- a/scripts/general/dblquad.m
+++ b/scripts/general/dblquad.m
@@ -17,30 +17,43 @@
 ## <http://www.gnu.org/licenses/>.
 
 ## -*- texinfo -*-
-## @deftypefn {Function File} {} dblquad (@var{f}, @var{xa}, @var{xb}, @var{ya}, @var{yb}, @var{tol}, @var{quadf}, @dots{})
-## Numerically evaluate a double integral.  The function over with to
-## integrate is defined by @code{@var{f}}, and the interval for the
-## integration is defined by @code{[@var{xa}, @var{xb}, @var{ya},
-## @var{yb}]}.  The function @var{f} must accept a vector @var{x} and a
-## scalar @var{y}, and return a vector of the same length as @var{x}.
+## @deftypefn  {Function File} {} dblquad (@var{f}, @var{xa}, @var{xb}, @var{ya}, @var{yb})
+## @deftypefnx {Function File} {} dblquad (@var{f}, @var{xa}, @var{xb}, @var{ya}, @var{yb}, @var{tol})
+## @deftypefnx {Function File} {} dblquad (@var{f}, @var{xa}, @var{xb}, @var{ya}, @var{yb}, @var{tol}, @var{quadf})
+## @deftypefnx {Function File} {} dblquad (@var{f}, @var{xa}, @var{xb}, @var{ya}, @var{yb}, @var{tol}, @var{quadf}, @dots{})
+## Numerically evaluate the double integral of @var{f}.
+## @var{f} is a function handle, inline function, or string
+## containing the name of the function to evaluate.  The function @var{f} must
+## have the form @math{z = f(x,y)} where @var{x} is a vector and @var{y} is a
+## scalar.  It should return a vector of the same length and orientation as
+## @var{x}.
 ##
-## If defined, @var{tol} defines the absolute tolerance to which to
-## which to integrate each sub-integral.
+## @var{xa}, @var{ya} and @var{xb}, @var{yb} are the lower and upper limits of
+## integration for x and y respectively.  The underlying integrator determines
+## whether infinite bounds are accepted.
+##
+## The optional argument @var{tol} defines the absolute tolerance used to 
+## integrate each sub-integral.  The default value is @math{1e^{-6}}.
+##
+## The optional argument @var{quadf} specifies which underlying integrator
+## function to use.  Any choice but @code{quad} is available and the default
+## is @code{quadcc}.
 ##
 ## Additional arguments, are passed directly to @var{f}.  To use the default
-## value for @var{tol} one may pass an empty matrix.
+## value for @var{tol} or @var{quadf} one may pass ':' or an empty matrix ([]).
 ## @seealso{triplequad, quad, quadv, quadl, quadgk, quadcc, trapz}
 ## @end deftypefn
 
-function q = dblquad(f, xa, xb, ya, yb, tol, quadf, varargin)
+function q = dblquad (f, xa, xb, ya, yb, tol = 1e-6, quadf = @quadcc, varargin)
+
   if (nargin < 5)
     print_usage ();
   endif
-  if (nargin < 6 || isempty (tol))
+  if (isempty (tol))
     tol = 1e-6;
   endif
-  if (nargin < 7 || isempty (quadf))
-    quadf = @quadgk;
+  if (isempty (quadf))
+    quadf = @quadcc;
   endif
 
   inner = @__dblquad_inner__;
@@ -60,9 +73,10 @@
   endfor
 endfunction
 
-%% Nasty integrand to show quadgk off
+%% Nasty integrand to show quadcc off
 %!assert (dblquad (@(x,y) 1 ./ (x+y), 0, 1, 0, 1), 2*log(2), 1e-6)
 
-%!assert (dblquad (@(x,y) exp(-x.^2 - y.^2) , -1, 1, -1, 1, [],  @quadgk), pi * erf(1).^2, 1e-6)
-%!assert (dblquad (@(x,y) exp(-x.^2 - y.^2) , -1, 1, -1, 1, [],  @quadl), pi * erf(1).^2, 1e-6)
-%!assert (dblquad (@(x,y) exp(-x.^2 - y.^2) , -1, 1, -1, 1, [],  @quadv), pi * erf(1).^2, 1e-6)
+%!assert (dblquad (@(x,y) exp(-x.^2 - y.^2) , -1, 1, -1, 1, 1e-6, @quadgk), pi * erf(1).^2, 1e-6)
+%!assert (dblquad (@(x,y) exp(-x.^2 - y.^2) , -1, 1, -1, 1, 1e-6, @quadl), pi * erf(1).^2, 1e-6)
+%!assert (dblquad (@(x,y) exp(-x.^2 - y.^2) , -1, 1, -1, 1, 1e-6, @quadv), pi * erf(1).^2, 1e-6)
+
--- a/scripts/general/quadgk.m
+++ b/scripts/general/quadgk.m
@@ -17,25 +17,29 @@
 ## <http://www.gnu.org/licenses/>.
 
 ## -*- texinfo -*-
-## @deftypefn  {Function File} {} quadgk (@var{f}, @var{a}, @var{b}, @var{abstol}, @var{trace})
-## @deftypefnx {Function File} {} quadgk (@var{f}, @var{a}, @var{b}, @var{prop}, @var{val}, @dots{})
+## @deftypefn  {Function File} {@var{q} =} quadgk (@var{f}, @var{a}, @var{b})
+## @deftypefnx {Function File} {@var{q} =} quadgk (@var{f}, @var{a}, @var{b}, @var{abstol})
+## @deftypefnx {Function File} {@var{q} =} quadgk (@var{f}, @var{a}, @var{b}, @var{abstol}, @var{trace})
+## @deftypefnx {Function File} {@var{q} =} quadgk (@var{f}, @var{a}, @var{b}, @var{prop}, @var{val}, @dots{})
 ## @deftypefnx {Function File} {[@var{q}, @var{err}] =} quadgk (@dots{})
-## Numerically evaluate integral using adaptive Gauss-Konrod quadrature.
+##
+## Numerically evaluate the integral of @var{f} from @var{a} to @var{b}
+## using adaptive Gauss-Konrod quadrature.
+## @var{f} is a function handle, inline function, or string
+## containing the name of the function to evaluate.
 ## The formulation is based on a proposal by L.F. Shampine,
 ## @cite{"Vectorized adaptive quadrature in @sc{matlab}", Journal of
 ## Computational and Applied Mathematics, pp131-140, Vol 211, Issue 2,
 ## Feb 2008} where all function evaluations at an iteration are
-## calculated with a single call to @var{f}.  Therefore the function
-## @var{f} must be of the form @code{@var{f} (@var{x})} and accept
-## vector values of @var{x} and return a vector of the same length
-## representing the function evaluations at the given values of @var{x}.
-## The function @var{f} can be defined in terms of a function handle,
-## inline function or string.
+## calculated with a single call to @var{f}.  Therefore, the function
+## @var{f} must be vectorized and must accept a vector of input values @var{x}
+## and return an output vector representing the function evaluations at the
+## given values of @var{x}.
 ##
-## The bounds of the quadrature @code{[@var{a}, @var{b}]} can be finite
-## or infinite and contain weak end singularities.  Variable
-## transformation will be used to treat infinite intervals and weaken
-## the singularities.  For example:
+## @var{a} and @var{b} are the lower and upper limits of integration.  Either
+## or both limits may be infinite or contain weak end singularities.
+## Variable transformation will be used to treat any infinite intervals and
+## weaken the singularities.  For example:
 ##
 ## @example
 ## quadgk(@@(x) 1 ./ (sqrt (x) .* (x + 1)), 0, Inf)
@@ -46,53 +50,57 @@
 ## element-by-element operator @code{./} and all user functions to
 ## @code{quadgk} should do the same.
 ##
-## The absolute tolerance can be passed as a fourth argument in a manner
-## compatible with @code{quadv}.  Equally the user can request that
-## information on the convergence can be printed is the fifth argument
-## is logically true.
+## The optional argument @var{tol} defines the absolute tolerance used to stop
+## the integration procedure.  The default value is @math{1e^{-10}}.
+## 
+## The algorithm used by @code{quadgk} involves subdividing the
+## integration interval and evaluating each subinterval.
+## If @var{trace} is true then after computing each of these partial
+## integrals display: (1) the number of subintervals at this step, 
+## (2) the current estimate of the error @var{err}, (3) the current estimate
+## for the integral @var{q}.
 ##
-## Alternatively, certain properties of @code{quadgk} can be passed as
-## pairs @code{@var{prop}, @var{val}}.  Valid properties are
+## Alternatively, properties of @code{quadgk} can be passed to the function as
+## pairs @code{"@var{prop}", @var{val}}.  Valid properties are
 ##
 ## @table @code
 ## @item AbsTol
-## Defines the absolute error tolerance for the quadrature.  The default
+## Define the absolute error tolerance for the quadrature.  The default
 ## absolute tolerance is 1e-10.
 ##
 ## @item RelTol
-## Defines the relative error tolerance for the quadrature.  The default
+## Define the relative error tolerance for the quadrature.  The default
 ## relative tolerance is 1e-5.
 ##
 ## @item MaxIntervalCount
 ## @code{quadgk} initially subdivides the interval on which to perform
-## the quadrature into 10 intervals.  Sub-intervals that have an
-## unacceptable error are sub-divided and re-evaluated.  If the number of
-## sub-intervals exceeds at any point 650 sub-intervals then a poor
+## the quadrature into 10 intervals.  Subintervals that have an
+## unacceptable error are subdivided and re-evaluated.  If the number of
+## subintervals exceeds 650 subintervals at any point then a poor
 ## convergence is signaled and the current estimate of the integral is
 ## returned.  The property 'MaxIntervalCount' can be used to alter the
-## number of sub-intervals that can exist before exiting.
+## number of subintervals that can exist before exiting.
 ##
 ## @item WayPoints
-## If there exists discontinuities in the first derivative of the
-## function to integrate, then these can be flagged with the
-## @code{"WayPoints"} property.  This forces the ends of a sub-interval
-## to fall on the breakpoints of the function and can result in
+## Discontinuities in the first derivative of the function to integrate can be
+## flagged with the  @code{"WayPoints"} property.  This forces the ends of
+## a subinterval to fall on the breakpoints of the function and can result in
 ## significantly improved estimation of the error in the integral, faster
-## computation or both.  For example,
+## computation, or both.  For example,
 ##
 ## @example
-## quadgk (@@(x) abs (1 - x .^ 2), 0, 2, 'Waypoints', 1)
+## quadgk (@@(x) abs (1 - x.^2), 0, 2, 'Waypoints', 1)
 ## @end example
 ##
 ## @noindent
 ## signals the breakpoint in the integrand at @code{@var{x} = 1}.
 ##
 ## @item Trace
-## If logically true, then @code{quadgk} prints information on the
+## If logically true @code{quadgk} prints information on the
 ## convergence of the quadrature at each iteration.
-##@end table
+## @end table
 ##
-## If any of @var{a}, @var{b} or @var{waypoints} is complex, then the
+## If any of @var{a}, @var{b}, or @var{waypoints} is complex then the
 ## quadrature is treated as a contour integral along a piecewise
 ## continuous path defined by the above.  In this case the integral is
 ## assumed to have no edge singularities.  For example,
@@ -108,9 +116,10 @@
 ## integrates @code{log (z)} along the square defined by @code{[1+1i,
 ##  1-1i, -1-1i, -1+1i]}
 ##
-## If two output arguments are requested, then @var{err} returns the
-## approximate bounds on the error in the integral @code{abs (@var{q} -
-## @var{i})}, where @var{i} is the exact value of the integral.
+## The result of the integration is returned in @var{q}.  
+## @var{err} is an approximate bound on the error in the integral 
+## @code{abs (@var{q} - @var{I})}, where @var{I} is the exact value of the
+## integral.
 ##
 ## @seealso{quad, quadv, quadl, quadcc, trapz, dblquad, triplequad}
 ## @end deftypefn
@@ -281,7 +290,7 @@
            3 .* (b - a) ./ 4 .* (1 - t.^2);
     endif
 
-    ## Split interval into at least 10 sub-interval with a 15 point
+    ## Split interval into at least 10 subinterval with a 15 point
     ## Gauss-Kronrod rule giving a minimum of 150 function evaluations
     while (length (subs) < 11)
       subs = [subs' ; subs(1:end-1)' + diff(subs') ./ 2, NaN](:)(1 : end - 1);
@@ -294,7 +303,7 @@
       ## Singularity will cause divide by zero warnings
       warning ("off", "Octave:divide-by-zero");
 
-      ## Initial evaluation of the integrand on the sub-intervals
+      ## Initial evaluation of the integrand on the subintervals
       [q_subs, q_errs] = __quadgk_eval__ (f, subs);
       q0 = sum (q_subs);
       err0 = sum (q_errs);
@@ -307,8 +316,8 @@
 
       first = true;
       while (true)
-        ## Check for sub-intervals that are too small. Test must be
-        ## performed in untransformed sub-intervals. What is a good
+        ## Check for subintervals that are too small. Test must be
+        ## performed in untransformed subintervals. What is a good
         ## value for this test. Shampine suggests 100*eps
         if (any (abs (diff (trans (subs), [], 2) / h0) < 100 * myeps))
           q = q0;
@@ -333,7 +342,7 @@
           break;
         endif
 
-        ## Accept the sub-intervals that meet the convergence criteria
+        ## Accept the subintervals that meet the convergence criteria
         idx = find (abs (q_errs) < tol .* abs(diff (subs, [], 2)) ./ h);
         if (first)
           q = sum (q_subs (idx));
@@ -347,7 +356,7 @@
         endif
         subs(idx,:) = [];
 
-        ## If no remaining sub-intervals exit
+        ## If no remaining subintervals exit
         if (rows (subs) == 0)
           break;
         endif
@@ -356,12 +365,12 @@
           disp([rows(subs), err, q0]);
         endif
 
-        ## Split remaining sub-intervals in two
+        ## Split remaining subintervals in two
         mid = (subs(:,2) + subs(:,1)) ./ 2;
         subs = [subs(:,1), mid; mid, subs(:,2)];
 
-        ## If the maximum sub-interval count is met accept remaining
-        ## sub-interval and exit
+        ## If the maximum subinterval count is met accept remaining
+        ## subinterval and exit
         if (rows (subs) > maxint)
           warning ("quadgk: maximum interval count (%d) met", maxint);
           q += sum (q_subs);
@@ -369,7 +378,7 @@
           break;
         endif
 
-        ## Evaluation of the integrand on the remaining sub-intervals
+        ## Evaluation of the integrand on the remaining subintervals
         [q_subs, q_errs] = __quadgk_eval__ (f, subs);
       endwhile
 
--- a/scripts/general/quadl.m
+++ b/scripts/general/quadl.m
@@ -22,21 +22,28 @@
 ## @deftypefnx {Function File} {@var{q} =} quadl (@var{f}, @var{a}, @var{b}, @var{tol}, @var{trace})
 ## @deftypefnx {Function File} {@var{q} =} quadl (@var{f}, @var{a}, @var{b}, @var{tol}, @var{trace}, @var{p1}, @var{p2}, @dots{})
 ##
-## Numerically evaluate integral using adaptive Lobatto rule.
-## @code{quadl (@var{f}, @var{a}, @var{b})} approximates the integral of
-## @code{@var{f}(@var{x})} to machine precision.  @var{f} is either a
-## function handle, inline function or string containing the name of
-## the function to evaluate.  The function @var{f} must return a vector
-## of output values if given a vector of input values.
+## Numerically evaluate the integral of @var{f} from @var{a} to @var{b}
+## using an adaptive Lobatto rule.
+## @var{f} is a function handle, inline function, or string
+## containing the name of the function to evaluate.
+## The function @var{f} must be vectorized and return a vector of output values
+## if given a vector of input values.
+##
+## @var{a} and @var{b} are the lower and upper limits of integration.  Both
+## limits must be finite.
 ##
-## If defined, @var{tol} defines the relative tolerance to which to
-## which to integrate @code{@var{f}(@var{x})}.  While if @var{trace} is
-## defined, displays the left end point of the current interval, the
-## interval length, and the partial integral.
+## The optional argument @var{tol} defines the relative tolerance with which
+## to perform the integration.  The default value is @code{eps}.
 ##
-## Additional arguments @var{p1}, etc., are passed directly to @var{f}.
-## To use default values for @var{tol} and @var{trace}, one may pass
-## empty matrices.
+## The algorithm used by @code{quadl} involves recursively subdividing the
+## integration interval.
+## If @var{trace} is defined then for each subinterval display: (1) the left
+## end of the subinterval, (2) the length of the subinterval, (3) the
+## approximation of the integral over the subinterval. 
+##
+## Additional arguments @var{p1}, etc., are passed directly to the function
+## @var{f}.  To use default values for @var{tol} and @var{trace}, one may pass
+## empty matrices ([]).
 ##
 ## Reference: W. Gander and W. Gautschi, @cite{Adaptive Quadrature -
 ## Revisited}, BIT Vol. 40, No. 1, March 2000, pp. 84--101.
@@ -55,7 +62,7 @@
 ##   * replace global variable terminate2 with local function need_warning
 ##   * add paper ref to docs
 
-function Q = quadl (f, a, b, tol, trace, varargin)
+function q = quadl (f, a, b, tol, trace, varargin)
   need_warning (1);
   if (nargin < 4)
     tol = [];
@@ -128,7 +135,7 @@
   if (is == 0)
     is = b-a;
   endif
-  Q = adaptlobstp (f, a, b, fa, fb, is, trace, varargin{:});
+  q = adaptlobstp (f, a, b, fa, fb, is, trace, varargin{:});
 endfunction
 
 ## ADAPTLOBSTP  Recursive function used by QUADL.
@@ -141,7 +148,7 @@
 ##
 ##   Walter Gautschi, 08/03/98
 
-function Q = adaptlobstp (f, a, b, fa, fb, is, trace, varargin)
+function q = adaptlobstp (f, a, b, fa, fb, is, trace, varargin)
   h = (b-a)/2;
   m = (a+b)/2;
   alpha = sqrt(2/3);
@@ -165,12 +172,12 @@
       warning ("quadl: required tolerance may not be met");
       need_warning (0);
     endif
-    Q = i1;
+    q = i1;
     if (trace)
-      disp ([a, b-a, Q]);
+      disp ([a, b-a, q]);
     endif
   else
-    Q = (adaptlobstp (f, a, mll, fa, fmll, is, trace, varargin{:})
+    q = (adaptlobstp (f, a, mll, fa, fmll, is, trace, varargin{:})
          + adaptlobstp (f, mll, ml, fmll, fml, is, trace, varargin{:})
          + adaptlobstp (f, ml, m, fml, fm, is, trace, varargin{:})
          + adaptlobstp (f, m, mr, fm, fmr, is, trace, varargin{:})
--- a/scripts/general/quadv.m
+++ b/scripts/general/quadv.m
@@ -21,32 +21,43 @@
 ## @deftypefnx {Function File} {@var{q} =} quadv (@var{f}, @var{a}, @var{b}, @var{tol})
 ## @deftypefnx {Function File} {@var{q} =} quadv (@var{f}, @var{a}, @var{b}, @var{tol}, @var{trace})
 ## @deftypefnx {Function File} {@var{q} =} quadv (@var{f}, @var{a}, @var{b}, @var{tol}, @var{trace}, @var{p1}, @var{p2}, @dots{})
-## @deftypefnx {Function File} {[@var{q}, @var{fcnt}] =} quadv (@dots{})
+## @deftypefnx {Function File} {[@var{q}, @var{nfun}] =} quadv (@dots{})
 ##
 ## Numerically evaluate the integral of @var{f} from @var{a} to @var{b}
-## using adaptive Simpson's rule.
-## @var{f} is either a function handle, inline function or string
+## using an adaptive Simpson's rule.
+## @var{f} is a function handle, inline function, or string
 ## containing the name of the function to evaluate.
-## The function defined by @var{f} may be a scalar, vector or array-valued.
+## @code{quadv} is a vectorized version of @code{quad} and the function
+## defined by @var{f} must accept a scalar or vector as input and return a
+## scalar, vector, or array as output.
 ##
-## If a value for @var{tol} is given, it defines the tolerance used to stop
-## the adaptation procedure, otherwise the default value of 1e-6 is used.
+## @var{a} and @var{b} are the lower and upper limits of integration.  Both
+## limits must be finite.
+## 
+## The optional argument @var{tol} defines the tolerance used to stop
+## the adaptation procedure.  The default value is @math{1e^{-6}}.
 ##
-## The algorithm used by @code{quadv}, involves recursively subdividing the
-## integration interval and applying Simpson's rule on each sub-interval.
-## If @var{trace} is @var{true}, after computing each of these partial
-## integrals, display the total number of function evaluations, the left end
-## of the sub-interval, the length of the sub-interval and the approximation
-## of the integral over the sub-interval.
+## The algorithm used by @code{quadv} involves recursively subdividing the
+## integration interval and applying Simpson's rule on each subinterval.
+## If @var{trace} is true then after computing each of these partial
+## integrals display: (1) the total number of function evaluations,
+## (2) the left end of the subinterval, (3) the length of the subinterval,
+## (4) the approximation of the integral over the subinterval.
 ##
-## Additional arguments @var{p1}, etc., are passed directly to @var{f}.
-## To use default values for @var{tol} and @var{trace}, one may pass
-## empty matrices.
+## Additional arguments @var{p1}, etc., are passed directly to the function
+## @var{f}.  To use default values for @var{tol} and @var{trace}, one may pass
+## empty matrices ([]).
 ##
-##@seealso{quad, quadl, quadgk, quadcc, trapz, dblquad, triplequad}
+## The result of the integration is returned in @var{q}.  @var{nfun} indicates
+## the number of function evaluations that were made.
+##
+## Note: @code{quadv} is written in Octave's scripting language and can be
+## used recursively in @code{dblquad} and @code{triplequad}, unlike the 
+## similar @code{quad} function.
+## @seealso{quad, quadl, quadgk, quadcc, trapz, dblquad, triplequad}
 ## @end deftypefn
 
-function [q, fcnt] = quadv (f, a, b, tol, trace, varargin)
+function [q, nfun] = quadv (f, a, b, tol, trace, varargin)
   if (nargin < 3)
     print_usage ();
   endif
@@ -73,7 +84,7 @@
   fa = feval (f, a, varargin{:});
   fc = feval (f, c, varargin{:});
   fb = feval (f, b, varargin{:});
-  fcnt = 3;
+  nfun = 3;
 
   ## If have edge singularities, move edge point by eps*(b-a) as
   ## discussed in Shampine paper used to implement quadgk
@@ -87,10 +98,10 @@
   h = (b - a);
   q = (b - a) / 6 * (fa + 4 * fc + fb);
 
-  [q, fcnt, hmin] = simpsonstp (f, a, b, c, fa, fb, fc, q, fcnt, abs (h),
+  [q, nfun, hmin] = simpsonstp (f, a, b, c, fa, fb, fc, q, nfun, abs (h),
                                 tol, trace, varargin{:});
 
-  if (fcnt > 10000)
+  if (nfun > 10000)
     warning ("maximum iteration count reached");
   elseif (isnan (q) || isinf (q))
     warning ("infinite or NaN function evaluations were returned");
@@ -99,16 +110,16 @@
   endif
 endfunction
 
-function [q, fcnt, hmin] = simpsonstp (f, a, b, c, fa, fb, fc, q0,
-                                       fcnt, hmin, tol, trace, varargin)
-  if (fcnt > 10000)
+function [q, nfun, hmin] = simpsonstp (f, a, b, c, fa, fb, fc, q0,
+                                       nfun, hmin, tol, trace, varargin)
+  if (nfun > 10000)
     q = q0;
   else
     d = (a + c) / 2;
     e = (c + b) / 2;
     fd = feval (f, d, varargin{:});
     fe = feval (f, e, varargin{:});
-    fcnt += 2;
+    nfun += 2;
     q1 = (c - a) / 6 * (fa + 4 * fd + fc);
     q2 = (b - c) / 6 * (fc + 4 * fe + fb);
     q = q1 + q2;
@@ -118,14 +129,14 @@
     endif
 
     if (trace)
-      disp ([fcnt, a, b-a, q]);
+      disp ([nfun, a, b-a, q]);
     endif
 
     ## Force at least one adpative step.
-    if (fcnt == 5 || abs (q - q0) > tol)
-      [q1, fcnt, hmin] = simpsonstp (f, a, c, d, fa, fc, fd, q1, fcnt, hmin,
+    if (nfun == 5 || abs (q - q0) > tol)
+      [q1, nfun, hmin] = simpsonstp (f, a, c, d, fa, fc, fd, q1, nfun, hmin,
                                     tol, trace, varargin{:});
-      [q2, fcnt, hmin] = simpsonstp (f, c, b, e, fc, fb, fe, q2, fcnt, hmin,
+      [q2, nfun, hmin] = simpsonstp (f, c, b, e, fc, fb, fe, q2, nfun, hmin,
                                      tol, trace, varargin{:});
       q = q1 + q2;
     endif
--- a/scripts/general/trapz.m
+++ b/scripts/general/trapz.m
@@ -17,15 +17,38 @@
 ## <http://www.gnu.org/licenses/>.
 
 ## -*- texinfo -*-
-## @deftypefn  {Function File} {@var{z} =} trapz (@var{y})
-## @deftypefnx {Function File} {@var{z} =} trapz (@var{x}, @var{y})
-## @deftypefnx {Function File} {@var{z} =} trapz (@dots{}, @var{dim})
+## @deftypefn  {Function File} {@var{q} =} trapz (@var{y})
+## @deftypefnx {Function File} {@var{q} =} trapz (@var{x}, @var{y})
+## @deftypefnx {Function File} {@var{q} =} trapz (@dots{}, @var{dim})
+##
+## Numerically evaluate the integral of points @var{y} using the trapezoidal
+## method.
+## @w{@code{trapz (@var{y})}} computes the integral of @var{y} along the first
+## non-singleton dimension.  When the argument @var{x} is omitted an
+## equally spaced @var{x} vector with unit spacing (1) is assumed.  
+## @code{trapz (@var{x}, @var{y})} evaluates the integral with respect
+## to the spacing in @var{x} and the values in @var{y}.  This is useful if
+## the points in @var{y} have been sampled unevenly.
+## If the optional @var{dim} argument is given, operate along this dimension.
 ##
-## Numerical integration using trapezoidal method.  @code{trapz
-## (@var{y})} computes the integral of the @var{y} along the first
-## non-singleton dimension.  If the argument @var{x} is omitted a
-## equally spaced vector is assumed.  @code{trapz (@var{x}, @var{y})}
-## evaluates the integral with respect to @var{x}.
+## If @var{x} is not specified then unit spacing will be used.  To scale
+## the integral to the correct value you must multiply by the actual spacing
+## value (deltaX).  As an example, the integral of @math{x^3} over the range
+## [0, 1] is @math{x^4/4} or 0.25.  The following code uses @code{trapz} to
+## calculate the integral in three different ways.
+##
+## @example
+## @group
+## x = 0:0.1:1;
+## y = x.^3;
+## q = trapz (y)
+##   @result{} q = 2.525   # No scaling
+## q * 0.1
+##   @result{} q = 0.2525  # Approximation to integral by scaling
+## trapz (x, y) 
+##   @result{} q = 0.2525  # Same result by specifying @var{x}
+## @end group
+## @end example
 ##
 ## @seealso{cumtrapz}
 ## @end deftypefn
--- a/scripts/general/triplequad.m
+++ b/scripts/general/triplequad.m
@@ -17,31 +17,45 @@
 ## <http://www.gnu.org/licenses/>.
 
 ## -*- texinfo -*-
-## @deftypefn {Function File} {} triplequad (@var{f}, @var{xa}, @var{xb}, @var{ya}, @var{yb}, @var{za}, @var{zb}, @var{tol}, @var{quadf}, @dots{})
-## Numerically evaluate a triple integral.  The function over which to
-## integrate is defined by @code{@var{f}}, and the interval for the
-## integration is defined by @code{[@var{xa}, @var{xb}, @var{ya},
-## @var{yb}, @var{za}, @var{zb}]}.  The function @var{f} must accept a
-## vector @var{x} and a scalar @var{y}, and return a vector of the same
-## length as @var{x}.
+## @deftypefn  {Function File} {} triplequad (@var{f}, @var{xa}, @var{xb}, @var{ya}, @var{yb}, @var{za}, @var{zb})
+## @deftypefnx {Function File} {} triplequad (@var{f}, @var{xa}, @var{xb}, @var{ya}, @var{yb}, @var{za}, @var{zb}, @var{tol})
+## @deftypefnx {Function File} {} triplequad (@var{f}, @var{xa}, @var{xb}, @var{ya}, @var{yb}, @var{za}, @var{zb}, @var{tol}, @var{quadf})
+## @deftypefnx {Function File} {} triplequad (@var{f}, @var{xa}, @var{xb}, @var{ya}, @var{yb}, @var{za}, @var{zb}, @var{tol}, @var{quadf}, @dots{})
+## Numerically evaluate the triple integral of @var{f}.
+## @var{f} is a function handle, inline function, or string
+## containing the name of the function to evaluate.  The function @var{f} must
+## have the form @math{w = f(x,y,z)} where either @var{x} or @var{y} is a
+## vector and the remaining inputs are scalars.  It should return a vector of
+## the same length and orientation as @var{x} or @var{y}.
 ##
-## If defined, @var{tol} defines the absolute tolerance to which to
-## which to integrate each sub-integral.
+## @var{xa}, @var{ya}, @var{za} and @var{xb}, @var{yb}, @var{zb} are the lower
+## and upper limits of integration for x, y, and z respectively.  The
+## underlying integrator determines whether infinite bounds are accepted.
+##
+## The optional argument @var{tol} defines the absolute tolerance used to 
+## integrate each sub-integral.  The default value is @math{1e^{-6}}.
+##
+## The optional argument @var{quadf} specifies which underlying integrator
+## function to use.  Any choice but @code{quad} is available and the default
+## is @code{quadcc}.
 ##
 ## Additional arguments, are passed directly to @var{f}.  To use the default
-## value for @var{tol} one may pass an empty matrix.
+## value for @var{tol} or @var{quadf} one may pass ':' or an empty matrix ([]).
 ## @seealso{dblquad, quad, quadv, quadl, quadgk, quadcc, trapz}
 ## @end deftypefn
 
-function Q = triplequad(f, xa, xb, ya, yb, za, zb, tol, quadf, varargin)
+function q = triplequad (f, xa, xb, ya, yb, za, zb, tol = 1e-6, quadf = @quadcc, varargin)
+
   if (nargin < 7)
     print_usage ();
   endif
-  if (nargin < 8 || isempty (tol))
+
+  ## Allow use of empty matrix ([]) to indicate default 
+  if (isempty (tol))
     tol = 1e-6;
   endif
-  if (nargin < 9 || isempty (quadf))
-    quadf = @quadgk;
+  if (isempty (quadf))
+    quadf = @quadcc;
   endif
 
   inner = @__triplequad_inner__;
@@ -50,17 +64,22 @@
     varargin = {};
   endif
 
-  Q = dblquad(@(y, z) inner (y, z, f, xa, xb, tol, quadf, varargin{:}),ya, yb, za, zb, tol);
+  q = dblquad (@(y, z) inner (y, z, f, xa, xb, tol, quadf, varargin{:}), ya, yb, za, zb, tol);
+
 endfunction
 
-function Q = __triplequad_inner__ (y, z, f, xa, xb, tol, quadf, varargin)
-  Q = zeros (size(y));
+function q = __triplequad_inner__ (y, z, f, xa, xb, tol, quadf, varargin)
+  q = zeros (size(y));
   for i = 1 : length (y)
-    Q(i) = feval (quadf, @(x) f (x, y(i), z, varargin{:}), xa, xb, tol);
+    q(i) = feval (quadf, @(x) f (x, y(i), z, varargin{:}), xa, xb, tol);
   endfor
 endfunction
 
-%% These tests are too expensive to run normally. Disable them
-% !#assert (triplequad (@(x,y,z) exp(-x.^2 - y.^2 - z.^2) , -1, 1, -1, 1, -1, 1, [],  @quadgk), pi ^ (3/2) * erf(1).^3, 1e-6)
-% !#assert (triplequad (@(x,y,z) exp(-x.^2 - y.^2 - z.^2) , -1, 1, -1, 1, -1, 1, [],  @quadl), pi ^ (3/2) * erf(1).^3, 1e-6)
-% !#assert (triplequad (@(x,y,z) exp(-x.^2 - y.^2 - z.^2) , -1, 1, -1, 1, -1, 1, [],  @quadv), pi ^ (3/2) * erf(1).^3, 1e-6)
+ 
+%!assert (triplequad (@(x,y,z) exp(-x.^2 - y.^2 - z.^2) , -1, 1, -1, 1, -1, 1, [],  @quadcc), pi ^ (3/2) * erf(1).^3, 1e-6)
+
+%% These tests are too expensive to run normally (~30 sec each).  Disable them
+#%!assert (triplequad (@(x,y,z) exp(-x.^2 - y.^2 - z.^2) , -1, 1, -1, 1, -1, 1, [],  @quadgk), pi ^ (3/2) * erf(1).^3, 1e-6)
+#%!#assert (triplequad (@(x,y,z) exp(-x.^2 - y.^2 - z.^2) , -1, 1, -1, 1, -1, 1, [],  @quadl), pi ^ (3/2) * erf(1).^3, 1e-6)
+#%!#assert (triplequad (@(x,y,z) exp(-x.^2 - y.^2 - z.^2) , -1, 1, -1, 1, -1, 1, [],  @quadv), pi ^ (3/2) * erf(1).^3, 1e-6)
+
--- a/scripts/image/image.m
+++ b/scripts/image/image.m
@@ -19,9 +19,9 @@
 ## -*- texinfo -*-
 ## @deftypefn  {Function File} {} image (@var{img})
 ## @deftypefnx {Function File} {} image (@var{x}, @var{y}, @var{img})
-## Display a matrix as a color image.  The elements of @var{x} are indices
+## Display a matrix as a color image.  The elements of @var{img} are indices
 ## into the current colormap, and the colormap will be scaled so that the
-## extremes of @var{x} are mapped to the extremes of the colormap.
+## extremes of @var{img} are mapped to the extremes of the colormap.
 ##
 ## The axis values corresponding to the matrix elements are specified in
 ## @var{x} and @var{y}.  If you're not using gnuplot 4.2 or later, these
new file mode 100644
--- /dev/null
+++ b/scripts/plot/colstyle.m
@@ -0,0 +1,89 @@
+## Copyright (C) 2011 David Bateman
+##
+## 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{style}, @var{color}, @var{marker}, @var{msg}] =} colstyle (@var{linespec})
+## Parse @var{linespec} and return the line style, color, and markers given.
+## In the case of an error, the string @var{msg} will return the text of the
+## error.
+## @end deftypefn
+
+function [l, c, m, msg] = colstyle (style)
+
+  if (nargin != 1)
+    print_usage ();
+  endif
+
+  if (! ischar (style))
+    error ("colstyle: STYLE must be a string");
+  endif
+
+  try
+    opt = __pltopt__ ("colstyle", style);
+    l = opt.linestyle;
+    c = opt.color;
+    m = opt.marker;
+    msg = [];
+    switch (c)
+      case [0 0 0]
+        c = "k";
+      case [1 0 0]
+        c = "r";
+      case [0 1 0]
+        c = "g";
+      case [0 0 1]
+        c = "b";
+      case [1 1 0]
+        c = "y";
+      case [1 0 1]
+        c = "m";
+      case [0 1 1]
+        c = "c";
+      case [0 1 1]
+        c = "w";
+    endswitch
+  catch
+    l = c = m = [];
+    msg = lasterr ();
+  end_try_catch
+
+endfunction
+
+%!test
+%! [l, c, m, msg] = colstyle ("r:x");
+%! assert (isempty (msg));
+%! assert (l, ":");
+%! assert (c, "r");
+%! assert (m, "x");
+
+%!test
+%! [l, c, m, msg] = colstyle (".");
+%! assert (isempty (msg));
+%! assert (l, "none");
+%! assert (c, []);
+%! assert (m, ".");
+
+%!test
+%! [l, c, m, msg] = colstyle ("~");
+%! assert (msg, "colstyle: unrecognized format character: `~'");
+
+%% Test input validation
+%!error colstyle ()
+%!error colstyle (1, 2)
+%!error colstyle (1.5)
+
--- a/scripts/plot/module.mk
+++ b/scripts/plot/module.mk
@@ -67,6 +67,7 @@
   plot/close.m \
   plot/closereq.m \
   plot/colorbar.m \
+  plot/colstyle.m \
   plot/comet.m \
   plot/comet3.m \
   plot/compass.m \
old mode 100755
new mode 100644
old mode 100755
new mode 100644
--- a/src/DLD-FUNCTIONS/__init_fltk__.cc
+++ b/src/DLD-FUNCTIONS/__init_fltk__.cc
@@ -648,13 +648,13 @@
     begin ();
     {
 
+      canvas = new
+        OpenGL_fltk (0, 0, ww , hh - status_h, number ());
+
       uimenu = new
         fltk_uimenu(0, 0, ww, menu_h);
       uimenu->hide ();
 
-      canvas = new
-        OpenGL_fltk (0, 0, ww , hh - status_h, number ());
-
       bottom = new
         Fl_Box (0,
                 hh - status_h,
--- a/src/DLD-FUNCTIONS/quad.cc
+++ b/src/DLD-FUNCTIONS/quad.cc
@@ -174,36 +174,31 @@
 
 DEFUN_DLD (quad, args, nargout,
   "-*- texinfo -*-\n\
-@deftypefn  {Loadable Function} {@var{v} =} quad (@var{f}, @var{a}, @var{b})\n\
-@deftypefnx {Loadable Function} {@var{v} =} quad (@var{f}, @var{a}, @var{b}, @var{tol})\n\
-@deftypefnx {Loadable Function} {@var{v} =} quad (@var{f}, @var{a}, @var{b}, @var{tol}, @var{sing})\n\
-@deftypefnx {Loadable Function} {[@var{v}, @var{ier}, @var{nfun}, @var{err}] =} quad (@dots{})\n\
-Integrate a nonlinear function of one variable using @sc{quadpack}.\n\
-The first argument is the name of the function, the function handle, or\n\
-the inline function to call to compute the value of the integrand.  It\n\
-must have the form\n\
+@deftypefn  {Loadable Function} {@var{q} =} quad (@var{f}, @var{a}, @var{b})\n\
+@deftypefnx {Loadable Function} {@var{q} =} quad (@var{f}, @var{a}, @var{b}, @var{tol})\n\
+@deftypefnx {Loadable Function} {@var{q} =} quad (@var{f}, @var{a}, @var{b}, @var{tol}, @var{sing})\n\
+@deftypefnx {Loadable Function} {[@var{q}, @var{ier}, @var{nfun}, @var{err}] =} quad (@dots{})\n\
+Numerically evaluate the integral of @var{f} from @var{a} to @var{b} using\n\
+Fortran routines from @w{@sc{quadpack}}.  @var{f} is a function handle, inline\n\
+function, or a string containing the name of the function to evaluate.\n\
+The function must have the form @code{y = f (x)} where @var{y} and @var{x}\n\
+are scalars.\n\
 \n\
-@example\n\
-y = f (x)\n\
-@end example\n\
-\n\
-@noindent\n\
-where @var{y} and @var{x} are scalars.\n\
-\n\
-The second and third arguments are limits of integration.  Either or\n\
-both may be infinite.\n\
+@var{a} and @var{b} are the lower and upper limits of integration.  Either\n\
+or both may be infinite.\n\
 \n\
 The optional argument @var{tol} is a vector that specifies the desired\n\
 accuracy of the result.  The first element of the vector is the desired\n\
 absolute tolerance, and the second element is the desired relative\n\
 tolerance.  To choose a relative test only, set the absolute\n\
 tolerance to zero.  To choose an absolute test only, set the relative\n\
-tolerance to zero.\n\
+tolerance to zero.  Both tolerances default to @code{sqrt(eps)} or\n\
+approximately @math{1.5e^{-8}}.\n\
 \n\
 The optional argument @var{sing} is a vector of values at which the\n\
 integrand is known to be singular.\n\
 \n\
-The result of the integration is returned in @var{v}.  @var{ier}\n\
+The result of the integration is returned in @var{q}.  @var{ier}\n\
 contains an integer error code (0 indicates a successful integration).\n\
 @var{nfun} indicates the number of function evaluations that were\n\
 made, and @var{err} contains an estimate of the error in the\n\
@@ -212,8 +207,9 @@
 The function @code{quad_options} can set other optional\n\
 parameters for @code{quad}.\n\
 \n\
-Note: because @code{quad} is written in Fortran it\n\
-cannot be called recursively.\n\
+Note: because @code{quad} is written in Fortran it cannot be called\n\
+recursively.  This prevents its use in integrating over more than one\n\
+variable by routines @code{dblquad} and @code{triplequad}.\n\
 @seealso{quad_options, quadv, quadl, quadgk, quadcc, trapz, dblquad, triplequad}\n\
 @end deftypefn")
 {
--- a/src/DLD-FUNCTIONS/quadcc.cc
+++ b/src/DLD-FUNCTIONS/quadcc.cc
@@ -50,8 +50,7 @@
   int depth, rdepth, ndiv;
 } cquad_ival;
 
-/* Some constants and matrices that we'll need.
-    */
+/* Some constants and matrices that we'll need.  */
 
 static const double xi[33] = {
   -1., -0.99518472667219688624, -0.98078528040323044912,
@@ -1473,17 +1472,63 @@
 }
 
 
-
-/* The actual integration routine.
-    */
+/* The actual integration routine.  */
 
 DEFUN_DLD (quadcc, args, nargout,
 "-*- texinfo -*-\n\
-@deftypefn  {Function File} {[@var{int}, @var{err}, @var{nr_points}] =} quadcc (@var{f}, @var{a}, @var{b}, @var{tol})\n\
-@deftypefnx {Function File} {[@var{int}, @var{err}, @var{nr_points}] =} quadcc (@var{f}, @var{a}, @var{b}, @var{tol}, @var{sing})\n\
-Numerically evaluate an integral using the doubly-adaptive\n\
-Clenshaw-Curtis quadrature described by P. Gonnet in @cite{Increasing the\n\
-Reliability of Adaptive Quadrature Using Explicit Interpolants}.\n\
+@deftypefn  {Function File} {@var{q} =} quadcc (@var{f}, @var{a}, @var{b})\n\
+@deftypefnx {Function File} {@var{q} =} quadcc (@var{f}, @var{a}, @var{b}, @var{tol})\n\
+@deftypefnx {Function File} {@var{q} =} quadcc (@var{f}, @var{a}, @var{b}, @var{tol}, @var{sing})\n\
+@deftypefnx {Function File} {[@var{q}, @var{err}, @var{nr_points}] =} quadcc (@dots{})\n\
+Numerically evaluate the integral of @var{f} from @var{a} to @var{b}\n\
+using the doubly-adaptive Clenshaw-Curtis quadrature described by P. Gonnet\n\
+in @cite{Increasing the Reliability of Adaptive Quadrature Using Explicit\n\
+Interpolants}.\n\
+@var{f} is a function handle, inline function, or string\n\
+containing the name of the function to evaluate.\n\
+The function @var{f} must be vectorized and must return a vector of output\n\
+values if given a vector of input values.  For example,\n\
+\n\
+@example\n\
+   f = @@(x) x .* sin (1./x) .* sqrt (abs (1 - x));\n\
+@end example\n\
+\n\
+@noindent\n\
+which uses the element-by-element `dot' form for all operators.\n\
+\n\
+@var{a} and @var{b} are the lower and upper limits of integration.  Either\n\
+or both limits may be infinite.  @code{quadcc} handles an inifinite limit\n\
+by substituting the variable of integration with @code{x=tan(pi/2*u)}.\n\
+\n\
+The optional argument @var{tol} defines the relative tolerance used to stop\n\
+the integration procedure.  The default value is @math{1e^{-6}}.\n\
+\n\
+The optional argument @var{sing} contains a list of points where the\n\
+integrand has known singularities, or discontinuities\n\
+in any of its derivatives, inside the integration interval.\n\
+For the example above, which has a discontinuity at x=1, the call to\n\
+@code{quadcc} would be as follows\n\
+\n\
+@example\n\
+   int = quadcc (f, a, b, 1.0e-6, [ 1 ]);\n\
+@end example\n\
+\n\
+The result of the integration is returned in @var{q}.\n\
+@var{err} is an estimate of the absolute integration error and\n\
+@var{nr_points} is the number of points at which the integrand was evaluated.\n\
+If the adaptive integration did not converge, the value of\n\
+@var{err} will be larger than the requested tolerance.  Therefore, it is\n\
+recommended to verify this value for difficult integrands.\n\
+\n\
+@code{quadcc} is capable of dealing with non-numeric\n\
+values of the integrand such as @code{NaN} or @code{Inf}.\n\
+If the integral diverges, and @code{quadcc} detects this, \n\
+then a warning is issued and @code{Inf} or @code{-Inf} is returned.\n\
+\n\
+Note: @code{quadcc} is a general purpose quadrature algorithm\n\
+and, as such, may be less efficient for a smooth or otherwise\n\
+well-behaved integrand than other methods such as @code{quadgk}.\n\
+\n\
 The algorithm uses Clenshaw-Curtis quadrature rules of increasing\n\
 degree in each interval and bisects the interval if either the\n\
 function does not appear to be smooth or a rule of maximum\n\
@@ -1491,63 +1536,13 @@
 L2-norm of the difference between two successive interpolations\n\
 of the integrand over the nodes of the respective quadrature rules.\n\
 \n\
-For example,\n\
-\n\
-@example\n\
-   int = quadcc (f, a, b, 1.0e-6);\n\
-@end example\n\
-\n\
-@noindent\n\
-computes the integral of a function @var{f} in the interval\n\
-[@var{a}, @var{b}] to the relative precision of six\n\
-decimal digits.\n\
-The integrand @var{f} should accept a vector argument and return a vector\n\
-result containing the integrand evaluated at each element of the\n\
-argument, for example:\n\
-\n\
-@example\n\
-   f = @@(x) x .* sin (1 ./ x) .* sqrt (abs (1 - x));\n\
-@end example\n\
-\n\
-If the integrand has known singularities or discontinuities\n\
-in any of its derivatives inside the interval,\n\
-as does the above example at x=1, these can be specified in\n\
-the additional argument @var{sing} as follows\n\
-\n\
-@example\n\
-   int = quadcc (f, a, b, 1.0e-6, [ 1 ]);\n\
-@end example\n\
-\n\
-The two additional output variables @var{err} and @var{nr_points}\n\
-return an estimate of the absolute integration error and\n\
-the number of points at which the integrand was evaluated\n\
-respectively.\n\
-If the adaptive integration did not converge, the value of\n\
-@var{err} will be larger than the requested tolerance.  It is\n\
-therefore recommended to verify this value for difficult\n\
-integrands.\n\
-\n\
-If either @var{a} or @var{b} are @code{+/-Inf}, @code{quadcc}\n\
-integrates @var{f} by substituting the variable of integration\n\
-with @code{x=tan(pi/2*u)}.\n\
-\n\
-@code{quadcc} is capable of dealing with non-numeric\n\
-values of the integrand such as @code{NaN} or @code{Inf}\n\
-, as in the above example at x=0.\n\
-If the integral diverges and @code{quadcc} detects this, \n\
-a warning is issued and @code{Inf} or @code{-Inf} is returned.\n\
-\n\
-Note that @code{quadcc} is a general purpose quadrature algorithm\n\
-and as such may be less efficient for smooth or otherwise\n\
-well-behaved integrand than other methods such as\n\
-@code{quadgk} or @code{trapz}.\n\
-\n\
 Reference: P. Gonnet, @cite{Increasing the Reliability of Adaptive\n\
 Quadrature Using Explicit Interpolants}, ACM Transactions on\n\
 Mathematical Software, Vol. 37, Issue 3, Article No. 3, 2010.\n\
 @seealso{quad, quadv, quadl, quadgk, trapz, dblquad, triplequad}\n\
 @end deftypefn")
 {
+  octave_value_list retval;
 
   /* Some constants that we will need. */
   static const int n[4] = { 4, 8, 16, 32 };
@@ -1566,11 +1561,11 @@
   double a, b, tol, iivals[cquad_heapsize], *sing;
 
   /* Variables needed for transforming the integrand. */
-  int wrap = 0;
+  bool wrap = false;
   double xw;
 
   /* Stuff we will need to call the integrand. */
-  octave_value_list fargs, retval;
+  octave_value_list fargs, fvals;
 
   /* Actual variables (as opposed to constants above). */
   double m, h, ml, hl, mr, hr, temp;
@@ -1583,48 +1578,49 @@
 
 
   /* Parse the input arguments. */
-  if (nargin < 1)
+  if (nargin < 3)
     {
-      error
-        ("quadcc: first argument (integrand) of type function handle required");
-      return octave_value_list ();
+      print_usage ();
+      return retval;
     }
+
+  if (args(0).is_function_handle () || args(0).is_inline_function ())
+    fcn = args(0).function_value ();
   else
     {
-      if (args (0).is_function_handle () || args (0).is_inline_function ())
-        fcn = args (0).function_value ();
-      else
-        {
-          error ("quadcc: first argument (integrand) must be a function handle or an inline function");
-          return octave_value_list();
-        }
+       std::string fcn_name = unique_symbol_name ("__quadcc_fcn_");
+       std::string fname = "function y = ";
+       fname.append (fcn_name);
+       fname.append ("(x) y = ");
+       fcn = extract_function (args(0), "quadcc", fcn_name, fname,
+                               "; endfunction");
     }
 
-  if (nargin < 2 || !args (1).is_real_scalar ())
+  if (!args(1).is_real_scalar ())
     {
-      error ("quadcc: second argument (left interval edge) must be a single real scalar");
-      return octave_value_list ();
+      error ("quadcc: lower limit of integration (A) must be a single real scalar");
+      return retval;
     }
   else
-    a = args (1).double_value ();
+    a = args(1).double_value ();
 
-  if (nargin < 3 || !args (2).is_real_scalar ())
+  if (!args(2).is_real_scalar ())
     {
-      error ("quadcc: third argument (right interval edge) must be a single real scalar");
-      return octave_value_list ();
+      error ("quadcc: upper limit of integration (B) must be a single real scalar");
+      return retval;
     }
   else
-    b = args (2).double_value ();
+    b = args(2).double_value ();
 
-  if (nargin < 4)
+  if (nargin < 4 || args(3).is_empty ())
     tol = 1.0e-6;
-  else if (!args (3).is_real_scalar ())
+  else if (!args(3).is_real_scalar () || args(3).double_value () <= 0)
     {
-      error ("quadcc: fourth argument (tolerance) must be a single real scalar");
-      return octave_value_list ();
+      error ("quadcc: tolerance (TOL) must be a single real scalar > 0");
+      return retval;
     }
   else
-    tol = args (3).double_value ();
+    tol = args(3).double_value ();
 
   if (nargin < 5)
     {
@@ -1632,20 +1628,21 @@
       iivals[0] = a;
       iivals[1] = b;
     }
-  else if (!(args (4).is_real_scalar () || args (4).is_real_matrix ()))
+  else if (!(args(4).is_real_scalar () || args(4).is_real_matrix ()))
     {
-      error ("quadcc: fifth argument (singularities) must be a vector of real values");
-      return octave_value_list ();
+      error ("quadcc: list of singularities (SING) must be a vector of real values");
+      return retval;
     }
   else
     {
-      nivals = 1 + args (4).length ();
-      if ( nivals > cquad_heapsize ) {
-        error ("quadcc: maximum number of singular points is limited to %i",
-               cquad_heapsize-1);
-        return octave_value_list();
+      nivals = 1 + args(4).length ();
+      if (nivals > cquad_heapsize) 
+        {
+          error ("quadcc: maximum number of singular points is limited to %i",
+                 cquad_heapsize-1);
+          return retval;
         }
-      sing = args (4).array_value ().fortran_vec ();
+      sing = args(4).array_value ().fortran_vec ();
       iivals[0] = a;
       for (i = 0; i < nivals - 2; i++)
         iivals[i + 1] = sing[i];
@@ -1655,7 +1652,7 @@
   /* If a or b are +/-Inf, transform the integral. */
   if (xisinf (a) || xisinf (b))
     {
-      wrap = 1;
+      wrap = true;
       for (i = 0; i <= nivals; i++)
         if (xisinf (iivals[i]))
           iivals[i] = copysign (1.0, iivals[i]);
@@ -1691,19 +1688,18 @@
           for (i = 0; i <= n[3]; i++)
             ex (i) = m + xi[i] * h;
         }
-      fargs (0) = ex;
-      retval = feval (fcn, fargs, 1);
-      if (retval.length () != 1 || !retval (0).is_real_matrix ())
+      fargs(0) = ex;
+      fvals = feval (fcn, fargs, 1);
+      if (fvals.length () != 1 || !fvals(0).is_real_matrix ())
         {
-          error
-            ("quadcc: integrand must return a single, real-valued vector");
-          return octave_value_list ();
+          error ("quadcc: integrand F must return a single, real-valued vector");
+          return retval;
         }
-      Matrix effex = retval (0).matrix_value ();
+      Matrix effex = fvals(0).matrix_value ();
       if (effex.length () != ex.length ())
         {
-          error ("quadcc: integrand must return a single, real-valued vector of the same size as the input");
-          return octave_value_list ();
+          error ("quadcc: integrand F must return a single, real-valued vector of the same size as the input");
+          return retval;
         }
       for (i = 0; i <= n[3]; i++)
         {
@@ -1812,18 +1808,18 @@
                 for (i = 0; i < n[d] / 2; i++)
                   ex (i) = m + xi[(2 * i + 1) * skip[d]] * h;
               }
-            fargs (0) = ex;
-            retval = feval (fcn, fargs, 1);
-            if (retval.length () != 1 || !retval (0).is_real_matrix ())
+            fargs(0) = ex;
+            fvals = feval (fcn, fargs, 1);
+            if (fvals.length () != 1 || !fvals(0).is_real_matrix ())
               {
-                error ("quadcc: integrand must return a single, real-valued vector");
-                return octave_value_list ();
+                error ("quadcc: integrand F must return a single, real-valued vector");
+                return retval;
               }
-            Matrix effex = retval (0).matrix_value ();
+            Matrix effex = fvals(0).matrix_value ();
             if (effex.length () != ex.length ())
               {
-                error ("quadcc: integrand must return a single, real-valued vector of the same size as the input");
-                return octave_value_list ();
+                error ("quadcc: integrand F must return a single, real-valued vector of the same size as the input");
+                return retval;
               }
             neval += effex.length ();
             for (i = 0; i < n[d] / 2; i++)
@@ -1960,18 +1956,18 @@
                 for (i = 0; i < n[0] - 1; i++)
                   ex (i) = ml + xi[(i + 1) * skip[0]] * hl;
               }
-            fargs (0) = ex;
-            retval = feval (fcn, fargs, 1);
-            if (retval.length () != 1 || !retval (0).is_real_matrix ())
+            fargs(0) = ex;
+            fvals = feval (fcn, fargs, 1);
+            if (fvals.length () != 1 || !fvals(0).is_real_matrix ())
               {
-                error ("quadcc: integrand must return a single, real-valued vector");
-                return octave_value_list ();
+                error ("quadcc: integrand F must return a single, real-valued vector");
+                return retval;
               }
-            Matrix effex = retval (0).matrix_value ();
+            Matrix effex = fvals(0).matrix_value ();
             if (effex.length () != ex.length ())
               {
-                error ("quadcc: integrand must return a single, real-valued vector of the same size as the input");
-                return octave_value_list ();
+                error ("quadcc: integrand F must return a single, real-valued vector of the same size as the input");
+                return retval;
               }
             neval += effex.length ();
             for (i = 0; i < n[0] - 1; i++)
@@ -2056,18 +2052,18 @@
                 for (i = 0; i < n[0] - 1; i++)
                   ex (i) = mr + xi[(i + 1) * skip[0]] * hr;
               }
-            fargs (0) = ex;
-            retval = feval (fcn, fargs, 1);
-            if (retval.length () != 1 || !retval (0).is_real_matrix ())
+            fargs(0) = ex;
+            fvals = feval (fcn, fargs, 1);
+            if (fvals.length () != 1 || !fvals(0).is_real_matrix ())
               {
-                error ("quadcc: integrand must return a single, real-valued vector");
-                return octave_value_list ();
+                error ("quadcc: integrand F must return a single, real-valued vector");
+                return retval;
               }
-            Matrix effex = retval (0).matrix_value ();
+            Matrix effex = fvals(0).matrix_value ();
             if (effex.length () != ex.length ())
               {
-                error ("quadcc: integrand must return a single, real-valued vector of the same size as the input");
-                return octave_value_list ();
+                error ("quadcc: integrand F must return a single, real-valued vector of the same size as the input");
+                return retval;
               }
             neval += effex.length ();
             for (i = 0; i < n[0] - 1; i++)
@@ -2237,11 +2233,39 @@
     }
 */
   /* Clean up and present the results. */
-  retval (0) = igral;
+  if (nargout > 2)
+    retval(2) = neval;
   if (nargout > 1)
-    retval (1) = err;
-  if (nargout > 2)
-    retval (2) = neval;
+    retval(1) = err;
+  retval(0) = igral;
   /* All is well that ends well. */
   return retval;
 }
+
+
+/* 
+
+%!assert (quadcc(@sin,-pi,pi), 0, 1e-6)
+%!assert (quadcc(inline('sin'),-pi,pi), 0, 1e-6)
+%!assert (quadcc('sin',-pi,pi), 0, 1e-6)
+
+%!assert (quadcc(@sin,-pi,0), -2, 1e-6)
+%!assert (quadcc(@sin,0,pi), 2, 1e-6)
+%!assert (quadcc(@(x) 1./sqrt(x), 0, 1), 2, 1e-6)
+%!assert (quadcc(@(x) 1./(sqrt(x).*(x+1)), 0, Inf), pi, 1e-6)
+
+%!assert (quadcc (@(x) exp(-x .^ 2), -Inf, Inf), sqrt(pi), 1e-6)
+%!assert (quadcc (@(x) exp(-x .^ 2), -Inf, 0), sqrt(pi)/2, 1e-6)
+
+%% Test input validation
+%!error (quadcc ())
+%!error (quadcc (@sin))
+%!error (quadcc (@sin, 0))
+%!error (quadcc (@sin, ones(2), pi))
+%!error (quadcc (@sin, -i, pi))
+%!error (quadcc (@sin, 0, ones(2)))
+%!error (quadcc (@sin, 0, i))
+%!error (quadcc (@sin, 0, pi, 0))
+%!error (quadcc (@sin, 0, pi, 1e-6, [ i ]))
+
+*/
--- a/src/help.cc
+++ b/src/help.cc
@@ -356,7 +356,7 @@
 
   pair_type ("function",
     "-*- texinfo -*-\n\
-@deftypefn {Keyword} {} function @var{outputs} = function (@var{input}, @dots{})\n\
+@deftypefn  {Keyword} {} function @var{outputs} = function (@var{input}, @dots{})\n\
 @deftypefnx {Keyword} {} function {} function (@var{input}, @dots{})\n\
 @deftypefnx {Keyword} {} function @var{outputs} = function\n\
 Begin a function body with @var{outputs} as results and @var{inputs} as\n\
@@ -382,7 +382,7 @@
 
   pair_type ("if",
     "-*- texinfo -*-\n\
-@deftypefn {Keyword} {} if (@var{cond}) @dots{} endif\n\
+@deftypefn  {Keyword} {} if (@var{cond}) @dots{} endif\n\
 @deftypefnx {Keyword} {} if (@var{cond}) @dots{} else @dots{} endif\n\
 @deftypefnx {Keyword} {} if (@var{cond}) @dots{} elseif (@var{cond}) @dots{} endif\n\
 @deftypefnx {Keyword} {} if (@var{cond}) @dots{} elseif (@var{cond}) @dots{} else @dots{} endif\n\
@@ -473,7 +473,7 @@
 execution will proceed after the catch block (though it is often\n\
 recommended to use the lasterr function to re-throw the error after cleanup\n\
 is completed).\n\
-@seealso{catch,unwind_protect}\n\
+@seealso{catch, unwind_protect}\n\
 @end deftypefn"),
 
   pair_type ("until",
@@ -494,7 +494,7 @@
 unwind_protect_cleanup block is still executed (in other words, the\n\
 unwind_protect_cleanup will be run with or without an error in the\n\
 unwind_protect block).\n\
-@seealso{unwind_protect_cleanup,try}\n\
+@seealso{unwind_protect_cleanup, try}\n\
 @end deftypefn"),
 
   pair_type ("unwind_protect_cleanup",
--- a/src/load-path.cc
+++ b/src/load-path.cc
@@ -2037,7 +2037,7 @@
 DEFUN (restoredefaultpath, , ,
     "-*- texinfo -*-\n\
 @deftypefn {Built-in Function} {} restoredefaultpath (@dots{})\n\
-Restore Octave's path to it's initial state at startup.\n\
+Restore Octave's path to its initial state at startup.\n\
 \n\
 @seealso{path, addpath, rmpath, genpath, pathdef, savepath, pathsep}\n\
 @end deftypefn")
--- a/src/ov-class.cc
+++ b/src/ov-class.cc
@@ -121,12 +121,15 @@
 {
   std::string retval = class_name ();
 
-  octave_function *fcn = octave_call_stack::current ();
+  if (nparents () > 0)
+    {
+      octave_function *fcn = octave_call_stack::current ();
 
-  // Here we are just looking to see if FCN is a method or constructor
-  // for any class, not specifically this one.
-  if (fcn->is_class_method () || fcn->is_class_constructor ())
-    retval = fcn->dispatch_class ();
+      // Here we are just looking to see if FCN is a method or constructor
+      // for any class, not specifically this one.
+      if (fcn && (fcn->is_class_method () || fcn->is_class_constructor ()))
+        retval = fcn->dispatch_class ();
+    }
 
   return retval;
 }
--- a/src/ov-class.h
+++ b/src/ov-class.h
@@ -127,10 +127,18 @@
   size_t nparents (void) const { return parent_list.size (); }
 
   octave_value reshape (const dim_vector& new_dims) const
-    { return map.reshape (new_dims); }
+    { 
+      octave_class retval = octave_class (*this);
+      retval.map = retval.map_value().reshape (new_dims);
+      return octave_value (new octave_class (retval));
+    }
 
   octave_value resize (const dim_vector& dv, bool = false) const
-    { octave_map tmap = map; tmap.resize (dv); return tmap; }
+    { 
+      octave_class retval = octave_class (*this);
+      retval.map.resize (dv);
+      return octave_value (new octave_class (retval));
+    }
 
   bool is_defined (void) const { return true; }