| blob | 89e4333a938af47b98f072d04d45b6ced22e126d |
1 @menu
2 * Introduction to Integration::
3 * Functions and Variables for Integration::
4 * Introduction to QUADPACK::
5 * Functions and Variables for QUADPACK::
6 @end menu
8 @c -----------------------------------------------------------------------------
9 @node Introduction to Integration, Functions and Variables for Integration, Integration, Integration
10 @section Introduction to Integration
11 @c -----------------------------------------------------------------------------
13 Maxima has several routines for handling integration.
14 The @mref{integrate} function makes use of most of them. There is also the
15 @mref{antid} package, which handles an unspecified function (and its
16 derivatives, of course). For numerical uses,
17 there is a set of adaptive integrators from QUADPACK, named @mrefcomma{quad_qag}
18 @mrefcomma{quad_qags} etc., which are described under the heading @code{QUADPACK}.
19 Hypergeometric functions are being worked on,
20 see @mref{specint} for details.
21 Generally speaking, Maxima only handles integrals which are
22 integrable in terms of the "elementary functions" (rational functions,
23 trigonometrics, logs, exponentials, radicals, etc.) and a few
24 extensions (error function, dilogarithm). It does not handle
25 integrals in terms of unknown functions such as @code{g(x)} and @code{h(x)}.
27 @c end concepts Integration
29 @c -----------------------------------------------------------------------------
30 @node Functions and Variables for Integration, Introduction to QUADPACK, Introduction to Integration, Integration
31 @section Functions and Variables for Integration
32 @c -----------------------------------------------------------------------------
34 @c NEEDS WORK
36 @c -----------------------------------------------------------------------------
37 @anchor{changevar}
38 @deffn {Function} changevar (@var{expr}, @var{f(x,y)}, @var{y}, @var{x})
40 Makes the change of variable given by @code{@var{f(x,y)} = 0} in all integrals
41 occurring in @var{expr} with integration with respect to @var{x}.
42 The new variable is @var{y}.
44 The change of variable can also be written @code{@var{f(x)} = @var{g(y)}}.
46 @c HMM, THIS EXAMPLE YIELDS A CORRECT BUT SLIGHTLY STRANGE RESULT...
47 @c ===beg===
48 @c assume(a > 0)$
49 @c 'integrate (%e**sqrt(a*y), y, 0, 4);
50 @c changevar (%, y-z^2/a, z, y);
51 @c ===end===
52 @example
53 (%i1) assume(a > 0)$
54 @group
55 (%i2) 'integrate (%e**sqrt(a*y), y, 0, 4);
56 4
57 /
58 [ sqrt(a) sqrt(y)
59 (%o2) I %e dy
60 ]
61 /
62 0
63 @end group
64 @group
65 (%i3) changevar (%, y-z^2/a, z, y);
66 0
67 /
68 [ abs(z)
69 2 I z %e dz
70 ]
71 /
72 - 2 sqrt(a)
73 (%o3) - ----------------------------
74 a
75 @end group
76 @end example
78 An expression containing a noun form, such as the instances of @code{'integrate}
79 above, may be evaluated by @code{ev} with the @code{nouns} flag.
80 For example, the expression returned by @code{changevar} above may be evaluated
81 by @code{ev (%o3, nouns)}.
83 @code{changevar} may also be used to make changes in the indices of a sum or
84 product. However, it must be realized that when a change is made in a
85 sum or product, this change must be a shift, i.e., @code{i = j+ ...}, not a
86 higher degree function. E.g.,
88 @c ===beg===
89 @c sum (a[i]*x^(i-2), i, 0, inf);
90 @c changevar (%, i-2-n, n, i);
91 @c ===end===
92 @example
93 @group
94 (%i4) sum (a[i]*x^(i-2), i, 0, inf);
95 inf
96 ====
97 \ i - 2
98 (%o4) > a x
99 / i
100 ====
101 i = 0
102 @end group
103 @group
104 (%i5) changevar (%, i-2-n, n, i);
105 inf
106 ====
107 \ n
108 (%o5) > a x
109 / n + 2
110 ====
111 n = - 2
112 @end group
113 @end example
115 @opencatbox{Categories:}
116 @category{Integral calculus}
117 @closecatbox
118 @end deffn
120 @c THIS ITEM IS A MESS, BUT DON'T BOTHER TO CLEAN IT UP:
121 @c THE GAUSS-KRONROD FUNCTIONS (QUADPACK) MAKE THIS OBSOLETE
123 @c -----------------------------------------------------------------------------
124 @anchor{dblint}
125 @deffn {Function} dblint (@var{f}, @var{r}, @var{s}, @var{a}, @var{b})
127 A double-integral routine which was written in
128 top-level Maxima and then translated and compiled to machine code.
129 Use @code{load ("dblint")} to access this package. It uses the Simpson's rule
130 method in both the x and y directions to calculate
132 @tex
133 $$\int_a^b \int_{r\left(x\right)}^{s\left(x\right)} f\left(x,y\right) \, dy \, dx.$$
134 @end tex
135 @ifnottex
136 @example
137 @group
138 /b /s(x)
139 | |
140 | | f(x,y) dy dx
141 | |
142 /a /r(x)
143 @end group
144 @end example
145 @end ifnottex
147 The function @var{f} must be a translated or compiled function of two variables,
148 and @var{r} and @var{s} must each be a translated or compiled function of one
149 variable, while @var{a} and @var{b} must be floating point numbers. The routine
150 has two global variables which determine the number of divisions of the x and y
151 intervals: @code{dblint_x} and @code{dblint_y}, both of which are initially 10,
152 and can be changed independently to other integer values (there are
153 @code{2*dblint_x+1} points computed in the x direction, and @code{2*dblint_y+1}
154 in the y direction). The routine subdivides the X axis and then for each value
155 of X it first computes @code{@var{r}(x)} and @code{@var{s}(x)}; then the Y axis
156 between @code{@var{r}(x)} and @code{@var{s}(x)} is subdivided and the integral
157 along the Y axis is performed using Simpson's rule; then the integral along the
158 X axis is done using Simpson's rule with the function values being the
159 Y-integrals. This procedure may be numerically unstable for a great variety of
160 reasons, but is reasonably fast: avoid using it on highly oscillatory functions
161 and functions with singularities (poles or branch points in the region). The Y
162 integrals depend on how far apart @code{@var{r}(x)} and @code{@var{s}(x)} are,
163 so if the distance @code{@var{s}(x) - @var{r}(x)} varies rapidly with X, there
164 may be substantial errors arising from truncation with different step-sizes in
165 the various Y integrals. One can increase @code{dblint_x} and @code{dblint_y}
166 in an effort to improve the coverage of the region, at the expense of
167 computation time. The function values are not saved, so if the function is very
168 time-consuming, you will have to wait for re-computation if you change anything
169 (sorry). It is required that the functions @var{f}, @var{r}, and @var{s} be
170 either translated or compiled prior to calling @code{dblint}. This will result
171 in orders of magnitude speed improvement over interpreted code in many cases!
173 @code{demo ("dblint")} executes a demonstration of @code{dblint} applied to an
174 example problem.
175 @c demo (dblint_1) FAILS WITH Could not find `fltdfnk.mc' -- DON'T BOTHER TO MENTION IT. !!!
176 @c @code{demo (dblint_1)} executes another demonstration.
178 @opencatbox{Categories:}
179 @category{Integral calculus}
180 @closecatbox
181 @end deffn
183 @c -----------------------------------------------------------------------------
184 @anchor{defint}
185 @deffn {Function} defint (@var{expr}, @var{x}, @var{a}, @var{b})
187 Attempts to compute a definite integral. @code{defint} is called by
188 @code{integrate} when limits of integration are specified, i.e., when
189 @code{integrate} is called as
190 @code{integrate (@var{expr}, @var{x}, @var{a}, @var{b})}.
191 Thus from the user's point of view, it is sufficient to call @code{integrate}.
192 @c SHOULD WE BOTHER TO DOCUMENT defint ??? NO FUNCTIONALITY HERE THAT IS NOT ALREADY PRESENT IN integrate !!!
194 @code{defint} returns a symbolic expression, either the computed integral or the
195 noun form of the integral. See @mref{quad_qag} and related functions for
196 numerical approximation of definite integrals.
198 @opencatbox{Categories:}
199 @category{Integral calculus}
200 @closecatbox
201 @end deffn
203 @c -----------------------------------------------------------------------------
204 @anchor{erfflag}
205 @defvr {Option variable} erfflag
206 Default value: @code{true}
208 When @code{erfflag} is @code{false}, prevents @code{risch} from introducing the
209 @code{erf} function in the answer if there were none in the integrand to
210 begin with.
212 @opencatbox{Categories:}
213 @category{Integral calculus}
214 @closecatbox
215 @end defvr
217 @c NEEDS WORK
219 @c -----------------------------------------------------------------------------
220 @anchor{ilt}
221 @deffn {Function} ilt (@var{expr}, @var{s}, @var{t})
223 Computes the inverse Laplace transform of @var{expr} with
224 respect to @var{s} and parameter @var{t}. @var{expr} must be a ratio of
225 polynomials whose denominator has only linear and quadratic factors;
226 there is an extension of @code{ilt}, called @mref{pwilt} (Piece-Wise
227 Inverse Laplace Transform) that handles several other cases where
228 @code{ilt} fails.
230 By using the functions @code{laplace} and @code{ilt} together with the
231 @code{solve} or @code{linsolve} functions the user can solve a single
232 differential or convolution integral equation or a set of them.
234 @c ===beg===
235 @c 'integrate (sinh(a*x)*f(t-x), x, 0, t) + b*f(t) = t**2;
236 @c laplace (%, t, s);
237 @c linsolve ([%], ['laplace(f(t), t, s)]);
238 @c ilt (rhs (first (%)), s, t);
239 @c input:pos;
240 @c ===end===
241 @example
242 @group
243 (%i1) 'integrate (sinh(a*x)*f(t-x), x, 0, t) + b*f(t) = t**2;
244 t
245 /
246 [ 2
247 (%o1) I f(t - x) sinh(a x) dx + b f(t) = t
248 ]
249 /
250 0
251 @end group
252 @group
253 (%i2) laplace (%, t, s);
254 a laplace(f(t), t, s) 2
255 (%o2) b laplace(f(t), t, s) + --------------------- = --
256 2 2 3
257 s - a s
258 @end group
259 @group
260 (%i3) linsolve ([%], ['laplace(f(t), t, s)]);
261 2 2
262 2 s - 2 a
263 (%o3) [laplace(f(t), t, s) = --------------------]
264 5 2 3
265 b s + (a - a b) s
266 @end group
267 @group
268 (%i4) ilt (rhs (first (%)), s, t);
269 Is a b (a b - 1) positive, negative, or zero?
271 pos;
272 sqrt(a b (a b - 1)) t
273 2 cosh(---------------------) 2
274 b a t
275 (%o4) - ----------------------------- + -------
276 3 2 2 a b - 1
277 a b - 2 a b + a
279 2
280 + ------------------
281 3 2 2
282 a b - 2 a b + a
283 @end group
284 @end example
286 @opencatbox{Categories:}
287 @category{Laplace transform}
288 @closecatbox
289 @end deffn
291 @c -----------------------------------------------------------------------------
292 @anchor{intanalysis}
293 @defvr {Option variable} intanalysis
294 Default value: @code{true}
296 When @code{true}, definite integration tries to find poles in the integrand in
297 the interval of integration. If there are, then the integral is evaluated
298 appropriately as a principal value integral. If intanalysis is @code{false},
299 this check is not performed and integration is done assuming there are no poles.
301 See also @mrefdot{ldefint}
303 Examples:
305 Maxima can solve the following integrals, when @mref{intanalysis} is set to
306 @code{false}:
308 @c ===beg===
309 @c integrate(1/(sqrt(x+1)+1),x,0,1);
310 @c integrate(1/(sqrt(x)+1),x,0,1),intanalysis:false;
311 @c integrate(cos(a)/sqrt((tan(a))^2+1),a,-%pi/2,%pi/2);
312 @c intanalysis:false$
313 @c integrate(cos(a)/sqrt((tan(a))^2 +1),a,-%pi/2,%pi/2);
314 @c ===end===
315 @example
316 (%i1) integrate(1/(sqrt(x)+1),x,0,1);
317 1
318 /
319 [ 1
320 (%o1) I ----------- dx
321 ] sqrt(x) + 1
322 /
323 0
325 (%i2) integrate(1/(sqrt(x)+1),x,0,1),intanalysis:false;
326 (%o2) 2 - 2 log(2)
328 (%i3) integrate(cos(a)/sqrt((tan(a))^2 +1),a,-%pi/2,%pi/2);
329 The number 1 isn't in the domain of atanh
330 -- an error. To debug this try: debugmode(true);
332 (%i4) intanalysis:false$
333 (%i5) integrate(cos(a)/sqrt((tan(a))^2+1),a,-%pi/2,%pi/2);
334 %pi
335 (%o5) ---
336 2
337 @end example
339 @opencatbox{Categories:}
340 @category{Integral calculus}
341 @closecatbox
342 @end defvr
344 @c -----------------------------------------------------------------------------
345 @anchor{integrate}
346 @deffn {Function} integrate @
347 @fname{integrate} (@var{expr}, @var{x}) @
348 @fname{integrate} (@var{expr}, @var{x}, @var{a}, @var{b})
350 Attempts to symbolically compute the integral of @var{expr} with respect to
351 @var{x}. @code{integrate (@var{expr}, @var{x})} is an indefinite integral,
352 while @code{integrate (@var{expr}, @var{x}, @var{a}, @var{b})} is a definite
353 integral, with limits of integration @var{a} and @var{b}. The limits should
354 not contain @var{x}, although @code{integrate} does not enforce this
355 restriction. @var{a} need not be less than @var{b}.
356 If @var{b} is equal to @var{a}, @code{integrate} returns zero.
358 See @mref{quad_qag} and related functions for numerical approximation of
359 definite integrals. See @mref{residue} for computation of residues
360 (complex integration). See @mref{antid} for an alternative means of computing
361 indefinite integrals.
363 The integral (an expression free of @code{integrate}) is returned if
364 @code{integrate} succeeds. Otherwise the return value is
365 the noun form of the integral (the quoted operator @code{'integrate})
366 or an expression containing one or more noun forms.
367 The noun form of @code{integrate} is displayed with an integral sign.
369 In some circumstances it is useful to construct a noun form by hand, by quoting
370 @code{integrate} with a single quote, e.g.,
371 @code{'integrate (@var{expr}, @var{x})}. For example, the integral may depend
372 on some parameters which are not yet computed.
373 The noun may be applied to its arguments by @code{ev (@var{i}, nouns)}
374 where @var{i} is the noun form of interest.
376 @c BEGIN EXPOSITION ON HEURISTICS
377 @code{integrate} handles definite integrals separately from indefinite, and
378 employs a range of heuristics to handle each case. Special cases of definite
379 integrals include limits of integration equal to zero or infinity (@mref{inf} or
380 @mref{minf}), trigonometric functions with limits of integration equal to zero
381 and @code{%pi} or @code{2 %pi}, rational functions, integrals related to the
382 definitions of the @mref{beta} and @mref{psi} functions, and some logarithmic
383 and trigonometric integrals. Processing rational functions may include
384 computation of residues. If an applicable special case is not found, an attempt
385 will be made to compute the indefinite integral and evaluate it at the limits of
386 integration. This may include taking a limit as a limit of integration goes to
387 infinity or negative infinity; see also @mrefdot{ldefint}
389 Special cases of indefinite integrals include trigonometric functions,
390 exponential and logarithmic functions,
391 and rational functions.
392 @code{integrate} may also make use of a short table of elementary integrals.
394 @code{integrate} may carry out a change of variable
395 if the integrand has the form @code{f(g(x)) * diff(g(x), x)}.
396 @code{integrate} attempts to find a subexpression @code{g(x)} such that
397 the derivative of @code{g(x)} divides the integrand.
398 This search may make use of derivatives defined by the @code{gradef} function.
399 See also @mref{changevar} and @mrefdot{antid}
401 If none of the preceding heuristics find the indefinite integral, the Risch
402 algorithm is executed. The flag @mref{risch} may be set as an @mrefcomma{evflag}
403 in a call to @code{ev} or on the command line, e.g.,
404 @code{ev (integrate (@var{expr}, @var{x}), risch)} or
405 @code{integrate (@var{expr}, @var{x}), risch}. If @code{risch} is present,
406 @code{integrate} calls the @mref{risch} function without attempting heuristics
407 first. See also @mrefdot{risch}
408 @c END EXPOSITION ON HEURISTICS
410 @code{integrate} works only with functional relations represented explicitly
411 with the @code{f(x)} notation. @code{integrate} does not respect implicit
412 dependencies established by the @mref{depends} function.
414 @code{integrate} may need to know some property of a parameter in the integrand.
415 @code{integrate} will first consult the @mref{assume} database,
416 and, if the variable of interest is not there,
417 @code{integrate} will ask the user.
418 Depending on the question,
419 suitable responses are @code{yes;} or @code{no;},
420 or @code{pos;}, @code{zero;}, or @code{neg;}.
422 @code{integrate} is not, by default, declared to be linear. See @code{declare}
423 and @code{linear}.
425 @code{integrate} attempts integration by parts only in a few special cases.
427 Examples:
429 @itemize @bullet
430 @item
431 Elementary indefinite and definite integrals.
433 @c ===beg===
434 @c integrate (sin(x)^3, x);
435 @c integrate (x/ sqrt (b^2 - x^2), x);
436 @c integrate (cos(x)^2 * exp(x), x, 0, %pi);
437 @c integrate (x^2 * exp(-x^2), x, minf, inf);
438 @c ===end===
439 @example
440 @group
441 (%i1) integrate (sin(x)^3, x);
442 3
443 cos (x)
444 (%o1) ------- - cos(x)
445 3
446 @end group
447 @group
448 (%i2) integrate (x/ sqrt (b^2 - x^2), x);
449 2 2
450 (%o2) - sqrt(b - x )
451 @end group
452 @group
453 (%i3) integrate (cos(x)^2 * exp(x), x, 0, %pi);
454 %pi
455 3 %e 3
456 (%o3) ------- - -
457 5 5
458 @end group
459 @group
460 (%i4) integrate (x^2 * exp(-x^2), x, minf, inf);
461 sqrt(%pi)
462 (%o4) ---------
463 2
464 @end group
465 @end example
467 @item
468 Use of @code{assume} and interactive query.
470 @c ===beg===
471 @c assume (a > 1)$
472 @c integrate (x**a/(x+1)**(5/2), x, 0, inf);
473 @c input:no;
474 @c input:neg;
475 @c ===end===
476 @example
477 (%i1) assume (a > 1)$
478 @group
479 (%i2) integrate (x**a/(x+1)**(5/2), x, 0, inf);
480 2 a + 2
481 Is ------- an integer?
482 5
484 no;
485 Is 2 a - 3 positive, negative, or zero?
487 neg;
488 3
489 (%o2) beta(a + 1, - - a)
490 2
491 @end group
492 @end example
494 @item
495 Change of variable. There are two changes of variable in this example:
496 one using a derivative established by @mrefcomma{gradef} and one using the
497 derivation @code{diff(r(x))} of an unspecified function @code{r(x)}.
499 @c ===beg===
500 @c gradef (q(x), sin(x**2));
501 @c diff (log (q (r (x))), x);
502 @c integrate (%, x);
503 @c ===end===
504 @example
505 @group
506 (%i3) gradef (q(x), sin(x**2));
507 (%o3) q(x)
508 @end group
509 @group
510 (%i4) diff (log (q (r (x))), x);
511 d 2
512 (-- (r(x))) sin(r (x))
513 dx
514 (%o4) ----------------------
515 q(r(x))
516 @end group
517 @group
518 (%i5) integrate (%, x);
519 (%o5) log(q(r(x)))
520 @end group
521 @end example
523 @item
524 Return value contains the @code{'integrate} noun form. In this example, Maxima
525 can extract one factor of the denominator of a rational function, but cannot
526 factor the remainder or otherwise find its integral. @mref{grind} shows the
527 noun form @code{'integrate} in the result. See also
528 @mref{integrate_use_rootsof} for more on integrals of rational functions.
530 @c ===beg===
531 @c expand ((x-4) * (x^3+2*x+1));
532 @c integrate (1/%, x);
533 @c grind (%);
534 @c ===end===
535 @example
536 @group
537 (%i1) expand ((x-4) * (x^3+2*x+1));
538 4 3 2
539 (%o1) x - 4 x + 2 x - 7 x - 4
540 @end group
541 @group
542 (%i2) integrate (1/%, x);
543 / 2
544 [ x + 4 x + 18
545 I ------------- dx
546 ] 3
547 log(x - 4) / x + 2 x + 1
548 (%o2) ---------- - ------------------
549 73 73
550 @end group
551 @group
552 (%i3) grind (%);
553 log(x-4)/73-('integrate((x^2+4*x+18)/(x^3+2*x+1),x))/73$
554 @end group
555 @end example
557 @item
558 Defining a function in terms of an integral. The body of a function is not
559 evaluated when the function is defined. Thus the body of @code{f_1} in this
560 example contains the noun form of @code{integrate}. The quote-quote operator
561 @code{'@w{}'} causes the integral to be evaluated, and the result becomes the
562 body of @code{f_2}.
564 @c ===beg===
565 @c f_1 (a) := integrate (x^3, x, 1, a);
566 @c ev (f_1 (7), nouns);
567 @c /* Note parentheses around integrate(...) here */ f_2 (a) := ''(integrate (x^3, x, 1, a));
568 @c f_2 (7);
569 @c ===end===
570 @example
571 @group
572 (%i1) f_1 (a) := integrate (x^3, x, 1, a);
573 3
574 (%o1) f_1(a) := integrate(x , x, 1, a)
575 @end group
576 @group
577 (%i2) ev (f_1 (7), nouns);
578 (%o2) 600
579 @end group
580 @group
581 (%i3) /* Note parentheses around integrate(...) here */
582 f_2 (a) := ''(integrate (x^3, x, 1, a));
583 4
584 a 1
585 (%o3) f_2(a) := -- - -
586 4 4
587 @end group
588 @group
589 (%i4) f_2 (7);
590 (%o4) 600
591 @end group
592 @end example
593 @end itemize
595 @opencatbox{Categories:}
596 @category{Integral calculus}
597 @closecatbox
598 @end deffn
600 @c -----------------------------------------------------------------------------
601 @anchor{integration_constant}
602 @defvr {System variable} integration_constant
603 Default value: @code{%c}
605 When a constant of integration is introduced by indefinite integration of an
606 equation, the name of the constant is constructed by concatenating
607 @code{integration_constant} and @code{integration_constant_counter}.
609 @code{integration_constant} may be assigned any symbol.
611 Examples:
613 @c ===beg===
614 @c integrate (x^2 = 1, x);
615 @c integration_constant : 'k;
616 @c integrate (x^2 = 1, x);
617 @c ===end===
618 @example
619 @group
620 (%i1) integrate (x^2 = 1, x);
621 3
622 x
623 (%o1) -- = x + %c1
624 3
625 @end group
626 @group
627 (%i2) integration_constant : 'k;
628 (%o2) k
629 @end group
630 @group
631 (%i3) integrate (x^2 = 1, x);
632 3
633 x
634 (%o3) -- = x + k2
635 3
636 @end group
637 @end example
639 @opencatbox{Categories:}
640 @category{Integral calculus}
641 @closecatbox
642 @end defvr
644 @c -----------------------------------------------------------------------------
645 @anchor{integration_constant_counter}
646 @defvr {System variable} integration_constant_counter
647 Default value: 0
649 When a constant of integration is introduced by indefinite integration of an
650 equation, the name of the constant is constructed by concatenating
651 @code{integration_constant} and @code{integration_constant_counter}.
653 @code{integration_constant_counter} is incremented before constructing the next
654 integration constant.
656 Examples:
658 @c ===beg===
659 @c integrate (x^2 = 1, x);
660 @c integrate (x^2 = 1, x);
661 @c integrate (x^2 = 1, x);
662 @c reset (integration_constant_counter);
663 @c integrate (x^2 = 1, x);
664 @c ===end===
665 @example
666 @group
667 (%i1) integrate (x^2 = 1, x);
668 3
669 x
670 (%o1) -- = x + %c1
671 3
672 @end group
673 @group
674 (%i2) integrate (x^2 = 1, x);
675 3
676 x
677 (%o2) -- = x + %c2
678 3
679 @end group
680 @group
681 (%i3) integrate (x^2 = 1, x);
682 3
683 x
684 (%o3) -- = x + %c3
685 3
686 @end group
687 @group
688 (%i4) reset (integration_constant_counter);
689 (%o4) [integration_constant_counter]
690 @end group
691 @group
692 (%i5) integrate (x^2 = 1, x);
693 3
694 x
695 (%o5) -- = x + %c1
696 3
697 @end group
698 @end example
700 @opencatbox{Categories:}
701 @category{Integral calculus}
702 @closecatbox
703 @end defvr
705 @c -----------------------------------------------------------------------------
706 @anchor{integrate_use_rootsof}
707 @defvr {Option variable} integrate_use_rootsof
708 Default value: @code{false}
710 When @code{integrate_use_rootsof} is @code{true} and the denominator of
711 a rational function cannot be factored, @mref{integrate} returns the integral
712 in a form which is a sum over the roots (not yet known) of the denominator.
714 For example, with @code{integrate_use_rootsof} set to @code{false},
715 @code{integrate} returns an unsolved integral of a rational function in noun
716 form:
718 @c ===beg===
719 @c integrate_use_rootsof: false$
720 @c integrate (1/(1+x+x^5), x);
721 @c ===end===
722 @example
723 (%i1) integrate_use_rootsof: false$
724 @group
725 (%i2) integrate (1/(1+x+x^5), x);
726 / 2
727 [ x - 4 x + 5
728 I ------------ dx 2 x + 1
729 ] 3 2 2 5 atan(-------)
730 / x - x + 1 log(x + x + 1) sqrt(3)
731 (%o2) ----------------- - --------------- + ---------------
732 7 14 7 sqrt(3)
733 @end group
734 @end example
736 Now we set the flag to be true and the unsolved part of the integral will be
737 expressed as a summation over the roots of the denominator of the rational
738 function:
740 @c ===beg===
741 @c integrate_use_rootsof: true$
742 @c integrate (1/(1+x+x^5), x);
743 @c ===end===
744 @example
745 (%i3) integrate_use_rootsof: true$
746 @group
747 (%i4) integrate (1/(1+x+x^5), x);
748 ==== 2
749 \ (%r4 - 4 %r4 + 5) log(x - %r4)
750 > -------------------------------
751 / 2
752 ==== 3 %r4 - 2 %r4
753 3 2
754 %r4 in rootsof(%r4 - %r4 + 1, %r4)
755 (%o4) ----------------------------------------------------------
756 7
758 2 x + 1
759 2 5 atan(-------)
760 log(x + x + 1) sqrt(3)
761 - --------------- + ---------------
762 14 7 sqrt(3)
763 @end group
764 @end example
766 Alternatively the user may compute the roots of the denominator separately,
767 and then express the integrand in terms of these roots, e.g.,
768 @code{1/((x - a)*(x - b)*(x - c))} or @code{1/((x^2 - (a+b)*x + a*b)*(x - c))}
769 if the denominator is a cubic polynomial.
770 Sometimes this will help Maxima obtain a more useful result.
772 @opencatbox{Categories:}
773 @category{Integral calculus}
774 @closecatbox
775 @end defvr
777 @c NEEDS EXAMPLES
779 @c -----------------------------------------------------------------------------
780 @anchor{ldefint}
781 @deffn {Function} ldefint (@var{expr}, @var{x}, @var{a}, @var{b})
783 Attempts to compute the definite integral of @var{expr} by using @mref{limit}
784 to evaluate the indefinite integral of @var{expr} with respect to @var{x}
785 at the upper limit @var{b} and at the lower limit @var{a}.
786 If it fails to compute the definite integral,
787 @code{ldefint} returns an expression containing limits as noun forms.
789 @code{ldefint} is not called from @mrefcomma{integrate} so executing
790 @code{ldefint (@var{expr}, @var{x}, @var{a}, @var{b})} may yield a different
791 result than @code{integrate (@var{expr}, @var{x}, @var{a}, @var{b})}.
792 @code{ldefint} always uses the same method to evaluate the definite integral,
793 while @code{integrate} may employ various heuristics and may recognize some
794 special cases.
796 @opencatbox{Categories:}
797 @category{Integral calculus}
798 @closecatbox
799 @end deffn
801 @c UMM, IS THERE SOME TEXT MISSING HERE ???
802 @c WHAT IS THIS ABOUT EXACTLY ??
804 @c -----------------------------------------------------------------------------
805 @anchor{pwilt}
806 @deffn {Function} pwilt (@var{expr}, @var{s}, @var{t})
808 Computes the inverse Laplace transform of @var{expr} with
809 respect to @var{s} and parameter @var{t}. Unlike @mrefcomma{ilt}
810 @code{pwilt} is able to return piece-wise and periodic functions
811 and can also handle some cases with polynomials of degree greater than 3
812 in the denominator.
814 Two examples where @code{ilt} fails:
815 @c ===beg===
816 @c pwilt (exp(-s)*s/(s^3-2*s-s+2), s, t);
817 @c pwilt ((s^2+2)/(s^2-1), s, t);
818 @c ===end===
819 @example
820 (%i1) pwilt (exp(-s)*s/(s^3-2*s-s+2), s, t);
821 t - 1 - 2 (t - 1)
822 (t - 1) %e 2 %e
823 (%o1) hstep(t - 1) (--------------- - ---------------)
824 3 9
826 (%i2) pwilt ((s^2+2)/(s^2-1), s, t);
827 t - t
828 3 %e 3 %e
829 (%o2) delta(t) + ----- - -------
830 2 2
831 @end example
833 @opencatbox{Categories:}
834 @category{Laplace transform}
835 @closecatbox
836 @end deffn
838 @c -----------------------------------------------------------------------------
839 @anchor{potential}
840 @deffn {Function} potential (@var{givengradient})
842 The calculation makes use of the global variable @code{potentialzeroloc[0]}
843 which must be @code{nonlist} or of the form
845 @example
846 [indeterminatej=expressionj, indeterminatek=expressionk, ...]
847 @end example
849 the former being equivalent to the nonlist expression for all right-hand
850 sides in the latter. The indicated right-hand sides are used as the
851 lower limit of integration. The success of the integrations may
852 depend upon their values and order. @code{potentialzeroloc} is initially set
853 to 0.
854 @end deffn
856 @c -----------------------------------------------------------------------------
857 @anchor{residue}
858 @deffn {Function} residue (@var{expr}, @var{z}, @var{z_0})
860 Computes the residue in the complex plane of the expression @var{expr} when the
861 variable @var{z} assumes the value @var{z_0}. The residue is the coefficient of
862 @code{(@var{z} - @var{z_0})^(-1)} in the Laurent series for @var{expr}.
864 @c ===beg===
865 @c residue (s/(s**2+a**2), s, a*%i);
866 @c residue (sin(a*x)/x**4, x, 0);
867 @c ===end===
868 @example
869 @group
870 (%i1) residue (s/(s**2+a**2), s, a*%i);
871 1
872 (%o1) -
873 2
874 @end group
875 @group
876 (%i2) residue (sin(a*x)/x**4, x, 0);
877 3
878 a
879 (%o2) - --
880 6
881 @end group
882 @end example
884 @opencatbox{Categories:}
885 @category{Integral calculus}
886 @category{Complex variables}
887 @closecatbox
888 @end deffn
890 @c -----------------------------------------------------------------------------
891 @anchor{risch}
892 @deffn {Function} risch (@var{expr}, @var{x})
894 Integrates @var{expr} with respect to @var{x} using the
895 transcendental case of the Risch algorithm. (The algebraic case of
896 the Risch algorithm has not been implemented.) This currently
897 handles the cases of nested exponentials and logarithms which the main
898 part of @code{integrate} can't do. @mref{integrate} will automatically apply
899 @code{risch} if given these cases.
901 @code{erfflag}, if @code{false}, prevents @code{risch} from introducing the
902 @code{erf} function in the answer if there were none in the integrand to begin
903 with.
905 @c ===beg===
906 @c risch (x^2*erf(x), x);
907 @c diff(%, x), ratsimp;
908 @c ===end===
909 @example
910 @group
911 (%i1) risch (x^2*erf(x), x);
912 2
913 3 2 - x
914 %pi x erf(x) + (sqrt(%pi) x + sqrt(%pi)) %e
915 (%o1) -------------------------------------------------
916 3 %pi
917 @end group
918 @group
919 (%i2) diff(%, x), ratsimp;
920 2
921 (%o2) x erf(x)
922 @end group
923 @end example
925 @opencatbox{Categories:}
926 @category{Integral calculus}
927 @closecatbox
928 @end deffn
930 @c NEEDS EXPANSION, CLARIFICATION, AND EXAMPLES
932 @c -----------------------------------------------------------------------------
933 @anchor{tldefint}
934 @deffn {Function} tldefint (@var{expr}, @var{x}, @var{a}, @var{b})
936 Equivalent to @code{ldefint} with @code{tlimswitch} set to @code{true}.
938 @opencatbox{Categories:}
939 @category{Integral calculus}
940 @closecatbox
941 @end deffn
943 @footnotestyle end
945 @c -----------------------------------------------------------------------------
946 @node Introduction to QUADPACK, Functions and Variables for QUADPACK, Functions and Variables for Integration, Integration
947 @section Introduction to QUADPACK
948 @c -----------------------------------------------------------------------------
950 @c FOLLOWING TEXT ADAPTED WITH HEAVY MODIFICATION FROM https://www.netlib.org/slatec/src/qpdoc.f
952 QUADPACK is a collection of functions for the numerical
953 computation of one-dimensional definite integrals.
954 It originated from a joint project of
955 R. Piessens @footnote{Applied Mathematics and Programming Division, K.U. Leuven},
956 E. de Doncker @footnote{Applied Mathematics and Programming Division, K.U. Leuven},
957 C. Ueberhuber @footnote{Institut f@"ur Mathematik, T.U. Wien},
958 and D. Kahaner @footnote{National Bureau of Standards, Washington, D.C., U.S.A}.
960 The QUADPACK library included in Maxima is an automatic translation (via the
961 program @code{f2cl}) of the Fortran source code of QUADPACK as it appears in
962 the SLATEC Common Mathematical Library, Version 4.1 @footnote{@url{https://www.netlib.org/slatec}}.
963 The SLATEC library is dated July 1993, but the QUADPACK functions
964 were written some years before.
965 There is another version of QUADPACK at Netlib @footnote{@url{https://www.netlib.org/quadpack}};
966 it is not clear how that version differs from the SLATEC version.
968 The QUADPACK functions included in Maxima are all automatic, in the sense that
969 these functions attempt to compute a result to a specified accuracy, requiring
970 an unspecified number of function evaluations. Maxima's Lisp translation of
971 QUADPACK also includes some non-automatic functions, but they are not exposed
972 at the Maxima level.
974 Further information about QUADPACK can be found in the QUADPACK book
975 @footnote{R. Piessens, E. de Doncker-Kapenga, C.W. Uberhuber, and D.K. Kahaner.
976 @i{QUADPACK: A Subroutine Package for Automatic Integration.}
977 Berlin: Springer-Verlag, 1983, ISBN 0387125531.}.
979 @c -----------------------------------------------------------------------------
980 @subsection Overview
981 @c -----------------------------------------------------------------------------
983 @table @code
984 @item @mref{quad_qag}
985 Integration of a general function over a finite interval.
986 @mref{quad_qag} implements a simple globally adaptive integrator using the
987 strategy of Aind (Piessens, 1973).
988 The caller may choose among 6 pairs of Gauss-Kronrod quadrature
989 formulae for the rule evaluation component.
990 The high-degree rules are suitable for strongly oscillating integrands.
992 @item @mref{quad_qags}
993 Integration of a general function over a finite interval.
994 @mref{quad_qags} implements globally adaptive interval subdivision with
995 extrapolation (de Doncker, 1978) by the Epsilon algorithm (Wynn, 1956).
997 @item @mref{quad_qagi}
998 Integration of a general function over an infinite or semi-infinite interval.
999 The interval is mapped onto a finite interval and
1000 then the same strategy as in @code{quad_qags} is applied.
1002 @item @mref{quad_qawo}
1004 Integration of m4_math(<<<\cos(\omega x) f(x)>>>,<<<@math{cos(omega x) f(x)}>>>) or m4_math(<<<\sin(\omega x) f(x)>>>,<<<@math{sin(omega x) f(x)}>>>) over a
1005 finite interval, where m4_math(\omega, @math{omega}) is a constant.
1006 The rule evaluation component is based on the modified Clenshaw-Curtis
1007 technique. @mref{quad_qawo} applies adaptive subdivision with extrapolation,
1008 similar to @mrefdot{quad_qags}
1010 @item @mref{quad_qawf}
1011 Calculates a Fourier cosine or Fourier sine transform on a semi-infinite
1012 interval. The same approach as in @mref{quad_qawo} is applied on successive
1013 finite intervals, and convergence acceleration by means of the Epsilon algorithm
1014 (Wynn, 1956) is applied to the series of the integral contributions.
1016 @item @mref{quad_qaws}
1017 Integration of m4_math(w(x)f(x), @math{w(x) f(x)}) over a finite interval @math{[a, b]}, where
1018 @math{w} is a function of the form m4_math((x-a)^\alpha (b-x)^\beta v(x), @math{(x - a)^alpha (b - x)^beta v(x)}) and
1019 @math{v(x)} is 1 or m4_math(\log(x-a), @math{log(x - a)}) or m4_math(\log(b-x), @math{log(b - x)}) or
1020 m4_math(\log(x-a)\log(b-x), @math{log(x - a) log(b - x)}), and m4_math(\alpha > -1, @math{alpha > -1}) and m4_math(\beta > -1, @math{beta > -1}).
1022 A globally adaptive subdivision strategy is applied, with modified
1023 Clenshaw-Curtis integration on the subintervals which contain @math{a}
1024 or @math{b}.
1026 @item @mref{quad_qawc}
1027 Computes the Cauchy principal value of @math{f(x)/(x - c)} over a finite
1028 interval @math{(a, b)} and specified @math{c}.
1029 The strategy is globally adaptive, and modified
1030 Clenshaw-Curtis integration is used on the subranges
1031 which contain the point @math{x = c}.
1033 @item @mref{quad_qagp}
1034 Basically the same as @mref{quad_qags} but points of singularity or
1035 discontinuity of the integrand must be supplied. This makes it easier
1036 for the integrator to produce a good solution.
1037 @end table
1040 @opencatbox{Categories:}
1041 @category{Integral calculus}
1042 @category{Numerical methods}
1043 @category{Share packages}
1044 @category{Package quadpack}
1045 @closecatbox
1047 @c -----------------------------------------------------------------------------
1048 @node Functions and Variables for QUADPACK, , Introduction to QUADPACK, Integration
1049 @section Functions and Variables for QUADPACK
1050 @c -----------------------------------------------------------------------------
1052 @c THERE ARE OPTIONAL ARGUMENTS WHICH MAKES LISTING THE VARIANTS A LITTLE TEDIOUS
1053 @c NEED A MORE CONVENIENT (AND NONAMBIGUOUS) NOTATION FOR OPTIONAL ARGUMENTS
1055 @c -----------------------------------------------------------------------------
1056 @anchor{quad_qag}
1057 @deffn {Function} quad_qag @
1058 @fname{quad_qag} (@var{f(x)}, @var{x}, @var{a}, @var{b}, @var{key}, [@var{epsrel}, @var{epsabs}, @var{limit}]) @
1059 @fname{quad_qag} (@var{f}, @var{x}, @var{a}, @var{b}, @var{key}, [@var{epsrel}, @var{epsabs}, @var{limit}])
1061 Integration of a general function over a finite interval. @code{quad_qag}
1062 implements a simple globally adaptive integrator using the strategy of Aind
1063 (Piessens, 1973). The caller may choose among 6 pairs of Gauss-Kronrod
1064 quadrature formulae for the rule evaluation component. The high-degree rules
1065 are suitable for strongly oscillating integrands.
1067 @code{quad_qag} computes the integral
1069 m4_displaymath(
1070 <<<\int_a^b f(x)\, dx>>>,
1071 @math{integrate (f(x), x, a, b)}>>>)
1073 The function to be integrated is @math{f(x)}, with dependent
1074 variable @math{x}, and the function is to be integrated between the
1075 limits @math{a} and @math{b}. @var{key} is the integrator to be used
1076 and should be an integer between 1 and 6, inclusive. The value of
1077 @var{key} selects the order of the Gauss-Kronrod integration rule.
1078 High-order rules are suitable for strongly oscillating integrands.
1080 The integrand may be specified as the name of a Maxima or Lisp function or
1081 operator, a Maxima lambda expression, or a general Maxima expression.
1083 The numerical integration is done adaptively by subdividing the
1084 integration region into sub-intervals until the desired accuracy is
1085 achieved.
1087 The keyword arguments are optional and may be specified in any order.
1088 They all take the form @code{key=val}. The keyword arguments are:
1090 @table @code
1091 @item epsrel
1092 Desired relative error of approximation. Default is 1d-8.
1093 @item epsabs
1094 Desired absolute error of approximation. Default is 0.
1095 @item limit
1096 Size of internal work array. @var{limit} is the
1097 maximum number of subintervals to use. Default is 200.
1098 @end table
1100 @code{quad_qag} returns a list of four elements:
1102 @itemize
1103 @item
1104 an approximation to the integral,
1105 @item
1106 the estimated absolute error of the approximation,
1107 @item
1108 the number integrand evaluations,
1109 @item
1110 an error code.
1111 @end itemize
1113 The error code (fourth element of the return value) can have the values:
1115 @table @code
1116 @item 0
1117 if no problems were encountered;
1118 @item 1
1119 if too many sub-intervals were done;
1120 @item 2
1121 if excessive roundoff error is detected;
1122 @item 3
1123 if extremely bad integrand behavior occurs;
1124 @item 6
1125 if the input is invalid.
1127 @end table
1129 @c NEED CROSS REFS HERE -- EITHER CROSS REF A QUADPACK OVERVIEW, OR CROSS REF EACH OF THE quad_* FUNCTIONS
1131 Examples:
1133 @c ===beg===
1134 @c quad_qag (x^(1/2)*log(1/x), x, 0, 1, 3, 'epsrel=5d-8);
1135 @c integrate (x^(1/2)*log(1/x), x, 0, 1);
1136 @c ===end===
1137 @example
1138 @group
1139 (%i1) quad_qag (x^(1/2)*log(1/x), x, 0, 1, 3, 'epsrel=5d-8);
1140 (%o1) [.4444444444492108, 3.1700968502883E-9, 961, 0]
1141 @end group
1142 @group
1143 (%i2) integrate (x^(1/2)*log(1/x), x, 0, 1);
1144 4
1145 (%o2) -
1146 9
1147 @end group
1148 @end example
1150 @opencatbox{Categories:}
1151 @category{Numerical methods}
1152 @category{Package quadpack}
1153 @closecatbox
1154 @end deffn
1156 @c THERE ARE OPTIONAL ARGUMENTS WHICH MAKES LISTING THE VARIANTS A LITTLE TEDIOUS
1157 @c NEED A MORE CONVENIENT (AND NONAMBIGUOUS) NOTATION FOR OPTIONAL ARGUMENTS
1159 @c -----------------------------------------------------------------------------
1160 @anchor{quad_qags}
1161 @deffn {Function} quad_qags @
1162 @fname{quad_qags} (@var{f(x)}, @var{x}, @var{a}, @var{b}, [@var{epsrel}, @var{epsabs}, @var{limit}]) @
1163 @fname{quad_qags} (@var{f}, @var{x}, @var{a}, @var{b}, [@var{epsrel}, @var{epsabs}, @var{limit}])
1165 Integration of a general function over a finite interval.
1166 @code{quad_qags} implements globally adaptive interval subdivision with
1167 extrapolation (de Doncker, 1978) by the Epsilon algorithm (Wynn, 1956).
1169 @code{quad_qags} computes the integral
1171 m4_displaymath(
1172 <<<\int_a^b f(x)\, dx>>>,
1173 @math{integrate (f(x), x, a, b)}>>>)
1175 The function to be integrated is @math{f(x)}, with
1176 dependent variable @math{x}, and the function is to be integrated
1177 between the limits @math{a} and @math{b}.
1179 The integrand may be specified as the name of a Maxima or Lisp function or
1180 operator, a Maxima lambda expression, or a general Maxima expression.
1182 The keyword arguments are optional and may be specified in any order.
1183 They all take the form @code{key=val}. The keyword arguments are:
1185 @table @code
1186 @item epsrel
1187 Desired relative error of approximation. Default is 1d-8.
1188 @item epsabs
1189 Desired absolute error of approximation. Default is 0.
1190 @item limit
1191 Size of internal work array. @var{limit} is the
1192 maximum number of subintervals to use. Default is 200.
1193 @end table
1195 @code{quad_qags} returns a list of four elements:
1197 @itemize
1198 @item
1199 an approximation to the integral,
1200 @item
1201 the estimated absolute error of the approximation,
1202 @item
1203 the number integrand evaluations,
1204 @item
1205 an error code.
1206 @end itemize
1208 The error code (fourth element of the return value) can have the values:
1210 @table @code
1211 @item 0
1212 no problems were encountered;
1213 @item 1
1214 too many sub-intervals were done;
1215 @item 2
1216 excessive roundoff error is detected;
1217 @item 3
1218 extremely bad integrand behavior occurs;
1219 @item 4
1220 failed to converge
1221 @item 5
1222 integral is probably divergent or slowly convergent
1223 @item 6
1224 if the input is invalid.
1225 @end table
1227 @c NEED CROSS REFS HERE -- EITHER CROSS REF A QUADPACK OVERVIEW, OR CROSS REF EACH OF THE quad_* FUNCTIONS
1229 Examples:
1231 @c ===beg===
1232 @c quad_qags (x^(1/2)*log(1/x), x, 0, 1, 'epsrel=1d-10);
1233 @c ===end===
1234 @example
1235 @group
1236 (%i1) quad_qags (x^(1/2)*log(1/x), x, 0, 1, 'epsrel=1d-10);
1237 (%o1) [.4444444444444448, 1.11022302462516E-15, 315, 0]
1238 @end group
1239 @end example
1241 Note that @code{quad_qags} is more accurate and efficient than @code{quad_qag} for this integrand.
1243 @opencatbox{Categories:}
1244 @category{Numerical methods}
1245 @category{Package quadpack}
1246 @closecatbox
1247 @end deffn
1249 @c THERE ARE OPTIONAL ARGUMENTS WHICH MAKES LISTING THE VARIANTS A LITTLE TEDIOUS
1250 @c NEED A MORE CONVENIENT (AND NONAMBIGUOUS) NOTATION FOR OPTIONAL ARGUMENTS
1252 @c -----------------------------------------------------------------------------
1253 @anchor{quad_qagi}
1254 @deffn {Function} quad_qagi @
1255 @fname{quad_qagi} (@var{f(x)}, @var{x}, @var{a}, @var{b}, [@var{epsrel}, @var{epsabs}, @var{limit}]) @
1256 @fname{quad_qagi} (@var{f}, @var{x}, @var{a}, @var{b}, [@var{epsrel}, @var{epsabs}, @var{limit}])
1258 Integration of a general function over an infinite or semi-infinite interval.
1259 The interval is mapped onto a finite interval and
1260 then the same strategy as in @code{quad_qags} is applied.
1262 @code{quad_qagi} evaluates one of the following integrals
1264 m4_displaymath(
1265 <<<\int_a^\infty f(x) \, dx>>>,
1266 <<<@math{integrate (f(x), x, a, inf)}>>>)
1268 m4_displaymath(
1269 <<<\int_\infty^a f(x) \, dx>>>,
1270 <<<@math{integrate (f(x), x, minf, a)}>>>)
1272 m4_displaymath(
1273 <<<\int_{-\infty}^\infty f(x) \, dx>>>,
1274 <<<@math{integrate (f(x), x, minf, inf)}>>>)
1276 using the Quadpack QAGI routine. The function to be integrated is
1277 @math{f(x)}, with dependent variable @math{x}, and the function is to
1278 be integrated over an infinite range.
1280 The integrand may be specified as the name of a Maxima or Lisp function or
1281 operator, a Maxima lambda expression, or a general Maxima expression.
1283 One of the limits of integration must be infinity. If not, then
1284 @code{quad_qagi} will just return the noun form.
1286 The keyword arguments are optional and may be specified in any order.
1287 They all take the form @code{key=val}. The keyword arguments are:
1289 @table @code
1290 @item epsrel
1291 Desired relative error of approximation. Default is 1d-8.
1292 @item epsabs
1293 Desired absolute error of approximation. Default is 0.
1294 @item limit
1295 Size of internal work array. @var{limit} is the
1296 maximum number of subintervals to use. Default is 200.
1297 @end table
1299 @code{quad_qagi} returns a list of four elements:
1301 @itemize
1302 @item
1303 an approximation to the integral,
1304 @item
1305 the estimated absolute error of the approximation,
1306 @item
1307 the number integrand evaluations,
1308 @item
1309 an error code.
1310 @end itemize
1312 The error code (fourth element of the return value) can have the values:
1314 @table @code
1315 @item 0
1316 no problems were encountered;
1317 @item 1
1318 too many sub-intervals were done;
1319 @item 2
1320 excessive roundoff error is detected;
1321 @item 3
1322 extremely bad integrand behavior occurs;
1323 @item 4
1324 failed to converge
1325 @item 5
1326 integral is probably divergent or slowly convergent
1327 @item 6
1328 if the input is invalid.
1330 @end table
1332 @c NEED CROSS REFS HERE -- EITHER CROSS REF A QUADPACK OVERVIEW, OR CROSS REF EACH OF THE quad_* FUNCTIONS
1334 Examples:
1336 @c ===beg===
1337 @c quad_qagi (x^2*exp(-4*x), x, 0, inf, 'epsrel=1d-8);
1338 @c integrate (x^2*exp(-4*x), x, 0, inf);
1339 @c ===end===
1340 @example
1341 @group
1342 (%i1) quad_qagi (x^2*exp(-4*x), x, 0, inf, 'epsrel=1d-8);
1343 (%o1) [0.03125, 2.95916102995002E-11, 105, 0]
1344 @end group
1345 @group
1346 (%i2) integrate (x^2*exp(-4*x), x, 0, inf);
1347 1
1348 (%o2) --
1349 32
1350 @end group
1351 @end example
1353 @opencatbox{Categories:}
1354 @category{Numerical methods}
1355 @category{Package quadpack}
1356 @closecatbox
1357 @end deffn
1359 @c THERE ARE OPTIONAL ARGUMENTS WHICH MAKES LISTING THE VARIANTS A LITTLE TEDIOUS
1360 @c NEED A MORE CONVENIENT (AND NONAMBIGUOUS) NOTATION FOR OPTIONAL ARGUMENTS
1362 @c -----------------------------------------------------------------------------
1363 @anchor{quad_qawc}
1364 @deffn {Function} quad_qawc @
1365 @fname{quad_qawc} (@var{f(x)}, @var{x}, @var{c}, @var{a}, @var{b}, [@var{epsrel}, @var{epsabs}, @var{limit}]) @
1366 @fname{quad_qawc} (@var{f}, @var{x}, @var{c}, @var{a}, @var{b}, [@var{epsrel}, @var{epsabs}, @var{limit}])
1368 Computes the Cauchy principal value of @math{f(x)/(x - c)} over a finite
1369 interval. The strategy is globally adaptive, and modified
1370 Clenshaw-Curtis integration is used on the subranges
1371 which contain the point @math{x = c}.
1373 @code{quad_qawc} computes the Cauchy principal value of
1375 m4_displaymath(
1376 <<<\int_{a}^{b}{{{f\left(x\right)}\over{x-c}}\>dx}>>>,
1377 <<<@math{integrate (f(x)/(x - c), x, a, b)}>>>)
1379 using the Quadpack QAWC routine. The function to be integrated is
1380 @math{f(x)/(x-c)}, with dependent variable @math{x}, and the
1381 function is to be integrated over the interval @math{a} to @math{b}.
1383 The integrand may be specified as the name of a Maxima or Lisp function or
1384 operator, a Maxima lambda expression, or a general Maxima expression.
1386 The keyword arguments are optional and may be specified in any order.
1387 They all take the form @code{key=val}. The keyword arguments are:
1389 @table @code
1390 @item epsrel
1391 Desired relative error of approximation. Default is 1d-8.
1392 @item epsabs
1393 Desired absolute error of approximation. Default is 0.
1394 @item limit
1395 Size of internal work array. @var{limit} is the
1396 maximum number of subintervals to use. Default is 200.
1397 @end table
1399 @code{quad_qawc} returns a list of four elements:
1401 @itemize
1402 @item
1403 an approximation to the integral,
1404 @item
1405 the estimated absolute error of the approximation,
1406 @item
1407 the number integrand evaluations,
1408 @item
1409 an error code.
1410 @end itemize
1412 The error code (fourth element of the return value) can have the values:
1414 @table @code
1415 @item 0
1416 no problems were encountered;
1417 @item 1
1418 too many sub-intervals were done;
1419 @item 2
1420 excessive roundoff error is detected;
1421 @item 3
1422 extremely bad integrand behavior occurs;
1423 @item 6
1424 if the input is invalid.
1426 @end table
1428 Examples:
1430 @c ===beg===
1431 @c quad_qawc (2^(-5)*((x-1)^2+4^(-5))^(-1), x, 2, 0, 5, 'epsrel=1d-7);
1432 @c integrate (2^(-alpha)*(((x-1)^2 + 4^(-alpha))*(x-2))^(-1), x, 0, 5);
1433 @c ev (%, alpha=5, numer);
1434 @c ===end===
1435 @example
1436 @group
1437 (%i1) quad_qawc (2^(-5)*((x-1)^2+4^(-5))^(-1), x, 2, 0, 5,
1438 'epsrel=1d-7);
1439 (%o1) [- 3.130120337415925, 1.306830140249558E-8, 495, 0]
1440 @end group
1441 @group
1442 (%i2) integrate (2^(-alpha)*(((x-1)^2 + 4^(-alpha))*(x-2))^(-1),
1443 x, 0, 5);
1444 Principal Value
1445 alpha
1446 alpha 9 4 9
1447 4 log(------------- + -------------)
1448 alpha alpha
1449 64 4 + 4 64 4 + 4
1450 (%o2) (-----------------------------------------
1451 alpha
1452 2 4 + 2
1454 3 alpha 3 alpha
1455 ------- -------
1456 2 alpha/2 2 alpha/2
1457 2 4 atan(4 4 ) 2 4 atan(4 ) alpha
1458 - --------------------------- - -------------------------)/2
1459 alpha alpha
1460 2 4 + 2 2 4 + 2
1461 @end group
1462 @group
1463 (%i3) ev (%, alpha=5, numer);
1464 (%o3) - 3.130120337415917
1465 @end group
1466 @end example
1468 @opencatbox{Categories:}
1469 @category{Numerical methods}
1470 @category{Package quadpack}
1471 @closecatbox
1472 @end deffn
1474 @c THERE ARE OPTIONAL ARGUMENTS WHICH MAKES LISTING THE VARIANTS A LITTLE TEDIOUS
1475 @c NEED A MORE CONVENIENT (AND NONAMBIGUOUS) NOTATION FOR OPTIONAL ARGUMENTS
1477 @c -----------------------------------------------------------------------------
1478 @anchor{quad_qawf}
1479 @deffn {Function} quad_qawf @
1480 @fname{quad_qawf} (@var{f(x)}, @var{x}, @var{a}, @var{omega}, @var{trig}, [@var{epsabs}, @var{limit}, @var{maxp1}, @var{limlst}]) @
1481 @fname{quad_qawf} (@var{f}, @var{x}, @var{a}, @var{omega}, @var{trig}, [@var{epsabs}, @var{limit}, @var{maxp1}, @var{limlst}])
1483 Calculates a Fourier cosine or Fourier sine transform on a semi-infinite
1484 interval using the Quadpack QAWF function. The same approach as in
1485 @code{quad_qawo} is applied on successive finite intervals, and convergence
1486 acceleration by means of the Epsilon algorithm (Wynn, 1956) is applied to the
1487 series of the integral contributions.
1489 @code{quad_qawf} computes the integral
1491 m4_displaymath(
1492 <<<\int_a^\infty f(x) \, w(x) \, dx>>>,
1493 <<<@math{integrate (f(x)*w(x), x, a, inf)}>>>)
1495 The weight function @math{w} is selected by @var{trig}:
1497 @table @code
1498 @item cos
1499 m4_math(w(x) = \cos\omega x, @math{w(x) = cos (omega x)})
1500 @item sin
1501 m4_math(w(x) = \sin\omega x, @math{w(x) = sin (omega x)})
1502 @end table
1504 The integrand may be specified as the name of a Maxima or Lisp function or
1505 operator, a Maxima lambda expression, or a general Maxima expression.
1507 The keyword arguments are optional and may be specified in any order.
1508 They all take the form @code{key=val}. The keyword arguments are:
1510 @table @code
1511 @item epsabs
1512 Desired absolute error of approximation. Default is 1d-10.
1513 @item limit
1514 Size of internal work array. (@var{limit} - @var{limlst})/2 is the
1515 maximum number of subintervals to use. Default is 200.
1516 @item maxp1
1517 Maximum number of Chebyshev moments. Must be greater than 0. Default
1518 is 100.
1519 @item limlst
1520 Upper bound on the number of cycles. Must be greater than or equal to
1521 3. Default is 10.
1522 @end table
1524 @code{quad_qawf} returns a list of four elements:
1526 @itemize
1527 @item
1528 an approximation to the integral,
1529 @item
1530 the estimated absolute error of the approximation,
1531 @item
1532 the number integrand evaluations,
1533 @item
1534 an error code.
1535 @end itemize
1537 The error code (fourth element of the return value) can have the values:
1539 @table @code
1540 @item 0
1541 no problems were encountered;
1542 @item 1
1543 too many sub-intervals were done;
1544 @item 2
1545 excessive roundoff error is detected;
1546 @item 3
1547 extremely bad integrand behavior occurs;
1548 @item 6
1549 if the input is invalid.
1551 @end table
1553 Examples:
1555 @c ===beg===
1556 @c quad_qawf (exp(-x^2), x, 0, 1, 'cos, 'epsabs=1d-9);
1557 @c integrate (exp(-x^2)*cos(x), x, 0, inf);
1558 @c ev (%, numer);
1559 @c ===end===
1560 @example
1561 @group
1562 (%i1) quad_qawf (exp(-x^2), x, 0, 1, 'cos, 'epsabs=1d-9);
1563 (%o1) [.6901942235215714, 2.84846300257552E-11, 215, 0]
1564 @end group
1565 @group
1566 (%i2) integrate (exp(-x^2)*cos(x), x, 0, inf);
1567 - 1/4
1568 %e sqrt(%pi)
1569 (%o2) -----------------
1570 2
1571 @end group
1572 @group
1573 (%i3) ev (%, numer);
1574 (%o3) .6901942235215714
1575 @end group
1576 @end example
1578 @opencatbox{Categories:}
1579 @category{Numerical methods}
1580 @category{Package quadpack}
1581 @closecatbox
1582 @end deffn
1584 @c THERE ARE OPTIONAL ARGUMENTS WHICH MAKES LISTING THE VARIANTS A LITTLE TEDIOUS
1585 @c NEED A MORE CONVENIENT (AND NONAMBIGUOUS) NOTATION FOR OPTIONAL ARGUMENTS
1587 @c -----------------------------------------------------------------------------
1588 @anchor{quad_qawo}
1589 @deffn {Function} quad_qawo @
1590 @fname{quad_qawo} (@var{f(x)}, @var{x}, @var{a}, @var{b}, @var{omega}, @var{trig}, [@var{epsrel}, @var{epsabs}, @var{limit}, @var{maxp1}, @var{limlst}]) @
1591 @fname{quad_qawo} (@var{f}, @var{x}, @var{a}, @var{b}, @var{omega}, @var{trig}, [@var{epsrel}, @var{epsabs}, @var{limit}, @var{maxp1}, @var{limlst}])
1593 Integration of m4_math(<<<\cos(\omega x) f(x)>>>,<<<@math{cos (omega x) f(x)}>>>) or m4_math(<<<\sin(\omega x)>>>, <<<@math{sin (omega x) f(x)}>>>) over a finite interval,
1594 where m4_math(\omega, @math{omega}) is a constant.
1595 The rule evaluation component is based on the modified
1596 Clenshaw-Curtis technique. @code{quad_qawo} applies adaptive subdivision with
1597 extrapolation, similar to @code{quad_qags}.
1599 @code{quad_qawo} computes the integral using the Quadpack QAWO
1600 routine:
1602 m4_displaymath(
1603 <<<\int_a^b f(x) \, w(x) \, dx>>>,
1604 <<<@math{integrate (f(x)*w(x), x, a, b)}>>>)
1606 The weight function @math{w} is selected by @var{trig}:
1608 @table @code
1609 @item cos
1610 m4_math(w(x) = \cos\omega x, @math{w(x) = cos (omega x)})
1611 @item sin
1612 m4_math(w(x) = \sin\omega x, @math{w(x) = sin (omega x)})
1613 @end table
1615 The integrand may be specified as the name of a Maxima or Lisp function or
1616 operator, a Maxima lambda expression, or a general Maxima expression.
1618 The keyword arguments are optional and may be specified in any order.
1619 They all take the form @code{key=val}. The keyword arguments are:
1621 @table @code
1622 @item epsrel
1623 Desired relative error of approximation. Default is 1d-8.
1624 @item epsabs
1625 Desired absolute error of approximation. Default is 0.
1626 @item limit
1627 Size of internal work array. @var{limit}/2 is the
1628 maximum number of subintervals to use. Default is 200.
1629 @item maxp1
1630 Maximum number of Chebyshev moments. Must be greater than 0. Default
1631 is 100.
1632 @item limlst
1633 Upper bound on the number of cycles. Must be greater than or equal to
1634 3. Default is 10.
1635 @end table
1637 @code{quad_qawo} returns a list of four elements:
1639 @itemize
1640 @item
1641 an approximation to the integral,
1642 @item
1643 the estimated absolute error of the approximation,
1644 @item
1645 the number integrand evaluations,
1646 @item
1647 an error code.
1648 @end itemize
1650 The error code (fourth element of the return value) can have the values:
1652 @table @code
1653 @item 0
1654 no problems were encountered;
1655 @item 1
1656 too many sub-intervals were done;
1657 @item 2
1658 excessive roundoff error is detected;
1659 @item 3
1660 extremely bad integrand behavior occurs;
1661 @item 6
1662 if the input is invalid.
1664 @end table
1666 Examples:
1668 @c ===beg===
1669 @c quad_qawo (x^(-1/2)*exp(-2^(-2)*x), x, 1d-8, 20*2^2, 1, cos);
1670 @c rectform (integrate (x^(-1/2)*exp(-2^(-alpha)*x) * cos(x), x, 0, inf));
1671 @c input:pos;
1672 @c ev (%, alpha=2, numer);
1673 @c ===end===
1674 @example
1675 @group
1676 (%i1) quad_qawo (x^(-1/2)*exp(-2^(-2)*x), x, 1d-8, 20*2^2, 1, cos);
1677 (%o1) [1.376043389877692, 4.72710759424899E-11, 765, 0]
1678 @end group
1679 @group
1680 (%i2) rectform (integrate (x^(-1/2)*exp(-2^(-alpha)*x) * cos(x),
1681 x, 0, inf));
1682 alpha/2 - 1/2 2 alpha
1683 sqrt(%pi) 2 sqrt(sqrt(2 + 1) + 1)
1684 (%o2) -----------------------------------------------------
1685 2 alpha
1686 sqrt(2 + 1)
1687 @end group
1688 @group
1689 (%i3) ev (%, alpha=2, numer);
1690 (%o3) 1.376043390090716
1691 @end group
1692 @end example
1694 @opencatbox{Categories:}
1695 @category{Numerical methods}
1696 @category{Package quadpack}
1697 @closecatbox
1698 @end deffn
1700 @c THERE ARE OPTIONAL ARGUMENTS WHICH MAKES LISTING THE VARIANTS A LITTLE TEDIOUS
1701 @c NEED A MORE CONVENIENT (AND NONAMBIGUOUS) NOTATION FOR OPTIONAL ARGUMENTS
1703 @c -----------------------------------------------------------------------------
1704 @anchor{quad_qaws}
1705 @deffn {Function} quad_qaws @
1706 @fname{quad_qaws} (@var{f(x)}, @var{x}, @var{a}, @var{b}, @var{alpha}, @var{beta}, @var{wfun}, [@var{epsrel}, @var{epsabs}, @var{limit}]) @
1707 @fname{quad_qaws} (@var{f}, @var{x}, @var{a}, @var{b}, @var{alpha}, @var{beta}, @var{wfun}, [@var{epsrel}, @var{epsabs}, @var{limit}])
1709 Integration of @math{w(x) f(x)} over a finite interval, where @math{w(x)} is a
1710 certain algebraic or logarithmic function. A globally adaptive subdivision
1711 strategy is applied, with modified Clenshaw-Curtis integration on the
1712 subintervals which contain the endpoints of the interval of integration.
1714 @code{quad_qaws} computes the integral using the Quadpack QAWS routine:
1716 m4_displaymath(
1717 <<<\int_a^b f(x) \, w(x) \, dx>>>,
1718 <<<@math{integrate (f(x)*w(x), x, a, b)}>>>)
1720 The weight function @math{w} is selected by @var{wfun}:
1722 @table @code
1723 @item 1
1724 m4_math(w(x) = (x - a)^\alpha (b - x)^\beta, @math{w(x) = (x - a)^alpha (b - x)^beta})
1725 @item 2
1726 m4_math(w(x) = (x - a)^\alpha (b - x)^\beta \log(x - a),
1727 @math{w(x) = (x - a)^alpha (b - x)^beta log(x - a)})
1728 @item 3
1729 m4_math(w(x) = (x - a)^\alpha (b - x)^\beta \log(b - x), @math{w(x) = (x - a)^alpha (b - x)^beta log(b - x)})
1730 @item 4
1731 m4_math(w(x) = (x - a)^\alpha (b - x)^\beta \log(x - a) \log(b - x),
1732 @math{w(x) = (x - a)^alpha (b - x)^beta log(x - a) log(b - x)})
1733 @end table
1735 The integrand may be specified as the name of a Maxima or Lisp function or
1736 operator, a Maxima lambda expression, or a general Maxima expression.
1738 The keyword arguments are optional and may be specified in any order.
1739 They all take the form @code{key=val}. The keyword arguments are:
1741 @table @code
1742 @item epsrel
1743 Desired relative error of approximation. Default is 1d-8.
1744 @item epsabs
1745 Desired absolute error of approximation. Default is 0.
1746 @item limit
1747 Size of internal work array. @var{limit}is the
1748 maximum number of subintervals to use. Default is 200.
1749 @end table
1751 @code{quad_qaws} returns a list of four elements:
1753 @itemize
1754 @item
1755 an approximation to the integral,
1756 @item
1757 the estimated absolute error of the approximation,
1758 @item
1759 the number integrand evaluations,
1760 @item
1761 an error code.
1762 @end itemize
1764 The error code (fourth element of the return value) can have the values:
1766 @table @code
1767 @item 0
1768 no problems were encountered;
1769 @item 1
1770 too many sub-intervals were done;
1771 @item 2
1772 excessive roundoff error is detected;
1773 @item 3
1774 extremely bad integrand behavior occurs;
1775 @item 6
1776 if the input is invalid.
1778 @end table
1780 Examples:
1782 @c ===beg===
1783 @c quad_qaws (1/(x+1+2^(-4)), x, -1, 1, -0.5, -0.5, 1, 'epsabs=1d-9);
1784 @c integrate ((1-x*x)^(-1/2)/(x+1+2^(-alpha)), x, -1, 1);
1785 @c input:pos;
1786 @c ev (%, alpha=4, numer);
1787 @c ===end===
1788 @example
1789 @group
1790 (%i1) quad_qaws (1/(x+1+2^(-4)), x, -1, 1, -0.5, -0.5, 1,
1791 'epsabs=1d-9);
1792 (%o1) [8.750097361672832, 1.24321522715422E-10, 170, 0]
1793 @end group
1794 @group
1795 (%i2) integrate ((1-x*x)^(-1/2)/(x+1+2^(-alpha)), x, -1, 1);
1796 alpha
1797 Is 4 2 - 1 positive, negative, or zero?
1799 pos;
1800 alpha alpha
1801 2 %pi 2 sqrt(2 2 + 1)
1802 (%o2) -------------------------------
1803 alpha
1804 4 2 + 2
1805 @end group
1806 @group
1807 (%i3) ev (%, alpha=4, numer);
1808 (%o3) 8.750097361672829
1809 @end group
1810 @end example
1812 @opencatbox{Categories:}
1813 @category{Numerical methods}
1814 @category{Package quadpack}
1815 @closecatbox
1816 @end deffn
1818 @c THERE ARE OPTIONAL ARGUMENTS WHICH MAKES LISTING THE VARIANTS A LITTLE TEDIOUS
1819 @c NEED A MORE CONVENIENT (AND NONAMBIGUOUS) NOTATION FOR OPTIONAL ARGUMENTS
1821 @c -----------------------------------------------------------------------------
1822 @anchor{quad_qagp}
1823 @deffn {Function} quad_qagp @
1824 @fname{quad_qagp} (@var{f(x)}, @var{x}, @var{a}, @var{b}, @var{points}, [@var{epsrel}, @var{epsabs}, @var{limit}]) @
1825 @fname{quad_qagp} (@var{f}, @var{x}, @var{a}, @var{b}, @var{points}, [@var{epsrel}, @var{epsabs}, @var{limit}])
1827 Integration of a general function over a finite interval.
1828 @code{quad_qagp} implements globally adaptive interval subdivision with
1829 extrapolation (de Doncker, 1978) by the Epsilon algorithm (Wynn, 1956).
1831 @code{quad_qagp} computes the integral
1833 m4_displaymath(
1834 <<<\int_a^b f(x) \, dx>>>,
1835 <<<@math{integrate (f(x), x, a, b)}>>>)
1837 The function to be integrated is @math{f(x)}, with
1838 dependent variable @math{x}, and the function is to be integrated
1839 between the limits @math{a} and @math{b}.
1841 The integrand may be specified as the name of a Maxima or Lisp function or
1842 operator, a Maxima lambda expression, or a general Maxima expression.
1844 To help the integrator, the user must supply a list of points where
1845 the integrand is singular or discontinuous.
1847 The keyword arguments are optional and may be specified in any order.
1848 They all take the form @code{key=val}. The keyword arguments are:
1850 @table @code
1851 @item epsrel
1852 Desired relative error of approximation. Default is 1d-8.
1853 @item epsabs
1854 Desired absolute error of approximation. Default is 0.
1855 @item limit
1856 Size of internal work array. @var{limit} is the
1857 maximum number of subintervals to use. Default is 200.
1858 @end table
1860 @code{quad_qagp} returns a list of four elements:
1862 @itemize
1863 @item
1864 an approximation to the integral,
1865 @item
1866 the estimated absolute error of the approximation,
1867 @item
1868 the number integrand evaluations,
1869 @item
1870 an error code.
1871 @end itemize
1873 The error code (fourth element of the return value) can have the values:
1875 @table @code
1876 @item 0
1877 no problems were encountered;
1878 @item 1
1879 too many sub-intervals were done;
1880 @item 2
1881 excessive roundoff error is detected;
1882 @item 3
1883 extremely bad integrand behavior occurs;
1884 @item 4
1885 failed to converge
1886 @item 5
1887 integral is probably divergent or slowly convergent
1888 @item 6
1889 if the input is invalid.
1890 @end table
1892 @c NEED CROSS REFS HERE -- EITHER CROSS REF A QUADPACK OVERVIEW, OR CROSS REF EACH OF THE quad_* FUNCTIONS
1894 Examples:
1896 @c ===beg===
1897 @c quad_qagp(x^3*log(abs((x^2-1)*(x^2-2))),x,0,3,[1,sqrt(2)]);
1898 @c quad_qags(x^3*log(abs((x^2-1)*(x^2-2))), x, 0, 3);
1899 @c ===end===
1900 @example
1901 @group
1902 (%i1) quad_qagp(x^3*log(abs((x^2-1)*(x^2-2))),x,0,3,[1,sqrt(2)]);
1903 (%o1) [52.74074838347143, 2.6247632689546663e-7, 1029, 0]
1904 @end group
1905 @group
1906 (%i2) quad_qags(x^3*log(abs((x^2-1)*(x^2-2))), x, 0, 3);
1907 (%o2) [52.74074847951494, 4.088443219529836e-7, 1869, 0]
1908 @end group
1909 @end example
1911 The integrand has singularities at @code{1} and @code{sqrt(2)} so we supply
1912 these points to @code{quad_qagp}. We also note that @code{quad_qagp} is
1913 more accurate and more efficient that @mrefdot{quad_qags}
1915 @opencatbox{Categories:}
1916 @category{Numerical methods}
1917 @category{Package quadpack}
1918 @closecatbox
1919 @end deffn
1921 @c -----------------------------------------------------------------------------
1922 @anchor{quad_control}
1923 @deffn {Function} quad_control (@var{parameter}, [@var{value}])
1925 Control error handling for quadpack. The parameter should be one of
1926 the following symbols:
1928 @table @code
1929 @item current_error
1930 The current error number
1931 @item control
1932 Controls if messages are printed or not. If it is set to zero or
1933 less, messages are suppressed.
1934 @item max_message
1935 The maximum number of times any message is to be printed.
1936 @end table
1938 If @var{value} is not given, then the current value of the
1939 @var{parameter} is returned. If @var{value} is given, the value of
1940 @var{parameter} is set to the given value.
1942 @opencatbox{Categories:}
1943 @category{Numerical methods}
1944 @category{Package quadpack}
1945 @closecatbox
1946 @end deffn