A benchmark of different languages that are potentially suitable to scientific computing with the adaptive trapezoid integration algorithm.
Currently implemented languages: Rust, C++, Scala, Haskell, Python
- rs: Rust
- cpp: C++
- hs: Haskell
- py: Python
- scala: Scala
The adaptive trapezoid quadrature method (i.e., the definite integration) works by dividing the integration interval iteratively (or in other words, recursively) and approximate the result by the summing areas of trapezoids of all the intervals.
The detailed algorithm is
- Setting up the function F to be integrated, setting the tolerant value
eps
. - Initialize a set (more than 1) of initial ticks, which defines the initial intervals. The x values of the initial ticks should be in increasing order. Say a set of points with X's=x_i (i=0,1,2,..n-1), then the corresponding initial intervals are
[x_0, x_1], [x_1, x_2], [x_2, x_3], ...[x_{n-2}, x_{n-1}]
- The result of any interval
[x_1, x_2]
is calculated asI(F, x_1, x_2)
, which is defined in following algorithm. - Sorting the result of each interval according to their
abs
in increasing order, and add them up.
The definition of I(F, x_1, x_2)
is
- Calculate the
diff=T(F, x_1, x_2)-T(F, x_1, (x_1+x_2)/2)-T(F, (x_1+x_2)/2, x_2)
, whereT(F, a, b)=(F(a)+F(b))*(b-a)/2
. - If
diff<eps/W*(x_2-x_1)
, go to 3, otherwise go to 4 whereW
is the width of the whole initial integration interval. - return
T(F, x_1, (x_1+x_2)/2)+T(F, (x_1+x_2)/2, x_2)
. - return
I(F, x_1, (x_1+x_2)/2)+I(F, (x_1+x_2)/2, x_2)
.