jwe/octave: scripts/general/quadgk.m annotate

annotate scripts/general/quadgk.m @ 31068:451fb63a10a0

update integral to call quadgk for 'ArrayValued' integrations (bug #62468) * integral.m: Modify integrator selection so that calls with â€˜ArrayValuedâ€™ go to quadgk instead of quadv. Remove error checks for previously incompatible parameter combinations. Update docstring to remove mention of quadv, point 'ArrayValued' reference to quadgk, change returned error parameter description to match current behavior, and remove parameter incompatibility note. Add BIST to verify combined parameter functionality, and change BISTs checking quadv err parameter. * quadgk.m: Correct parameter name in docstring. * NEWS.8.md: Under General Improvements add note about quadgk now accepting 'ArrayValued' parameter and update integral improvement description of optional returned error parameter. Under Matlab Compatibility add note about integral now accepting all parameter combinations.

author	Nicholas R. Jankowski <jankowski.nicholas@gmail.com>
date	Thu, 02 Jun 2022 19:56:15 -0400
parents	12f8fb75fc30
children

rev	line source
27928 bd51beb6205e update formatting of copyright notices John W. Eaton <jwe@octave.org> parents: 27924 diff changeset	1 ########################################################################
bd51beb6205e update formatting of copyright notices John W. Eaton <jwe@octave.org> parents: 27924 diff changeset	2 ##
30569 796f54d4ddbf update Octave Project Developers copyright for the new year John W. Eaton <jwe@octave.org> parents: 30384 diff changeset	3 ## Copyright (C) 2008-2022 The Octave Project Developers
27923 b442ec6dda5c use centralized file for copyright info for individual contributors John W. Eaton <jwe@octave.org> parents: 27805 diff changeset	4 ##
27928 bd51beb6205e update formatting of copyright notices John W. Eaton <jwe@octave.org> parents: 27924 diff changeset	5 ## See the file COPYRIGHT.md in the top-level directory of this
bd51beb6205e update formatting of copyright notices John W. Eaton <jwe@octave.org> parents: 27924 diff changeset	6 ## distribution or <https://octave.org/copyright/>.
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	7 ##
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	8 ## This file is part of Octave.
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	9 ##
24534 194eb4bd202b maint: Update punctuation for GPL v3 license text. Rik <rik@octave.org> parents: 24208 diff changeset	10 ## Octave is free software: you can redistribute it and/or modify it
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	11 ## under the terms of the GNU General Public License as published by
24534 194eb4bd202b maint: Update punctuation for GPL v3 license text. Rik <rik@octave.org> parents: 24208 diff changeset	12 ## the Free Software Foundation, either version 3 of the License, or
22755 3a2b891d0b33 maint: Standardize Copyright formatting. Rik <rik@octave.org> parents: 22474 diff changeset	13 ## (at your option) any later version.
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	14 ##
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	15 ## Octave is distributed in the hope that it will be useful, but
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	16 ## WITHOUT ANY WARRANTY; without even the implied warranty of
22755 3a2b891d0b33 maint: Standardize Copyright formatting. Rik <rik@octave.org> parents: 22474 diff changeset	17 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
3a2b891d0b33 maint: Standardize Copyright formatting. Rik <rik@octave.org> parents: 22474 diff changeset	18 ## GNU General Public License for more details.
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	19 ##
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	20 ## You should have received a copy of the GNU General Public License
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	21 ## along with Octave; see the file COPYING. If not, see
24534 194eb4bd202b maint: Update punctuation for GPL v3 license text. Rik <rik@octave.org> parents: 24208 diff changeset	22 ## <https://www.gnu.org/licenses/>.
27928 bd51beb6205e update formatting of copyright notices John W. Eaton <jwe@octave.org> parents: 27924 diff changeset	23 ##
bd51beb6205e update formatting of copyright notices John W. Eaton <jwe@octave.org> parents: 27924 diff changeset	24 ########################################################################
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	25
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	26 ## -- texinfo --
20852 516bb87ea72e 2015 Code Sprint: remove class of function from docstring for all m-files. Rik <rik@octave.org> parents: 20715 diff changeset	27 ## @deftypefn {} {@var{q} =} quadgk (@var{f}, @var{a}, @var{b})
516bb87ea72e 2015 Code Sprint: remove class of function from docstring for all m-files. Rik <rik@octave.org> parents: 20715 diff changeset	28 ## @deftypefnx {} {@var{q} =} quadgk (@var{f}, @var{a}, @var{b}, @var{abstol})
516bb87ea72e 2015 Code Sprint: remove class of function from docstring for all m-files. Rik <rik@octave.org> parents: 20715 diff changeset	29 ## @deftypefnx {} {@var{q} =} quadgk (@var{f}, @var{a}, @var{b}, @var{abstol}, @var{trace})
31063 12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	30 ## @deftypefnx {} {@var{q} =} quadgk (@var{f}, @var{a}, @var{b}, "@var{prop}", @var{val}, @dots{})
20852 516bb87ea72e 2015 Code Sprint: remove class of function from docstring for all m-files. Rik <rik@octave.org> parents: 20715 diff changeset	31 ## @deftypefnx {} {[@var{q}, @var{err}] =} quadgk (@dots{})
12612 16cca721117b doc: Update all documentation for chapter on Numerical Integration Rik <octave@nomad.inbox5.com> parents: 12595 diff changeset	32 ##
16cca721117b doc: Update all documentation for chapter on Numerical Integration Rik <octave@nomad.inbox5.com> parents: 12595 diff changeset	33 ## Numerically evaluate the integral of @var{f} from @var{a} to @var{b}
27951 14efec874546 quadgk.m: Fix direction of complex path integrals with "Waypoints" property (bug #57614) Rik <rik@octave.org> parents: 27928 diff changeset	34 ## using adaptive @nospell{Gauss-Kronrod} quadrature.
20158 7503499a252b doc: Update docstrings to have one sentence summary as first line. Rik <rik@octave.org> parents: 19833 diff changeset	35 ##
7503499a252b doc: Update docstrings to have one sentence summary as first line. Rik <rik@octave.org> parents: 19833 diff changeset	36 ## @var{f} is a function handle, inline function, or string containing the name
7503499a252b doc: Update docstrings to have one sentence summary as first line. Rik <rik@octave.org> parents: 19833 diff changeset	37 ## of the function to evaluate. The function @var{f} must be vectorized and
31063 12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	38 ## return a vector of output values when given a vector of input values (See
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	39 ## property @qcode{"ArrayValued"} for an exception to this rule).
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	40 ##
12612 16cca721117b doc: Update all documentation for chapter on Numerical Integration Rik <octave@nomad.inbox5.com> parents: 12595 diff changeset	41 ## @var{a} and @var{b} are the lower and upper limits of integration. Either
20158 7503499a252b doc: Update docstrings to have one sentence summary as first line. Rik <rik@octave.org> parents: 19833 diff changeset	42 ## or both limits may be infinite or contain weak end singularities. Variable
7503499a252b doc: Update docstrings to have one sentence summary as first line. Rik <rik@octave.org> parents: 19833 diff changeset	43 ## transformation will be used to treat any infinite intervals and weaken the
7503499a252b doc: Update docstrings to have one sentence summary as first line. Rik <rik@octave.org> parents: 19833 diff changeset	44 ## singularities. For example:
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	45 ##
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	46 ## @example
14327 4d917a6a858b doc: Use Octave coding conventions in @example blocks of docstrings. Rik <octave@nomad.inbox5.com> parents: 14138 diff changeset	47 ## quadgk (@@(x) 1 ./ (sqrt (x) .* (x + 1)), 0, Inf)
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	48 ## @end example
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	49 ##
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	50 ## @noindent
20158 7503499a252b doc: Update docstrings to have one sentence summary as first line. Rik <rik@octave.org> parents: 19833 diff changeset	51 ## Note that the formulation of the integrand uses the element-by-element
7503499a252b doc: Update docstrings to have one sentence summary as first line. Rik <rik@octave.org> parents: 19833 diff changeset	52 ## operator @code{./} and all user functions to @code{quadgk} should do the
7503499a252b doc: Update docstrings to have one sentence summary as first line. Rik <rik@octave.org> parents: 19833 diff changeset	53 ## same.
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	54 ##
31068 451fb63a10a0 update integral to call quadgk for 'ArrayValued' integrations (bug #62468) Nicholas R. Jankowski <jankowski.nicholas@gmail.com> parents: 31063 diff changeset	55 ## The optional argument @var{abstol} defines the absolute tolerance used to stop
24102 c723faa56ab4 Add missing functions integral2.m and integral3.m (bug #52074). Nicholas R. Jankowski <jankowskin@asme.org> parents: 24066 diff changeset	56 ## the integration procedure. The default value is 1e-10 (1e-5 for single).
12642 f96b9b9f141b doc: Periodic grammarcheck and spellcheck of documentation. Rik <octave@nomad.inbox5.com> parents: 12612 diff changeset	57 ##
20158 7503499a252b doc: Update docstrings to have one sentence summary as first line. Rik <rik@octave.org> parents: 19833 diff changeset	58 ## The algorithm used by @code{quadgk} involves subdividing the integration
7503499a252b doc: Update docstrings to have one sentence summary as first line. Rik <rik@octave.org> parents: 19833 diff changeset	59 ## interval and evaluating each subinterval. If @var{trace} is true then after
7503499a252b doc: Update docstrings to have one sentence summary as first line. Rik <rik@octave.org> parents: 19833 diff changeset	60 ## computing each of these partial integrals display: (1) the number of
7503499a252b doc: Update docstrings to have one sentence summary as first line. Rik <rik@octave.org> parents: 19833 diff changeset	61 ## subintervals at this step, (2) the current estimate of the error @var{err},
7503499a252b doc: Update docstrings to have one sentence summary as first line. Rik <rik@octave.org> parents: 19833 diff changeset	62 ## (3) the current estimate for the integral @var{q}.
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	63 ##
31063 12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	64 ## The behavior of the algorithm can be configured by passing arguments to
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	65 ## @code{quadgk} as pairs @qcode{"@var{prop}", @var{val}}. Valid properties
27951 14efec874546 quadgk.m: Fix direction of complex path integrals with "Waypoints" property (bug #57614) Rik <rik@octave.org> parents: 27928 diff changeset	66 ## are
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	67 ##
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	68 ## @table @code
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	69 ## @item AbsTol
12612 16cca721117b doc: Update all documentation for chapter on Numerical Integration Rik <octave@nomad.inbox5.com> parents: 12595 diff changeset	70 ## Define the absolute error tolerance for the quadrature. The default
20699 e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	71 ## absolute tolerance is 1e-10 (1e-5 for single).
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	72 ##
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	73 ## @item RelTol
12612 16cca721117b doc: Update all documentation for chapter on Numerical Integration Rik <octave@nomad.inbox5.com> parents: 12595 diff changeset	74 ## Define the relative error tolerance for the quadrature. The default
20699 e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	75 ## relative tolerance is 1e-6 (1e-4 for single).
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	76 ##
31063 12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	77 ## @item ArrayValued
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	78 ## When set to true, the function @var{f} produces an array output for a scalar
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	79 ## input. The default is false which requires that @var{f} produce an output
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	80 ## that is the same size as the input. For example,
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	81 ##
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	82 ## @example
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	83 ## quadgk (@@(x) x .^ (1:5), 0, 2, "ArrayValued", 1)
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	84 ## @end example
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	85 ##
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	86 ## will integrate @code{[x.^1, x.^2, x.^3, x.^4, x.^5]} in one function call
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	87 ## rather than having to repeatedly define a single anonymous function and
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	88 ## use a normal invocation of @code{quadgk}.
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	89 ##
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	90 ## @item WayPoints
31063 12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	91 ## Specify points which will become endpoints for subintervals in the
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	92 ## algorithm which can result in significantly improved estimation of the error
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	93 ## in the integral, faster computation, or both. It can be useful to specify
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	94 ## more subintervals around a region where the integrand is rapidly changing or
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	95 ## to flag locations where there is a discontinuity in the first derivative
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	96 ## of the function. For example, the signum function has a discontinuity at
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	97 ## @code{x == 0} and by specifying a waypoint
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	98 ##
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	99 ## @example
31063 12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	100 ## quadgk (@@(x) sign (x), -0.5, 1, "Waypoints", [0])
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	101 ## @end example
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	102 ##
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	103 ## @noindent
31063 12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	104 ## the error bound is reduced from 4e-7 to 1e-13.
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	105 ##
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	106 ## If the function has @strong{singularities} within the region of integration
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	107 ## those should not be addressed with waypoints. Instead, the overall integral
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	108 ## should be decomposed into a sum of several smaller integrals such that the
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	109 ## singularity occurs as one of the bounds of integration in the call to
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	110 ## @code{quadgk}.
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	111 ##
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	112 ## If any of the waypoints are complex then contour integration is performed as
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	113 ## documented below.
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	114 ##
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	115 ## @item MaxIntervalCount
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	116 ## @code{quadgk} initially subdivides the interval on which to perform the
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	117 ## quadrature into 10 intervals or, if WayPoints are given, at each waypoint.
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	118 ## Subintervals that have an unacceptable error are subdivided and
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	119 ## re-evaluated. If the number of subintervals exceeds 650 subintervals at any
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	120 ## point then a poor convergence is signaled and the current estimate of the
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	121 ## integral is returned. The property @qcode{"MaxIntervalCount"} can be used
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	122 ## to alter the number of subintervals that can exist before exiting.
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	123 ##
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	124 ## @item Trace
20158 7503499a252b doc: Update docstrings to have one sentence summary as first line. Rik <rik@octave.org> parents: 19833 diff changeset	125 ## If logically true @code{quadgk} prints information on the convergence of the
7503499a252b doc: Update docstrings to have one sentence summary as first line. Rik <rik@octave.org> parents: 19833 diff changeset	126 ## quadrature at each iteration.
12612 16cca721117b doc: Update all documentation for chapter on Numerical Integration Rik <octave@nomad.inbox5.com> parents: 12595 diff changeset	127 ## @end table
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	128 ##
12612 16cca721117b doc: Update all documentation for chapter on Numerical Integration Rik <octave@nomad.inbox5.com> parents: 12595 diff changeset	129 ## If any of @var{a}, @var{b}, or @var{waypoints} is complex then the
31063 12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	130 ## quadrature is treated as a contour integral along a piecewise linear
27951 14efec874546 quadgk.m: Fix direction of complex path integrals with "Waypoints" property (bug #57614) Rik <rik@octave.org> parents: 27928 diff changeset	131 ## path defined by
14efec874546 quadgk.m: Fix direction of complex path integrals with "Waypoints" property (bug #57614) Rik <rik@octave.org> parents: 27928 diff changeset	132 ## @code{[@var{a}, @var{waypoints}(1), @var{waypoints}(2), @dots{}, @var{b}]}.
14efec874546 quadgk.m: Fix direction of complex path integrals with "Waypoints" property (bug #57614) Rik <rik@octave.org> parents: 27928 diff changeset	133 ## In this case the integral is assumed to have no edge singularities. For
14efec874546 quadgk.m: Fix direction of complex path integrals with "Waypoints" property (bug #57614) Rik <rik@octave.org> parents: 27928 diff changeset	134 ## example,
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	135 ##
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	136 ## @example
9051 1bf0ce0930be Grammar check TexInfo in all .m files Rik <rdrider0-list@yahoo.com> parents: 8920 diff changeset	137 ## @group
8516 e2a179415bac doc fixes John W. Eaton <jwe@octave.org> parents: 8507 diff changeset	138 ## quadgk (@@(z) log (z), 1+1i, 1+1i, "WayPoints",
27951 14efec874546 quadgk.m: Fix direction of complex path integrals with "Waypoints" property (bug #57614) Rik <rik@octave.org> parents: 27928 diff changeset	139 ## [-1+1i, -1-1i, +1-1i])
9051 1bf0ce0930be Grammar check TexInfo in all .m files Rik <rdrider0-list@yahoo.com> parents: 8920 diff changeset	140 ## @end group
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	141 ## @end example
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	142 ##
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	143 ## @noindent
20158 7503499a252b doc: Update docstrings to have one sentence summary as first line. Rik <rik@octave.org> parents: 19833 diff changeset	144 ## integrates @code{log (z)} along the square defined by
27951 14efec874546 quadgk.m: Fix direction of complex path integrals with "Waypoints" property (bug #57614) Rik <rik@octave.org> parents: 27928 diff changeset	145 ## @code{[1+1i, -1+1i, -1-1i, +1-1i]}.
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	146 ##
12642 f96b9b9f141b doc: Periodic grammarcheck and spellcheck of documentation. Rik <octave@nomad.inbox5.com> parents: 12612 diff changeset	147 ## The result of the integration is returned in @var{q}.
20158 7503499a252b doc: Update docstrings to have one sentence summary as first line. Rik <rik@octave.org> parents: 19833 diff changeset	148 ##
12642 f96b9b9f141b doc: Periodic grammarcheck and spellcheck of documentation. Rik <octave@nomad.inbox5.com> parents: 12612 diff changeset	149 ## @var{err} is an approximate bound on the error in the integral
27951 14efec874546 quadgk.m: Fix direction of complex path integrals with "Waypoints" property (bug #57614) Rik <rik@octave.org> parents: 27928 diff changeset	150 ## @w{@code{abs (@var{q} - @var{I})}}, where @var{I} is the exact value of the
31000 1fcfe9443917 quadgk.m: Don't warn about failing to converge if second output ERR is requested (bug #62412) Rik <rik@octave.org> parents: 30569 diff changeset	151 ## integral. If the adaptive integration did not converge, the value of
1fcfe9443917 quadgk.m: Don't warn about failing to converge if second output ERR is requested (bug #62412) Rik <rik@octave.org> parents: 30569 diff changeset	152 ## @var{err} will be larger than the requested tolerance. If only a single
1fcfe9443917 quadgk.m: Don't warn about failing to converge if second output ERR is requested (bug #62412) Rik <rik@octave.org> parents: 30569 diff changeset	153 ## output is requested then a warning will be emitted when the requested
1fcfe9443917 quadgk.m: Don't warn about failing to converge if second output ERR is requested (bug #62412) Rik <rik@octave.org> parents: 30569 diff changeset	154 ## tolerance is not met. If the second output @var{err} is requested then no
1fcfe9443917 quadgk.m: Don't warn about failing to converge if second output ERR is requested (bug #62412) Rik <rik@octave.org> parents: 30569 diff changeset	155 ## warning is issued and it is the responsibility of the programmer to inspect
1fcfe9443917 quadgk.m: Don't warn about failing to converge if second output ERR is requested (bug #62412) Rik <rik@octave.org> parents: 30569 diff changeset	156 ## and determine whether the results are satisfactory.
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	157 ##
20158 7503499a252b doc: Update docstrings to have one sentence summary as first line. Rik <rik@octave.org> parents: 19833 diff changeset	158 ## Reference: @nospell{L.F. Shampine},
7503499a252b doc: Update docstrings to have one sentence summary as first line. Rik <rik@octave.org> parents: 19833 diff changeset	159 ## @cite{"Vectorized adaptive quadrature in @sc{matlab}"}, Journal of
27805 5a6a19a4e3da doc: Use Texinfo non-sentence ending periods in citations. Rik <rik@octave.org> parents: 26949 diff changeset	160 ## Computational and Applied Mathematics, pp.@: 131--140, Vol 211, Issue 2,
20158 7503499a252b doc: Update docstrings to have one sentence summary as first line. Rik <rik@octave.org> parents: 19833 diff changeset	161 ## Feb 2008.
7503499a252b doc: Update docstrings to have one sentence summary as first line. Rik <rik@octave.org> parents: 19833 diff changeset	162 ##
24208 eec262017c6a maint: Strip trailing whitespace from source files. John W. Eaton <jwe@octave.org> parents: 24102 diff changeset	163 ## @seealso{quad, quadv, quadl, quadcc, trapz, dblquad, triplequad, integral,
31063 12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	164 ## integral2, integral3}
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	165 ## @end deftypefn
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	166
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	167 function [q, err] = quadgk (f, a, b, varargin)
20699 e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	168
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	169 if (nargin < 3)
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	170 print_usage ();
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	171 endif
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	172
20699 e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	173 abstol = [];
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	174 reltol = [];
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	175 waypoints = [];
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	176 maxint = 650;
31063 12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	177 arrayvalued = false;
20699 e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	178 trace = false;
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	179
20699 e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	180 ## Parse options if present.
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	181 if (nargin > 3)
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	182 if (! ischar (varargin{1}))
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	183 if (! isempty (varargin{1}))
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	184 abstol = varargin{1};
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	185 reltol = 0;
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	186 endif
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	187 if (nargin > 4)
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	188 trace = varargin{2};
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	189 endif
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	190 if (nargin > 5)
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	191 error ("quadgk: can not pass additional arguments to user function");
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	192 endif
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	193 else
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	194 if (mod (nargin - 3, 2) != 0)
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	195 error ("quadgk: property/value options must occur in pairs");
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	196 endif
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	197
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	198 idx = 1;
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	199 while (idx < nargin - 3)
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	200 if (! ischar (varargin{idx}))
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	201 error ("quadgk: property PROP must be a string");
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	202 endif
31063 12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	203 prop = varargin{idx++};
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	204 switch (tolower (prop))
20699 e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	205 case "reltol"
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	206 reltol = varargin{idx++};
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	207 case "abstol"
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	208 abstol = varargin{idx++};
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	209 case "waypoints"
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	210 waypoints = varargin{idx++}(:);
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	211 case "maxintervalcount"
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	212 maxint = varargin{idx++};
31063 12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	213 case "arrayvalued"
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	214 arrayvalued = varargin{idx++};
20699 e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	215 case "trace"
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	216 trace = varargin{idx++};
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	217 otherwise
31063 12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	218 error ("quadgk: unknown property '%s'", prop);
20699 e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	219 endswitch
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	220 endwhile
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	221 endif
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	222 endif
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	223
31063 12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	224 reverse = 1;
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	225 contour = iscomplex (a) \|\| iscomplex (b) \|\| iscomplex (waypoints);
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	226 if ((b < a) && ! contour)
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	227 ## Reverse integration
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	228 [b, a] = deal (a, b);
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	229 waypoints = sort (waypoints(waypoints > a & waypoints < b));
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	230 reverse = -1;
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	231 endif
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	232
20699 e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	233 issingle = (isa (a, "single") \|\| isa (b, "single")
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	234 \|\| isa (waypoints, "single"));
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	235
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	236 if (isempty (abstol))
20715 5b7643257978 Remove trailing whitespace at end of lines. Rik <rik@octave.org> parents: 20699 diff changeset	237 abstol = ifelse (issingle, 1e-5, 1e-10);
20699 e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	238 elseif (! isscalar (abstol) \|\| abstol < 0)
23324 336f2e8331fa quadgk.m: Correct error messages which point to quadv (bug #50604). Nicholas R. Jankowski <jankowskin@asme.org> parents: 23219 diff changeset	239 error ("quadgk: ABSTOL must be a scalar >=0");
20699 e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	240 endif
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	241
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	242 if (isempty (reltol))
20715 5b7643257978 Remove trailing whitespace at end of lines. Rik <rik@octave.org> parents: 20699 diff changeset	243 reltol = ifelse (issingle, 1e-4, 1e-6);
20699 e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	244 elseif (! isscalar (reltol) \|\| reltol < 0)
23324 336f2e8331fa quadgk.m: Correct error messages which point to quadv (bug #50604). Nicholas R. Jankowski <jankowskin@asme.org> parents: 23219 diff changeset	245 error ("quadgk: RELTOL must be a scalar >=0");
20699 e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	246 endif
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	247
31063 12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	248 ## FIXME: No validation of inputs MaxIntervalCount, Waypoints, ArrayValued,
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	249 ## Trace.
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	250
20699 e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	251 ## Convert function given as a string to a function handle
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	252 if (ischar (f))
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	253 f = @(x) feval (f, x);
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	254 endif
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	255
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	256 ## Use variable substitution to weaken endpoint singularities and
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	257 ## to perform integration with endpoints at infinity.
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	258 ## No transform for contour integrals.
31063 12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	259 if (contour)
20699 e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	260 ## contour integral, no transform
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	261 subs = [a; waypoints; b];
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	262 h = sum (abs (diff (subs)));
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	263 trans = @(t) t;
31063 12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	264 ## Ensure f is always vectorized even if specified as, e.g., f = @(x) 1;
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	265 f = @(t) f (t) + 0*t;
20699 e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	266 elseif (isinf (a) && isinf (b))
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	267 ## Standard infinite to finite integral transformation.
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	268 ## \int_{-\infinity_^\infinity f(x) dx = \int_-1^1 f (g(t)) g'(t) dt
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	269 ## where
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	270 ## g(t) = t / (1 - t^2)
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	271 ## g'(t) = (1 + t^2) / (1 - t^2) ^ 2
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	272 ## waypoint transform is then
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	273 ## t = (2 * g(t)) ./ (1 + sqrt(1 + 4 * g(t) .^ 2))
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	274 if (! isempty (waypoints))
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	275 trans = @(x) (2 * x) ./ (1 + sqrt (1 + 4 * x .^ 2));
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	276 subs = [-1; trans(waypoints); 1];
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	277 else
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	278 subs = linspace (-1, 1, 11)';
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	279 endif
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	280 h = 2;
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	281 trans = @(t) t ./ (1 - t.^2);
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	282 f = @(t) f (t ./ (1 - t .^ 2)) .* (1 + t .^ 2) ./ ((1 - t .^ 2) .^ 2);
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	283 elseif (isinf (a))
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	284 ## Formula defined in Shampine paper as two separate steps.
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	285 ## One to weaken singularity at finite end, then a second to transform to
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	286 ## a finite interval. The singularity weakening transform is
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	287 ## \int_{-\infinity}^b f(x) dx =
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	288 ## - \int_{-\infinity}^0 f (b - t^2) 2 t dt
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	289 ## (note minus sign) and the finite interval transform is
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	290 ## \int_{-\infinity}^0 f(b - t^2) 2 t dt =
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	291 ## \int_{-1}^0 f (b - g(s) ^ 2) 2 g(s) g'(s) ds
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	292 ## where
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	293 ## g(s) = s / (1 + s)
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	294 ## g'(s) = 1 / (1 + s) ^ 2
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	295 ## waypoint transform is then
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	296 ## t = sqrt (b - x)
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	297 ## s = - t / (t + 1)
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	298 if (! isempty (waypoints))
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	299 tmp = sqrt (b - waypoints);
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	300 trans = @(x) - x ./ (x + 1);
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	301 subs = [-1; trans(tmp); 0];
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	302 else
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	303 subs = linspace (-1, 0, 11)';
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	304 endif
20699 e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	305 h = 1;
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	306 trans = @(t) b - (t ./ (1 + t)).^2;
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	307 f = @(s) - 2 * s .* f (b - (s ./ (1 + s)) .^ 2) ./ ((1 + s) .^ 3);
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	308 elseif (isinf (b))
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	309 ## Formula defined in Shampine paper as two separate steps.
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	310 ## One to weaken singularity at finite end, then a second to transform to
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	311 ## a finite interval. The singularity weakening transform is
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	312 ## \int_a^\infinity f(x) dx = \int_0^\infinity f (a + t^2) 2 t dt
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	313 ## and the finite interval transform is
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	314 ## \int_0^\infinity f(a + t^2) 2 t dt =
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	315 ## \int_0^1 f (a + g(s) ^ 2) 2 g(s) g'(s) ds
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	316 ## where
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	317 ## g(s) = s / (1 - s)
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	318 ## g'(s) = 1 / (1 - s) ^ 2
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	319 ## waypoint transform is then
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	320 ## t = sqrt (x - a)
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	321 ## s = t / (t + 1)
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	322 if (! isempty (waypoints))
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	323 tmp = sqrt (waypoints - a);
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	324 trans = @(x) x ./ (x + 1);
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	325 subs = [0; trans(tmp); 1];
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	326 else
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	327 subs = linspace (0, 1, 11)';
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	328 endif
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	329 h = 1;
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	330 trans = @(t) a + (t ./ (1 - t)).^2;
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	331 f = @(s) 2 * s .* f (a + (s ./ (1 - s)) .^ 2) ./ ((1 - s) .^ 3);
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	332 else
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	333 ## Davis, Rabinowitz, "Methods of Numerical Integration" p441 2ed.
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	334 ## Presented in section 5 of the Shampine paper as
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	335 ## g(t) = ((b - a) / 2) * (t / 2 * (3 - t^2)) + (b + a) / 2
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	336 ## g'(t) = ((b-a)/4) * (3 - 3t^2);
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	337 ## waypoint transform can then be found by solving for t with
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	338 ## Maxima (solve (c + 3*t - 3^3, t);). This gives 3 roots, two of
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	339 ## which are complex for values between a and b and so can be ignored.
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	340 ## The third is
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	341 ## c = (-4x + 2(b+a)) / (b-a);
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	342 ## k = ((sqrt(c^2 - 4) + c)/2)^(1/3);
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	343 ## t = (sqrt(3)* 1i * (1 - k^2) - (1 + k^2)) / 2 / k;
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	344 if (! isempty (waypoints))
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	345 trans = @__quadgk_finite_waypoint__;
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	346 subs = [-1; trans(waypoints, a, b); 1];
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	347 else
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	348 subs = linspace (-1, 1, 11)';
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	349 endif
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	350 h = 2;
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	351 trans = @(t) ((b - a) ./ 4) * t .* (3 - t.^2) + (b + a) ./ 2;
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	352 f = @(t) f((b - a) ./ 4 .* t .* (3 - t.^2) + (b + a) ./ 2) .* ...
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	353 3 .* (b - a) ./ 4 .* (1 - t.^2);
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	354 endif
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	355
20699 e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	356 ## Split interval into at least 10 subinterval with a 15 point
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	357 ## Gauss-Kronrod rule giving a minimum of 150 function evaluations.
31063 12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	358 while (numel (subs) < 11)
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	359 subs = [subs.' ; subs(1:end-1).' + diff(subs.') ./ 2, NaN](:)(1:end-1);
20699 e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	360 endwhile
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	361 subs = [subs(1:end-1), subs(2:end)];
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	362
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	363 warn_id = "Octave:quadgk:warning-termination";
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	364
31063 12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	365 if (! arrayvalued)
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	366 ## Initial evaluation of the integrand on the subintervals.
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	367 [q_subs, q_errs] = __quadgk_eval__ (f, subs, trans);
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	368 q0 = sum (q_subs);
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	369 err0 = sum (q_errs);
22474 ed1722f70fad Fix early exit from quadgk for infinite intervals (bug #37613). Marco Caliari <marco.caliari@univr.it> parents: 22323 diff changeset	370
31063 12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	371 first = true;
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	372 while (true)
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	373 ## Quit if any evaluations are not finite (Inf or NaN).
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	374 if (any (! isfinite (q_subs)))
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	375 warning (warn_id, "quadgk: non-finite integrand encountered");
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	376 q = q0;
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	377 err = err0;
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	378 break;
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	379 endif
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	380
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	381 tol = max (abstol, reltol .* abs (q0));
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	382
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	383 ## If the global error estimate is met then exit.
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	384 if (err0 <= tol)
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	385 q = q0;
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	386 err = err0;
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	387 break;
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	388 endif
20699 e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	389
31063 12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	390 ## Accept the subintervals that meet the convergence criteria.
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	391 idx = find (abs (q_errs) < tol .* abs (diff (subs, 1, 2)) ./ h);
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	392 if (first)
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	393 q = sum (q_subs(idx));
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	394 err = sum (q_errs(idx));
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	395 first = false;
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	396 else
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	397 q0 = q + sum (q_subs);
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	398 err0 = err + sum (q_errs);
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	399 q += sum (q_subs(idx));
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	400 err += sum (q_errs(idx));
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	401 endif
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	402 subs(idx,:) = [];
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	403
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	404 ## If no remaining subintervals then exit.
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	405 if (isempty (subs))
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	406 break;
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	407 endif
20699 e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	408
31063 12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	409 if (trace)
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	410 disp ([rows(subs), err, q0]);
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	411 endif
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	412
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	413 ## Split remaining subintervals in two
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	414 mid = (subs(:,1) + subs(:,2)) ./ 2;
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	415 subs = [subs(:,1), mid; mid, subs(:,2)];
20699 e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	416
31063 12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	417 ## If the maximum subinterval count is met, then
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	418 ## accept remaining subinterval and exit.
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	419 if (rows (subs) > maxint)
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	420 warning (warn_id, "quadgk: maximum interval count (%d) exceeded", maxint);
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	421 q += sum (q_subs);
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	422 err += sum (q_errs);
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	423 break;
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	424 endif
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	425
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	426 ## Evaluation of the integrand on the remaining subintervals
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	427 [q_subs, q_errs] = __quadgk_eval__ (f, subs, trans);
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	428 endwhile
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	429
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	430 if (nargout < 2 && err > max (abstol, reltol * abs (q)))
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	431 warning (warn_id,
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	432 "quadgk: Error tolerance not met. Estimated error %g", err);
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	433 endif
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	434
31063 12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	435 ## Reverse integral if necessary.
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	436 q = reverse * q;
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	437
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	438 else
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	439 ## f is array-valued
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	440 sz = size (f (subs(1)));
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	441
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	442 ## Initial evaluation of the integrand on the subintervals
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	443 [q_subs, q_errs] = __quadgk_eval_array__ (f, subs, trans, prod (sz));
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	444 q0 = sum (q_subs, 1);
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	445 err0 = sum (q_errs, 1);
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	446
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	447 first = true;
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	448 while (true)
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	449 ## Quit if any evaluations are not finite (Inf or NaN).
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	450 if (any (! isfinite (q_subs)(:)))
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	451 warning (warn_id, "quadgk: non-finite integrand encountered");
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	452 q = q0;
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	453 err = err0;
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	454 break;
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	455 endif
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	456
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	457 tol = max (abstol, reltol .* abs (q0));
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	458
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	459 ## If the global error estimate is met then exit
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	460 if (err0 <= tol)
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	461 q = q0;
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	462 err = err0;
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	463 break;
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	464 endif
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	465
31063 12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	466 ## Accept subintervals that meet the convergence criteria in all entries.
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	467 idx = find (all (abs (q_errs) < tol .* abs (diff (subs, 1, 2)) ./ h, 2));
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	468 if (first)
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	469 q = sum (q_subs(idx,:), 1);
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	470 err = sum (q_errs(idx,:), 1);
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	471 first = false;
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	472 else
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	473 q0 = q + sum (q_subs, 1);
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	474 err0 = err + sum (q_errs, 1);
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	475 q += sum (q_subs(idx,:), 1);
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	476 err += sum (q_errs(idx,:), 1);
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	477 endif
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	478 subs(idx,:) = [];
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	479
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	480 ## If no remaining subintervals exit
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	481 if (isempty (subs))
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	482 break;
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	483 endif
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	484
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	485 if (trace)
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	486 disp ([rows(subs), err(1, 1), q0(1, 1)]); # print only first entry
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	487 endif
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	488
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	489 ## Split remaining subintervals in two
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	490 mid = (subs(:,1) + subs(:,2)) ./ 2;
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	491 subs = [subs(:,1), mid; mid, subs(:,2)];
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	492
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	493 ## If the maximum subinterval count is met accept remaining subinterval
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	494 ## and exit
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	495 if (rows (subs) > maxint)
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	496 warning (warn_id, "quadgk: maximum interval count (%d) exceeded", maxint);
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	497 q += sum (q_subs, 1);
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	498 err += sum (q_errs, 1);
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	499 break;
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	500 endif
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	501
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	502 ## Evaluation of the integrand on the remaining subintervals
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	503 [q_subs, q_errs] = __quadgk_eval_array__ (f, subs, trans, prod (sz));
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	504 endwhile
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	505
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	506 i = find (err > max (abstol, reltol * abs (q)), 1);
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	507 if (nargout < 2 && length (i) > 0)
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	508 ## like ind2sub, only as vector.
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	509 j = mod (floor ((i-1)./cumprod ([1 sz(1:end-1)])),sz)+1;
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	510 s = ["(" sprintf("%d,",j)(1:end-1) ")"];
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	511 warning (warn_id,
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	512 "quadgk: Error tolerance not met. First entry at index %s with estimated error %g", s, err(i));
20699 e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	513 endif
15579 5fb80374c881 quadgk.m: cleanup unwind_protect cruft Jordi Gutiérrez Hermoso <jordigh@octave.org> parents: 15578 diff changeset	514
31063 12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	515 q = reverse * reshape (q, sz);
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	516 err = reshape (err, sz);
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	517 endif
20699 e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	518
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	519 endfunction
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	520
31063 12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	521 function [q, err] = __quadgk_eval__ (f, subs, trans)
30384 363fb10055df maint: Style check m-files ahead of 7.1 release. Rik <rik@octave.org> parents: 29364 diff changeset	522
22474 ed1722f70fad Fix early exit from quadgk for infinite intervals (bug #37613). Marco Caliari <marco.caliari@univr.it> parents: 22323 diff changeset	523 ## A (15,7) point pair of Gauss-Kronrod quadrature rules.
20699 e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	524 ## The abscissa and weights are copied directly from dqk15w.f from quadpack.
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	525
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	526 persistent abscissa = [-0.9914553711208126e+00, -0.9491079123427585e+00, ...
10549 95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 9590 diff changeset	527 -0.8648644233597691e+00, -0.7415311855993944e+00, ...
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 9590 diff changeset	528 -0.5860872354676911e+00, -0.4058451513773972e+00, ...
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 9590 diff changeset	529 -0.2077849550078985e+00, 0.0000000000000000e+00, ...
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 9590 diff changeset	530 0.2077849550078985e+00, 0.4058451513773972e+00, ...
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 9590 diff changeset	531 0.5860872354676911e+00, 0.7415311855993944e+00, ...
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 9590 diff changeset	532 0.8648644233597691e+00, 0.9491079123427585e+00, ...
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 9590 diff changeset	533 0.9914553711208126e+00];
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	534
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	535 persistent weights15 = ...
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	536 diag ([0.2293532201052922e-01, 0.6309209262997855e-01, ...
10549 95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 9590 diff changeset	537 0.1047900103222502e+00, 0.1406532597155259e+00, ...
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 9590 diff changeset	538 0.1690047266392679e+00, 0.1903505780647854e+00, ...
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 9590 diff changeset	539 0.2044329400752989e+00, 0.2094821410847278e+00, ...
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 9590 diff changeset	540 0.2044329400752989e+00, 0.1903505780647854e+00, ...
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 9590 diff changeset	541 0.1690047266392679e+00, 0.1406532597155259e+00, ...
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 9590 diff changeset	542 0.1047900103222502e+00, 0.6309209262997855e-01, ...
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 9590 diff changeset	543 0.2293532201052922e-01]);
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	544
21568 3d60ed163b70 maint: Eliminate bad spacing around '='. Rik <rik@octave.org> parents: 20852 diff changeset	545 persistent weights7 = ...
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	546 diag ([0.1294849661688697e+00, 0.2797053914892767e+00, ...
10549 95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 9590 diff changeset	547 0.3818300505051889e+00, 0.4179591836734694e+00, ...
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 9590 diff changeset	548 0.3818300505051889e+00, 0.2797053914892767e+00, ...
95c3e38098bf Untabify .m scripts Rik <code@nomad.inbox5.com> parents: 9590 diff changeset	549 0.1294849661688697e+00]);
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	550
31063 12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	551 halfwidth = diff (subs, 1, 2) ./ 2;
20699 e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	552 center = sum (subs, 2) ./ 2;
22474 ed1722f70fad Fix early exit from quadgk for infinite intervals (bug #37613). Marco Caliari <marco.caliari@univr.it> parents: 22323 diff changeset	553 t = (halfwidth * abscissa) + center;
ed1722f70fad Fix early exit from quadgk for infinite intervals (bug #37613). Marco Caliari <marco.caliari@univr.it> parents: 22323 diff changeset	554 x = trans ([t(:,1), t(:,end)]);
ed1722f70fad Fix early exit from quadgk for infinite intervals (bug #37613). Marco Caliari <marco.caliari@univr.it> parents: 22323 diff changeset	555
ed1722f70fad Fix early exit from quadgk for infinite intervals (bug #37613). Marco Caliari <marco.caliari@univr.it> parents: 22323 diff changeset	556 y = reshape (f (t(:)), size (t));
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	557
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	558 ## This is faster than using bsxfun as the * operator can use a
20699 e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	559 ## single BLAS call, rather than rows (sub) calls to the @times function.
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	560 q = sum (y * weights15, 2) .* halfwidth;
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	561 err = abs (sum (y(:,2:2:end) * weights7, 2) .* halfwidth - q);
21758 ffad2baa90f7 maint: Use newlines to make code more readable. Rik <rik@octave.org> parents: 21568 diff changeset	562
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	563 endfunction
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	564
31063 12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	565 function [q, err] = __quadgk_eval_array__ (f, subs, trans, nel)
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	566
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	567 ## A (15,7) point pair of Gauss-Kronrod quadrature rules.
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	568 ## The abscissa and weights are copied directly from dqk15w.f from quadpack.
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	569
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	570 persistent abscissa = [-0.9914553711208126e+00, -0.9491079123427585e+00, ...
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	571 -0.8648644233597691e+00, -0.7415311855993944e+00, ...
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	572 -0.5860872354676911e+00, -0.4058451513773972e+00, ...
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	573 -0.2077849550078985e+00, 0.0000000000000000e+00, ...
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	574 0.2077849550078985e+00, 0.4058451513773972e+00, ...
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	575 0.5860872354676911e+00, 0.7415311855993944e+00, ...
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	576 0.8648644233597691e+00, 0.9491079123427585e+00, ...
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	577 0.9914553711208126e+00];
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	578
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	579 persistent weights15 = ...
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	580 [0.2293532201052922e-01, 0.6309209262997855e-01, ...
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	581 0.1047900103222502e+00, 0.1406532597155259e+00, ...
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	582 0.1690047266392679e+00, 0.1903505780647854e+00, ...
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	583 0.2044329400752989e+00, 0.2094821410847278e+00, ...
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	584 0.2044329400752989e+00, 0.1903505780647854e+00, ...
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	585 0.1690047266392679e+00, 0.1406532597155259e+00, ...
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	586 0.1047900103222502e+00, 0.6309209262997855e-01, ...
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	587 0.2293532201052922e-01];
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	588
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	589 persistent weights7 = ...
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	590 [0.1294849661688697e+00, 0.2797053914892767e+00, ...
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	591 0.3818300505051889e+00, 0.4179591836734694e+00, ...
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	592 0.3818300505051889e+00, 0.2797053914892767e+00, ...
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	593 0.1294849661688697e+00];
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	594
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	595 halfwidth = diff (subs, 1, 2) ./ 2;
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	596 center = sum (subs, 2) ./ 2;
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	597 t = (halfwidth * abscissa) + center;
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	598 x = trans ([t(:,1), t(:,end)]);
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	599
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	600 y = zeros (nel, columns(t), rows(t));
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	601 for i = 1:rows (t)
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	602 for j = 1:columns(t)
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	603 y(:,j,i) = f (t(i,j))(:);
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	604 endfor
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	605 endfor
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	606 y = permute (y, [2 3 1]);
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	607
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	608 q = reshape (weights15 * y(:,:), [rows(t), nel]) .* halfwidth;
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	609 err = abs (reshape (weights7 * y(2:2:end,:), rows (t), nel) .* halfwidth - q);
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	610
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	611 endfunction
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	612
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	613 function t = __quadgk_finite_waypoint__ (x, a, b)
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	614 c = (-4 .* x + 2.* (b + a)) ./ (b - a);
14868 5d3a684236b0 maint: Use Octave coding conventions for cuddling parentheses in scripts directory Rik <octave@nomad.inbox5.com> parents: 14363 diff changeset	615 k = ((sqrt (c .^ 2 - 4) + c) ./ 2) .^ (1/3);
28917 0de38a6ef693 maint: Use Octave convention of space after function name in scripts dir. Rik <rik@octave.org> parents: 28665 diff changeset	616 t = real ((sqrt (3) .* 1i * (1 - k .^ 2) - (1 + k .^ 2)) ./ 2 ./ k);
7771 680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	617 endfunction
680631e787aa Add quadv, quadgk, dblquad and triplequad functions David Bateman <dbateman@free.fr> parents: diff changeset	618
14363 f3d52523cde1 Use Octave coding conventions in all m-file %!test blocks Rik <octave@nomad.inbox5.com> parents: 14359 diff changeset	619
24066 d9b0d8ae734f Update tolerances in BIST tests for quad functions. Rik <rik@octave.org> parents: 23960 diff changeset	620 %!assert (quadgk (@sin,-pi,pi), 0, 1e-10)
31000 1fcfe9443917 quadgk.m: Don't warn about failing to converge if second output ERR is requested (bug #62412) Rik <rik@octave.org> parents: 30569 diff changeset	621 %!test
1fcfe9443917 quadgk.m: Don't warn about failing to converge if second output ERR is requested (bug #62412) Rik <rik@octave.org> parents: 30569 diff changeset	622 %! warning ("off", "Octave:legacy-function", "local");
1fcfe9443917 quadgk.m: Don't warn about failing to converge if second output ERR is requested (bug #62412) Rik <rik@octave.org> parents: 30569 diff changeset	623 %! assert (quadgk (inline ("sin"), -pi, pi), 0, 1e-10);
24066 d9b0d8ae734f Update tolerances in BIST tests for quad functions. Rik <rik@octave.org> parents: 23960 diff changeset	624 %!assert (quadgk ("sin",-pi,pi), 0, 1e-10)
d9b0d8ae734f Update tolerances in BIST tests for quad functions. Rik <rik@octave.org> parents: 23960 diff changeset	625 %!assert (quadgk (@sin,-pi,pi, "WayPoints", 0, "MaxIntervalCount", 100,
d9b0d8ae734f Update tolerances in BIST tests for quad functions. Rik <rik@octave.org> parents: 23960 diff changeset	626 %! "RelTol", 1e-3, "AbsTol", 1e-6, "trace", false), 0, 1e-6)
d9b0d8ae734f Update tolerances in BIST tests for quad functions. Rik <rik@octave.org> parents: 23960 diff changeset	627 %!assert (quadgk (@sin,-pi,pi, 1e-6, false), 0, 1e-6)
26949 04e5cb5e2cb3 update bug status in tests John W. Eaton <jwe@octave.org> parents: 26599 diff changeset	628 %!assert <*51867> (quadgk (@(x) x, 0, 0), 0, 0)
14363 f3d52523cde1 Use Octave coding conventions in all m-file %!test blocks Rik <octave@nomad.inbox5.com> parents: 14359 diff changeset	629
24066 d9b0d8ae734f Update tolerances in BIST tests for quad functions. Rik <rik@octave.org> parents: 23960 diff changeset	630 %!assert (quadgk (@sin,-pi,0), -2, 1e-10)
d9b0d8ae734f Update tolerances in BIST tests for quad functions. Rik <rik@octave.org> parents: 23960 diff changeset	631 %!assert (quadgk (@sin,0,pi), 2, 1e-10)
d9b0d8ae734f Update tolerances in BIST tests for quad functions. Rik <rik@octave.org> parents: 23960 diff changeset	632 %!assert (quadgk (@(x) 1./sqrt (x),0,1), 2, 1e-10)
d9b0d8ae734f Update tolerances in BIST tests for quad functions. Rik <rik@octave.org> parents: 23960 diff changeset	633 %!assert (quadgk (@(x) abs (1 - x.^2),0,2, "Waypoints", 1), 2, 1e-10)
d9b0d8ae734f Update tolerances in BIST tests for quad functions. Rik <rik@octave.org> parents: 23960 diff changeset	634 %!assert (quadgk (@(x) 1./(sqrt (x) .* (x+1)),0,Inf), pi, 1e-10)
28917 0de38a6ef693 maint: Use Octave convention of space after function name in scripts dir. Rik <rik@octave.org> parents: 28665 diff changeset	635 %!assert <*57614> (quadgk (@(z) exp (z)./z, 1, 1,
27951 14efec874546 quadgk.m: Fix direction of complex path integrals with "Waypoints" property (bug #57614) Rik <rik@octave.org> parents: 27928 diff changeset	636 %! "Waypoints", [1+i, -1+i, -1-i, 1-i]),
28665 dc80e087df4b maint: Strip trailing spaces from files. Rik <rik@octave.org> parents: 27951 diff changeset	637 %! complex (0, 2*pi), 1e-10)
28917 0de38a6ef693 maint: Use Octave convention of space after function name in scripts dir. Rik <rik@octave.org> parents: 28665 diff changeset	638 %!assert <*57614> (quadgk (@(z) exp (z)./z, 1, 1,
27951 14efec874546 quadgk.m: Fix direction of complex path integrals with "Waypoints" property (bug #57614) Rik <rik@octave.org> parents: 27928 diff changeset	639 %! "Waypoints", [1-i, -1-i, -1+i, 1+i]),
28665 dc80e087df4b maint: Strip trailing spaces from files. Rik <rik@octave.org> parents: 27951 diff changeset	640 %! complex (0, -2*pi), 1e-10)
24066 d9b0d8ae734f Update tolerances in BIST tests for quad functions. Rik <rik@octave.org> parents: 23960 diff changeset	641 %!assert (quadgk (@(z) log (z),1+1i,1+1i, "WayPoints", [1-1i, -1,-1i, -1+1i]),
27951 14efec874546 quadgk.m: Fix direction of complex path integrals with "Waypoints" property (bug #57614) Rik <rik@octave.org> parents: 27928 diff changeset	642 %! complex (0, pi), 1e-10)
24066 d9b0d8ae734f Update tolerances in BIST tests for quad functions. Rik <rik@octave.org> parents: 23960 diff changeset	643 %!assert (quadgk (@(x) exp (-x .^ 2),-Inf,Inf), sqrt (pi), -1e-6)
d9b0d8ae734f Update tolerances in BIST tests for quad functions. Rik <rik@octave.org> parents: 23960 diff changeset	644 %!assert (quadgk (@(x) exp (-x .^ 2),-Inf,0), sqrt (pi)/2, -1e-6)
22474 ed1722f70fad Fix early exit from quadgk for infinite intervals (bug #37613). Marco Caliari <marco.caliari@univr.it> parents: 22323 diff changeset	645 %!test
ed1722f70fad Fix early exit from quadgk for infinite intervals (bug #37613). Marco Caliari <marco.caliari@univr.it> parents: 22323 diff changeset	646 %! f = @(x) x .^ 5 .* exp (-x) .* sin (x);
27951 14efec874546 quadgk.m: Fix direction of complex path integrals with "Waypoints" property (bug #57614) Rik <rik@octave.org> parents: 27928 diff changeset	647 %! assert (quadgk (f, 0, Inf, "RelTol", 1e-8, "AbsTol", 1e-12), -15, -1e-8);
14363 f3d52523cde1 Use Octave coding conventions in all m-file %!test blocks Rik <octave@nomad.inbox5.com> parents: 14359 diff changeset	648
31063 12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	649 ## Test vector-valued functions
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	650 %!assert (quadgk (@(x) [(sin (x)), (sin (2 * x))], 0, pi, "arrayvalued", 1),
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	651 %! [2, 0], 1e-6)
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	652
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	653 ## Test matrix-valued functions
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	654 %!assert (quadgk (@(x) [ x,x,x; x,1./sqrt(x),x; x,x,x ], 0, 1, "arrayvalued",1),
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	655 %! [0.5,0.5,0.5; 0.5,2,0.5; 0.5,0.5,0.5], 15*1e-6);
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	656
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	657 ## Bug #62412
12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	658 %!warning <Error tolerance not met>
31000 1fcfe9443917 quadgk.m: Don't warn about failing to converge if second output ERR is requested (bug #62412) Rik <rik@octave.org> parents: 30569 diff changeset	659 %! f = @(t) -1 ./ t.^1.1;
31063 12f8fb75fc30 quadgk.m: Overhaul function to add new "ArrayValued" option (bug #62468) Michael Leitner <michael.leitner@frm2.tum.de>, Rik <rik@octave.org> parents: 31001 diff changeset	660 %! quadgk (f, 1, Inf);
31000 1fcfe9443917 quadgk.m: Don't warn about failing to converge if second output ERR is requested (bug #62412) Rik <rik@octave.org> parents: 30569 diff changeset	661
20699 e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	662 ## Test input validation
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	663 %!error quadgk (@sin)
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	664 %!error quadgk (@sin, 0)
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	665 %!error <can not pass additional arguments> quadgk (@sin, 0, 1, 1e-6, true, 4)
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	666 %!error <options must occur in pairs> quadgk (@sin, 0, 1, "DummyArg")
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	667 %!error <PROP must be a string> quadgk (@sin, 0, 1, "AbsTol", 1e-6, 2, 3)
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	668 %!error <unknown property 'foo'> quadgk (@sin, 0, 1, "foo", 3)
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	669 %!error <ABSTOL must be a scalar> quadgk (@sin, 0, 1, ones (2,2))
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	670 %!error <ABSTOL must be a scalar .=0> quadgk (@sin, 0, 1, -1)
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	671 %!error <RELTOL must be a scalar> quadgk (@sin, 0, 1, "RelTol", ones (2,2))
e15b7963746b quadgk.m: Overhaul function and increase Matlab compatibility. Rik <rik@octave.org> parents: 20617 diff changeset	672 %!error <RELTOL must be a scalar> quadgk (@sin, 0, 1, "RelTol", -1)