title |
---|
quadrature |
[TOC]
Experimental
Returns the trapezoidal rule integral of an array y
representing discrete samples of a function. The integral is computed assuming either equidistant abscissas with spacing dx
or arbitrary abscissas x
.
result =
[[stdlib_quadrature(module):trapz(interface)]] (y, x)
result =
[[stdlib_quadrature(module):trapz(interface)]] (y, dx)
y
: Shall be a rank-one array of type real
.
x
: Shall be a rank-one array of type real
having the same kind and size as y
.
dx
: Shall be a scalar of type real
having the same kind as y
.
The result is a scalar of type real
having the same kind as y
.
If the size of y
is zero or one, the result is zero.
{!example/quadrature/example_trapz.f90!}
Experimental
Given an array of abscissas x
, computes the array of weights w
such that if y
represented function values tabulated at x
, then sum(w*y)
produces a trapezoidal rule approximation to the integral.
result =
[[stdlib_quadrature(module):trapz_weights(interface)]] (x)
x
: Shall be a rank-one array of type real
.
The result is a real
array with the same size and kind as x
.
If the size of x
is one, then the sole element of the result is zero.
{!example/quadrature/example_trapz_weights.f90!}
Experimental
Returns the Simpson's rule integral of an array y
representing discrete samples of a function. The integral is computed assuming either equidistant abscissas with spacing dx
or arbitrary abscissas x
.
Simpson's ordinary ("1/3") rule is used for odd-length arrays. For even-length arrays, Simpson's 3/8 rule is also utilized in a way that depends on the value of even
. If even
is negative (positive), the 3/8 rule is used at the beginning (end) of the array. If even
is zero or not present, the result is as if the 3/8 rule were first used at the beginning of the array, then at the end of the array, and these two results were averaged.
result =
[[stdlib_quadrature(module):simps(interface)]] (y, x [, even])
result =
[[stdlib_quadrature(module):simps(interface)]] (y, dx [, even])
y
: Shall be a rank-one array of type real
.
x
: Shall be a rank-one array of type real
having the same kind and size as y
.
dx
: Shall be a scalar of type real
having the same kind as y
.
even
: (Optional) Shall be a default-kind integer
.
The result is a scalar of type real
having the same kind as y
.
If the size of y
is zero or one, the result is zero.
If the size of y
is two, the result is the same as if trapz
had been called instead.
{!example/quadrature/example_simps.f90!}
Experimental
Given an array of abscissas x
, computes the array of weights w
such that if y
represented function values tabulated at x
, then sum(w*y)
produces a Simpson's rule approximation to the integral.
Simpson's ordinary ("1/3") rule is used for odd-length arrays. For even-length arrays, Simpson's 3/8 rule is also utilized in a way that depends on the value of even
. If even
is negative (positive), the 3/8 rule is used at the beginning (end) of the array and the 1/3 rule used elsewhere. If even
is zero or not present, the result is as if the 3/8 rule were first used at the beginning of the array, then at the end of the array, and then these two results were averaged.
result =
[[stdlib_quadrature(module):simps_weights(interface)]] (x [, even])
x
: Shall be a rank-one array of type real
.
even
: (Optional) Shall be a default-kind integer
.
The result is a real
array with the same size and kind as x
.
If the size of x
is one, then the sole element of the result is zero.
If the size of x
is two, then the result is the same as if trapz_weights
had been called instead.
{!example/quadrature/example_simps_weights.f90!}
Experimental
Computes Gauss-Legendre quadrature (also known as simply Gaussian quadrature) nodes and weights,
for any N
(number of nodes).
Using the nodes x
and weights w
, you can compute the integral of some function f
as follows:
integral = sum(f(x) * w)
.
Only double precision is supported - if lower precision is required, you must do the appropriate conversion yourself. Accuracy has been validated up to N=64 by comparing computed results to tablulated values known to be accurate to machine precision (maximum difference from those values is 2 epsilon).
subroutine
[[stdlib_quadrature(module):gauss_legendre(interface)]] (x, w[, interval])
x
: Shall be a rank-one array of type real(real64)
. It is an output argument, representing the quadrature nodes.
w
: Shall be a rank-one array of type real(real64)
, with the same dimension as x
.
It is an output argument, representing the quadrature weights.
interval
: (Optional) Shall be a two-element array of type real(real64)
.
If present, the nodes and weigts are calculated for integration from interval(1)
to interval(2)
.
If not specified, the default integral is -1 to 1.
{!example/quadrature/example_gauss_legendre.f90!}
Experimental
Computes Gauss-Legendre-Lobatto quadrature nodes and weights,
for any N
(number of nodes).
Using the nodes x
and weights w
, you can compute the integral of some function f
as follows:
integral = sum(f(x) * w)
.
Only double precision is supported - if lower precision is required, you must do the appropriate conversion yourself. Accuracy has been validated up to N=64 by comparing computed results to tablulated values known to be accurate to machine precision (maximum difference from those values is 2 epsilon).
subroutine
[[stdlib_quadrature(module):gauss_legendre_lobatto(interface)]] (x, w[, interval])
x
: Shall be a rank-one array of type real(real64)
. It is an output argument, representing the quadrature nodes.
w
: Shall be a rank-one array of type real(real64)
, with the same dimension as x
.
It is an output argument, representing the quadrature weights.
interval
: (Optional) Shall be a two-element array of type real(real64)
.
If present, the nodes and weigts are calculated for integration from interval(1)
to interval(2)
.
If not specified, the default integral is -1 to 1.
{!example/quadrature/example_gauss_legendre_lobatto.f90!}