# HG changeset patch
# User Nicholas R. Jankowski <jankowskin@asme.org>
# Date 1508263999 14400
#      Tue Oct 17 14:13:19 2017 -0400
# Node ID 907d1ef5a6c21ba200fb4ecd2fcbc3cecda5272e
# Parent  ea06ac5003480d790cb2bb0e3cb687d2a817acdf
integral.m: update to use recently changed quadcc and preserve type single (bug #42073)

* scripts/general/integral.m: integral will now make general use of quadcc for all but arrayvalued and waypoints calls

diff -r ea06ac500348 -r 907d1ef5a6c2 scripts/general/integral.m
--- a/scripts/general/integral.m	Tue Oct 17 08:58:58 2017 -0700
+++ b/scripts/general/integral.m	Tue Oct 17 14:13:19 2017 -0400
@@ -23,10 +23,11 @@
 ## Numerically evaluate the integral of @var{f} from @var{a} to @var{b} using
 ## adaptive quadrature.
 ##
-## @code{integral} is a wrapper for @code{quadgk} (for scalar integrands) and
-## @code{quadv} (for array-valued integrands) intended to provide @sc{matlab}
-## compatibility.  More control of the numerical integration may be achievable
-## by calling the various quadrature functions directly.
+## @code{integral} is a wrapper for @code{quadcc} (for general scalar 
+## integrands), @code{quadgk} (for integrals with specified integratio paths), 
+## and @code{quadv} (for array-valued integrands) intended to provide 
+## @sc{matlab} compatibility.  More control of the numerical integration may be 
+## achievable by calling the various quadrature functions directly.
 ##
 ## @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
@@ -52,7 +53,8 @@
 ## @code{integral} expects @var{f} to return a scalar value unless
 ## @var{arrayvalued} is specified as true.  This option will cause
 ## @code{integral} to perform the integration over the entire array and return
-## @var{q} with the same dimensions as returned by @var{f}.
+## @var{q} with the same dimensions as returned by @var{f}. For more detail see
+## @code{quadv}.
 ##
 ## @item AbsTol
 ## Define the absolute error tolerance for the quadrature.  The default
@@ -90,9 +92,10 @@
 ## inputs in place regardless of the class of the integration limits.
 ##
 ## @item
-## As a consequence of using @code{quadgk} and @code{quadv}, certain option
-## combinations are not supported.  Currently, @qcode{"ArrayValued"} cannot be
-## combined with @qcode{"RelTol"} or @qcode{"Waypoints"}.
+## As a consequence of using @code{quadcc}, @code{quadgk} and @code{quadv}, 
+## certain option combinations are not supported.  Currently, 
+## @qcode{"ArrayValued"} cannot be combined with @qcode{"RelTol"} or 
+## @qcode{"Waypoints"}. Any @qcode{"ArrayValued"} call will simply 
 ## @end enumerate
 ##
 ## @seealso{integral2, integral3, quad, quadgk, quadv, quadl, quadcc, trapz,
@@ -106,9 +109,10 @@
   endif
 
   if (nargin == 3)
+    issingle = isa (a, "single") || isa (b, "single");
     ## Pass the simplest case directly to general integrator.
-    ## Let quadgk function handle input checks on function and limits.
-    q = quadgk (f, a, b);
+    ## Let quadcc function handle input checks on function and limits.
+    q = quadcc (f, a, b);
   else
     ## Parse options to determine how to call integrator.
     abstol = [];
@@ -136,34 +140,51 @@
           error ("integral: unknown property '%s'", prop);
       endswitch
     endwhile
-
-    if (! arrayvalued)
-      ## quadgk will accept empty values of optional parameters.
-      q = quadgk (f, a, b,
-                  "AbsTol", abstol, "RelTol", reltol, "WayPoints", waypoints);
-    else
-      ## FIXME: Replace warning with array compatible call to quadgk
+    
+    issingle = isa (a, "single") || isa (b, "single") ...
+                 || isa (waypoints, "single");
+        
+    if (arrayvalued)
+      ## pass vector valued function to quadgv, checking for conflicting params
+      
+      ## FIXME: Replace warning when have array compatible call with waypoints
       if (! isempty (waypoints))
         warning(["integral: array-valued quadrature routine currently ", ...
                  "unable to handle WayPoints.  WayPoints are ignored."]);
       endif
 
-      ## FIXME: Remove warning once we have reltol compatible arrayval'd
-      ##        quad or arrayfun call to quadgk.
+      ## FIXME: Remove warning once we have reltol compatible arrayval'd quadfn
       if (! isempty (reltol))
         warning(["integral: array-valued quadrature only accepts AbsTol.", ...
                  "  RelTol ignored."]);
       endif
-
       if (isempty (abstol))
-        issingle = isa (a, "single") || isa (b, "single");
         abstol = ifelse (issingle, 1e-5, 1e-10);
       endif
 
       q = quadv (f, a, b, abstol);
+    
+    else
+      if (isempty (abstol))
+        abstol = ifelse (issingle, 1e-5, 1e-10);
+      endif
+
+      if (isempty (reltol))
+        reltol = ifelse (issingle, 1e-4, 1e-6);
+      endif
+
+      if (! isempty (waypoints))
+        
+        q = quadgk (f, a, b,
+                      "AbsTol", abstol, "RelTol", reltol, "WayPoints", waypoints);
+      else
+        q = quadcc (f, a, b, [abstol reltol]);
+      endif
     endif
   endif
 
+  q = ifelse(issingle, single (q), q); #preserve type single for output
+  
 endfunction
 
 
@@ -192,6 +213,16 @@
 %! assert (integral (f, 0, 1, "ArrayValued", true), 1./[1:5]-cos(1:5)./[1:5],
 %!         1e-10);
 
+%!test  # test single input/output
+%! assert (integral (@sin, 0, 1), cos(0)-cos(1), 1e-10);
+%! assert (class(integral (@sin, single (0), 1)), 'single');
+%! assert (class(integral (@sin, 0, single (1))), 'single');
+%! assert (class(integral (@sin, single (0), single (1))), 'single');
+%! assert (integral (@sin, 0, 1, "Waypoints", 0.5), cos(0)-cos(1), 1e-10);
+%! assert (class (integral (@sin, 0, 1, "Waypoints", single (0.5))), 'single');
+%! assert (class (integral (@sin, single (0), 1, "Waypoints", 0.5)), 'single');
+%! assert (class (integral (@sin, 0, single (1), "Waypoints", 0.5)), 'single');
+
 %!test
 %! f = @(x) x.^5 .* exp (-x) .* sin (x);
 %! assert (integral (f, 0, inf, "RelTol", 1e-8, "AbsTol", 1e-12), -15, -1e-8);
