QSIMP
The QSIMP function performs numerical integration of a function over the closed interval [A, B] using Simpson’s rule.
Examples
To integrate the SIMPSON function (listed above) over the interval [0, π/2] and print the result:
; Define lower limit of integration:
A = 0.0
; Define upper limit of integration:
B = !PI/2.0
PRINT, QSIMP('simpson', A, B)
IDL prints:
-0.479158
The exact solution can be found using the integration-by-parts formula:
FB = 4.*B*(B^2-7.)*SIN(B) - (B^4-14.*B^2+28.)*COS(B)
FA = 4.*A*(A^2-7.)*SIN(A) - (A^4-14.*A^2+28.)*COS(A)
exact = FB - FA
PRINT, exact
IDL prints:
-0.479156
Syntax
Result = QSIMP( Func, A, B [, /DOUBLE] [, EPS=value] [, JMAX=value] )
Return Value
The result will have the same structure as the smaller of A and B, and the resulting type will be single- or double-precision floating, depending on the input types.
Arguments
Func
A scalar string specifying the name of a user-supplied IDL function to be integrated. This function must accept a single scalar argument X and return a scalar result. It must be defined over the closed interval [A, B].
For example, if we wish to integrate the fourth-order polynomial
y = (x4 - 2x2) sin(x)
we define a function SIMPSON to express this relationship in the IDL language:
FUNCTION simpson, X
RETURN, (X^4 - 2.0 * X^2) * SIN(X)
END
Note: If QSIMP is complex then only the real part is used for the computation.
A
The lower limit of the integration. A can be either a scalar or an array.
B
The upper limit of the integration. B can be either a scalar or an array.
Note: If arrays are specified for A and B, then QSIMP integrates the user-supplied function over the interval [Ai, Bi] for each i. If either A or B is a scalar and the other an array, the scalar is paired with each array element in turn.
Keywords
DOUBLE
Set this keyword to force the computation to be done in double-precision arithmetic.
EPS
The desired fractional accuracy. For single-precision calculations, the default value is 1.0 x 10-6. For double-precision calculations, the default value is 1.0 x 10-12.
JMAX
2(JMAX - 1) is the maximum allowed number of steps. If not specified, a default of 20 is used.
Version History
4.0 |
Introduced |
Resources and References
QSIMP is based on the routine qsimp
described in section 4.2 of Numerical Recipes in C: The Art of Scientific Computing (Second Edition), published by Cambridge University Press, and is used by permission.