simpadpt | R Documentation |
Numerically evaluate an integral using adaptive Simpson's rule.
simpadpt(f, a, b, tol = 1e-6, ...)
f |
univariate function, the integrand. |
a , b |
lower limits of integration; must be finite. |
tol |
relative tolerance |
... |
additional arguments to be passed to |
Approximates the integral of the function f
from a to b to within
an error of tol
using recursive adaptive Simpson quadrature.
A numerical value or vector, the computed integral.
Based on code from the book by Quarteroni et al., with some tricks borrowed from Matlab and Octave.
Quarteroni, A., R. Sacco, and F. Saleri (2007). Numerical Mathematics. Second Edition, Springer-Verlag, Berlin Heidelberg.
quad
, simpson2d
myf <- function(x, n) 1/(x+n) # 0.0953101798043249 , log((n+1)/n) for n=10
simpadpt(myf, 0, 1, n = 10) # 0.095310179804535
## Dilogarithm function
flog <- function(t) log(1-t) / t # singularity at t=1, almost at t=0
dilog <- function(x) simpadpt(flog, x, 0, tol = 1e-12)
dilog(1) # 1.64493406685615
# 1.64493406684823 = pi^2/6
## Not run:
N <- 51
xs <- seq(-5, 1, length.out = N)
ys <- numeric(N)
for (i in 1:N) ys[i] <- dilog(xs[i])
plot(xs, ys, type = "l", col = "blue",
main = "Dilogarithm function")
grid()
## End(Not run)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.