Function: dblint (f, r, s, a, b)
A double-integral routine which was
written in top-level Maxima and then translated and compiled to machine code.
Use load (dblint) to access this package. It uses the Simpson's rule method in both the
x and y directions to calculate
|
|
/b /s(x)| || | f(x,y) dy dx| |/a /r(x) |
The function f must be a
translated or compiled function of two variables, and r and s
must each be a translated or compiled function of one variable, while a
and b must be floating point numbers. The routine has two global
variables which determine the number of divisions of the x and y intervals: dblint_x and dblint_y, both of
which are initially 10, and can be changed independently to other integer
values (there are 2*dblint_x+1 points computed in the x direction, and 2*dblint_y+1 in the y
direction). The routine subdivides the X axis and then for each value of X it
first computes r(x) and s(x); then the Y axis between r(x) and s(x) is subdivided
and the integral along the Y axis is performed using Simpson's rule; then the
integral along the X axis is done using Simpson's rule with the function values
being the Y-integrals. This procedure may be numerically unstable for a great
variety of reasons, but is reasonably fast: avoid using it on highly
oscillatory functions and functions with singularities (poles or branch points
in the region). The Y integrals depend on how far apart r(x) and s(x) are, so if the
distance s(x) - r(x) varies rapidly
with X, there may be substantial errors arising from truncation with different
step-sizes in the various Y integrals. One can increase dblint_x and dblint_y in an
effort to improve the coverage of the region, at the expense of computation
time. The function values are not saved, so if the function is very time-consuming,
you will have to wait for re-computation if you change anything (sorry). It is
required that the functions f, r, and s be
either translated or compiled prior to calling dblint. This will result
in orders of magnitude speed improvement over interpreted code in many cases!
demo
(dblint) executes a demonstration of dblint applied to
an example problem.
Function: defint (expr, x, a, b)