Mercurial > octave

annotate scripts/general/quadl.m @ 22323:bac0d6f07a3e

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
maint: Update copyright notices for 2016.
author John W. Eaton <jwe@octave.org>
date Wed, 17 Aug 2016 01:05:19 -0400
parents 516bb87ea72e
children 3a2b891d0b33 e9a0469dedd9
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
22323
bac0d6f07a3e maint: Update copyright notices for 2016.
John W. Eaton <jwe@octave.org>
parents: 20852
diff changeset
1 ## Copyright (C) 1998-2016 Walter Gautschi
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
2 ##
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
3 ## This file is part of Octave.
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
4 ##
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
5 ## Octave is free software; you can redistribute it and/or modify it
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
6 ## under the terms of the GNU General Public License as published by
7016
93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe]
jwe
parents: 5838
diff changeset
7 ## the Free Software Foundation; either version 3 of the License, or (at
93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe]
jwe
parents: 5838
diff changeset
8 ## your option) any later version.
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
9 ##
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
10 ## Octave is distributed in the hope that it will be useful, but
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
11 ## WITHOUT ANY WARRANTY; without even the implied warranty of
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
12 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
13 ## General Public License for more details.
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
14 ##
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
15 ## You should have received a copy of the GNU General Public License
7016
93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe]
jwe
parents: 5838
diff changeset
16 ## along with Octave; see the file COPYING. If not, see
93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe]
jwe
parents: 5838
diff changeset
17 ## <http://www.gnu.org/licenses/>.
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
18
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
19 ## -*- texinfo -*-
20852
516bb87ea72e 2015 Code Sprint: remove class of function from docstring for all m-files.
Rik <rik@octave.org>
parents: 20689
diff changeset
20 ## @deftypefn {} {@var{q} =} quadl (@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: 20689
diff changeset
21 ## @deftypefnx {} {@var{q} =} quadl (@var{f}, @var{a}, @var{b}, @var{tol})
516bb87ea72e 2015 Code Sprint: remove class of function from docstring for all m-files.
Rik <rik@octave.org>
parents: 20689
diff changeset
22 ## @deftypefnx {} {@var{q} =} quadl (@var{f}, @var{a}, @var{b}, @var{tol}, @var{trace})
516bb87ea72e 2015 Code Sprint: remove class of function from docstring for all m-files.
Rik <rik@octave.org>
parents: 20689
diff changeset
23 ## @deftypefnx {} {@var{q} =} quadl (@var{f}, @var{a}, @var{b}, @var{tol}, @var{trace}, @var{p1}, @var{p2}, @dots{})
516bb87ea72e 2015 Code Sprint: remove class of function from docstring for all m-files.
Rik <rik@octave.org>
parents: 20689
diff changeset
24 ## @deftypefnx {} {[@var{q}, @var{nfun}] =} quadl (@dots{})
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
25 ##
20158
7503499a252b doc: Update docstrings to have one sentence summary as first line.
Rik <rik@octave.org>
parents: 19697
diff changeset
26 ## Numerically evaluate the integral of @var{f} from @var{a} to @var{b} using
7503499a252b doc: Update docstrings to have one sentence summary as first line.
Rik <rik@octave.org>
parents: 19697
diff changeset
27 ## an adaptive Lobatto rule.
7503499a252b doc: Update docstrings to have one sentence summary as first line.
Rik <rik@octave.org>
parents: 19697
diff changeset
28 ##
7503499a252b doc: Update docstrings to have one sentence summary as first line.
Rik <rik@octave.org>
parents: 19697
diff changeset
29 ## @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: 19697
diff changeset
30 ## of the function to evaluate. The function @var{f} must be vectorized and
7503499a252b doc: Update docstrings to have one sentence summary as first line.
Rik <rik@octave.org>
parents: 19697
diff changeset
31 ## return a vector of output values when given a vector of input values.
12612
16cca721117b doc: Update all documentation for chapter on Numerical Integration
Rik <octave@nomad.inbox5.com>
parents: 12575
diff changeset
32 ##
16cca721117b doc: Update all documentation for chapter on Numerical Integration
Rik <octave@nomad.inbox5.com>
parents: 12575
diff changeset
33 ## @var{a} and @var{b} are the lower and upper limits of integration. Both
16cca721117b doc: Update all documentation for chapter on Numerical Integration
Rik <octave@nomad.inbox5.com>
parents: 12575
diff changeset
34 ## limits must be finite.
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
35 ##
20689
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
36 ## The optional argument @var{tol} defines the absolute tolerance with which
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
37 ## to perform the integration. The default value is @math{1e-6}.
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
38 ##
12612
16cca721117b doc: Update all documentation for chapter on Numerical Integration
Rik <octave@nomad.inbox5.com>
parents: 12575
diff changeset
39 ## The algorithm used by @code{quadl} involves recursively subdividing the
20158
7503499a252b doc: Update docstrings to have one sentence summary as first line.
Rik <rik@octave.org>
parents: 19697
diff changeset
40 ## integration interval. If @var{trace} is defined then for each subinterval
20689
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
41 ## display: (1) the total number of function evaluations, (2) the left end of
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
42 ## the subinterval, (3) the length of the subinterval, (4) the approximation of
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
43 ## the integral over the subinterval.
12612
16cca721117b doc: Update all documentation for chapter on Numerical Integration
Rik <octave@nomad.inbox5.com>
parents: 12575
diff changeset
44 ##
16cca721117b doc: Update all documentation for chapter on Numerical Integration
Rik <octave@nomad.inbox5.com>
parents: 12575
diff changeset
45 ## Additional arguments @var{p1}, etc., are passed directly to the function
16cca721117b doc: Update all documentation for chapter on Numerical Integration
Rik <octave@nomad.inbox5.com>
parents: 12575
diff changeset
46 ## @var{f}. To use default values for @var{tol} and @var{trace}, one may pass
16cca721117b doc: Update all documentation for chapter on Numerical Integration
Rik <octave@nomad.inbox5.com>
parents: 12575
diff changeset
47 ## empty matrices ([]).
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
48 ##
20689
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
49 ## The result of the integration is returned in @var{q}.
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
50 ##
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
51 ## The optional output @var{nfun} indicates the total number of function
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
52 ## evaluations performed.
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
53 ##
19040
0850b5212619 doc: Add @nospell macro around proper names in documentation.
Rik <rik@octave.org>
parents: 17744
diff changeset
54 ## Reference: @nospell{W. Gander and W. Gautschi}, @cite{Adaptive Quadrature -
10791
3140cb7a05a1 Add spellchecker scripts for Octave and run spellcheck of documentation
Rik <octave@nomad.inbox5.com>
parents: 10549
diff changeset
55 ## Revisited}, BIT Vol. 40, No. 1, March 2000, pp. 84--101.
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
56 ## @url{http://www.inf.ethz.ch/personal/gander/}
12575
d0b799dafede Grammarcheck files for 3.4.1 release.
Rik <octave@nomad.inbox5.com>
parents: 11587
diff changeset
57 ## @seealso{quad, quadv, quadgk, quadcc, trapz, dblquad, triplequad}
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
58 ## @end deftypefn
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
59
20689
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
60 ## Original Author: Walter Gautschi
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
61 ## Date: 08/03/98
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
62 ## Reference: Gander, Computermathematik, Birkhaeuser, 1992.
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
63
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
64 ## 2003-08-05 Shai Ayal
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
65 ## * permission from author to release as GPL
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
66
20689
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
67 function [q, nfun] = quadl (f, a, b, tol = [], trace = false, varargin)
13008
85dac13a911b quadl.m: Fix integration with large error tolerances (Bug #33792)
Rik <octave@nomad.inbox5.com>
parents: 12802
diff changeset
68
85dac13a911b quadl.m: Fix integration with large error tolerances (Bug #33792)
Rik <octave@nomad.inbox5.com>
parents: 12802
diff changeset
69 if (nargin < 3)
85dac13a911b quadl.m: Fix integration with large error tolerances (Bug #33792)
Rik <octave@nomad.inbox5.com>
parents: 12802
diff changeset
70 print_usage ();
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
71 endif
13008
85dac13a911b quadl.m: Fix integration with large error tolerances (Bug #33792)
Rik <octave@nomad.inbox5.com>
parents: 12802
diff changeset
72
7795
df9519e9990c Handle single precision eps values
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
73 if (isa (a, "single") || isa (b, "single"))
20689
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
74 eps = eps ("single");
7795
df9519e9990c Handle single precision eps values
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
75 else
20689
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
76 eps = eps ("double");
7795
df9519e9990c Handle single precision eps values
David Bateman <dbateman@free.fr>
parents: 7017
diff changeset
77 endif
5838
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
78 if (isempty (tol))
20689
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
79 tol = 1e-6;
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
80 elseif (! isscalar (tol) || tol < 0)
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
81 error ("quadl: TOL must be a scalar >=0");
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
82 elseif (tol < eps)
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
83 tol = eps;
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
84 endif
5838
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
85 if (isempty (trace))
13008
85dac13a911b quadl.m: Fix integration with large error tolerances (Bug #33792)
Rik <octave@nomad.inbox5.com>
parents: 12802
diff changeset
86 trace = false;
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
87 endif
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
88
20689
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
89 y = feval (f, [a, b], varargin{:});
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
90 nfun = 1;
5838
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
91
11587
c792872f8942 all script files: untabify and strip trailing whitespace
John W. Eaton <jwe@octave.org>
parents: 11562
diff changeset
92 fa = y(1);
20689
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
93 fb = y(2);
5838
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
94
20689
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
95 h = b - a;
5838
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
96
20689
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
97 [q, nfun, hmin] = adaptlobstp (f, a, b, fa, fb, Inf, nfun, abs (h),
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
98 tol, trace, varargin{:});
5838
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
99
20689
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
100 if (nfun > 10_000)
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
101 warning ("quadl: maximum iteration count reached -- possible singular integral");
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
102 elseif (any (! isfinite (q(:))))
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
103 warning ("quadl: infinite or NaN function evaluations were returned");
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
104 elseif (hmin < (b - a) * eps)
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
105 warning ("quadl: minimum step size reached -- possible singular integral");
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
106 endif
13008
85dac13a911b quadl.m: Fix integration with large error tolerances (Bug #33792)
Rik <octave@nomad.inbox5.com>
parents: 12802
diff changeset
107
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
108 endfunction
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
109
20689
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
110 function [q, nfun, hmin] = adaptlobstp (f, a, b, fa, fb, q0, nfun, hmin,
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
111 tol, trace, varargin)
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
112 persistent alpha = sqrt (2/3);
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
113 persistent beta = 1 / sqrt (5);
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
114
20689
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
115 if (nfun > 10_000)
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
116 q = q0;
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
117 return;
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
118 endif
13008
85dac13a911b quadl.m: Fix integration with large error tolerances (Bug #33792)
Rik <octave@nomad.inbox5.com>
parents: 12802
diff changeset
119
20689
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
120 h = (b - a) / 2;
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
121 m = (a + b) / 2;
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
122 mll = m - alpha*h;
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
123 ml = m - beta*h;
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
124 mr = m + beta*h;
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
125 mrr = m + alpha*h;
5838
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
126 x = [mll, ml, m, mr, mrr];
13008
85dac13a911b quadl.m: Fix integration with large error tolerances (Bug #33792)
Rik <octave@nomad.inbox5.com>
parents: 12802
diff changeset
127 y = feval (f, x, varargin{:});
20689
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
128 nfun += 1;
11587
c792872f8942 all script files: untabify and strip trailing whitespace
John W. Eaton <jwe@octave.org>
parents: 11562
diff changeset
129 fmll = y(1);
13008
85dac13a911b quadl.m: Fix integration with large error tolerances (Bug #33792)
Rik <octave@nomad.inbox5.com>
parents: 12802
diff changeset
130 fml = y(2);
85dac13a911b quadl.m: Fix integration with large error tolerances (Bug #33792)
Rik <octave@nomad.inbox5.com>
parents: 12802
diff changeset
131 fm = y(3);
85dac13a911b quadl.m: Fix integration with large error tolerances (Bug #33792)
Rik <octave@nomad.inbox5.com>
parents: 12802
diff changeset
132 fmr = y(4);
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
133 fmrr = y(5);
5838
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
134 i2 = (h/6)*(fa + fb + 5*(fml+fmr));
376e02b2ce70 [project @ 2006-06-01 20:23:53 by jwe]
jwe
parents: 5837
diff changeset
135 i1 = (h/1470)*(77*(fa+fb) + 432*(fmll+fmrr) + 625*(fml+fmr) + 672*fm);
20689
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
136
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
137 if (abs (b - a) < hmin)
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
138 hmin = abs (b - a);
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
139 endif
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
140
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
141 if (trace)
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
142 disp ([nfun, a, b-a, i1]);
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
143 endif
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
144
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
145 ## Force at least one adaptive step (nfun > 2 test).
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
146 if ((abs (i1-i2) < tol || mll <= a || b <= mrr) && nfun > 2)
12612
16cca721117b doc: Update all documentation for chapter on Numerical Integration
Rik <octave@nomad.inbox5.com>
parents: 12575
diff changeset
147 q = i1;
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
148 else
20689
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
149 q = zeros (6, 1);
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
150 [q(1), nfun, hmin] = adaptlobstp (f, a , mll, fa , fmll, q0/6, nfun, hmin,
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
151 tol, trace, varargin{:});
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
152 [q(2), nfun, hmin] = adaptlobstp (f, mll, ml , fmll, fml , q0/6, nfun, hmin,
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
153 tol, trace, varargin{:});
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
154 [q(3), nfun, hmin] = adaptlobstp (f, ml , m , fml , fm , q0/6, nfun, hmin,
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
155 tol, trace, varargin{:});
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
156 [q(4), nfun, hmin] = adaptlobstp (f, m , mr , fm , fmr , q0/6, nfun, hmin,
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
157 tol, trace, varargin{:});
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
158 [q(5), nfun, hmin] = adaptlobstp (f, mr , mrr, fmr , fmrr, q0/6, nfun, hmin,
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
159 tol, trace, varargin{:});
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
160 [q(6), nfun, hmin] = adaptlobstp (f, mrr, b , fmrr, fb , q0/6, nfun, hmin,
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
161 tol, trace, varargin{:});
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
162 q = sum (q);
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
163 endif
20689
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
164
5837
55404f3b0da1 [project @ 2006-06-01 19:05:31 by jwe]
jwe
parents:
diff changeset
165 endfunction
12802
412882f498b4 codesprint: Wrote 5 tests for quadl.m
David Wells <drwells@vt.edu>
parents: 12642
diff changeset
166
412882f498b4 codesprint: Wrote 5 tests for quadl.m
David Wells <drwells@vt.edu>
parents: 12642
diff changeset
167
412882f498b4 codesprint: Wrote 5 tests for quadl.m
David Wells <drwells@vt.edu>
parents: 12642
diff changeset
168 ## basic functionality
20689
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
169 %!assert (quadl (@(x) sin (x), 0, pi), 2, 5e-15)
12802
412882f498b4 codesprint: Wrote 5 tests for quadl.m
David Wells <drwells@vt.edu>
parents: 12642
diff changeset
170
412882f498b4 codesprint: Wrote 5 tests for quadl.m
David Wells <drwells@vt.edu>
parents: 12642
diff changeset
171 ## the values here are very high so it may be unavoidable that this fails
20689
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
172 %!assert (quadl (@(x) sin (3*x).*cosh (x).*sinh (x),10,15, 1e-3),
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
173 %! 2.588424538641647e+10, 1e-3)
12802
412882f498b4 codesprint: Wrote 5 tests for quadl.m
David Wells <drwells@vt.edu>
parents: 12642
diff changeset
174
412882f498b4 codesprint: Wrote 5 tests for quadl.m
David Wells <drwells@vt.edu>
parents: 12642
diff changeset
175 ## extra parameters
412882f498b4 codesprint: Wrote 5 tests for quadl.m
David Wells <drwells@vt.edu>
parents: 12642
diff changeset
176 %!assert (quadl (@(x,a,b) sin (a + b*x), 0, 1, [], [], 2, 3),
20689
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
177 %! cos(2)/3 - cos(5)/3, 1e-15)
12802
412882f498b4 codesprint: Wrote 5 tests for quadl.m
David Wells <drwells@vt.edu>
parents: 12642
diff changeset
178
13141
e81ddf9cacd5 maint: untabify and remove trailing whitespace from source files
John W. Eaton <jwe@octave.org>
parents: 13008
diff changeset
179 ## test different tolerances.
20689
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
180 %!test
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
181 %! [q, nfun1] = quadl (@(x) sin (2 + 3*x).^2, 0, 10, 0.5, []);
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
182 %! assert (q, (60 + sin(4) - sin(64))/12, 0.5);
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
183 %! [q, nfun2] = quadl (@(x) sin (2 + 3*x).^2, 0, 10, 0.1, []);
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
184 %! assert (q, (60 + sin(4) - sin(64))/12, 0.1);
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
185 %! assert (nfun2 > nfun1);
12802
412882f498b4 codesprint: Wrote 5 tests for quadl.m
David Wells <drwells@vt.edu>
parents: 12642
diff changeset
186
20689
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
187 ## Test input validation
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
188 %!error quadl ()
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
189 %!error quadl (@sin)
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
190 %!error quadl (@sin,1)
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
191 %!error <TOL must be a scalar> quadl (@sin, 0, 1, ones (2,2))
655816377845 quadl.m: Overhal function and switch to absolute tolerance.
Rik <rik@octave.org>
parents: 20231
diff changeset
192