Nothing
### == expint: Exponential Integral and Incomplete Gamma Function ==
###
### Tests for the incomplete gamma function
###
### G(a,x) = int_x^infty t^{a-1} exp(-t) dt
###
### for a *real* and x >= 0.
###
### AUTHOR: Vincent Goulet <vincent.goulet@act.ulaval.ca>
## Load the package
library(expint)
## a > 0; direct link to the standard incomplete gamma function
x <- c(0.2, 2.5, 5, 8, 10)
a <- 1.2
stopifnot(exprs = {
identical(gammainc(a, x),
gamma(a) * pgamma(x, a, 1, lower = FALSE))
})
## a = 0; direct link to the exponential integral
x <- c(0.2, 2.5, 5, 8, 10)
a <- 0
stopifnot(exprs = {
identical(gammainc(a, x), expint(x))
identical(gammainc(a, x), expint_E1(x))
})
## a < 0; compare to the recursive formula
x <- c(0.2, 2.5, 5, 8, 10)
a <- c(-0.25, -1.2, -2)
stopifnot(exprs = {
all.equal(gammainc(a[1], x),
-(x^a[1] * exp(-x))/a[1] +
gamma(a[1] + 1) * pgamma(x, a[1] + 1, 1, lower = FALSE)/a[1])
all.equal(gammainc(a[2], x),
-(x^a[2] * exp(-x))/a[2] +
(-(x^(a[2] + 1) * exp(-x))/(a[2] + 1) +
gamma(a[2] + 2) * pgamma(x, a[2] + 2, 1, lower = FALSE)/(a[2] + 1))/a[2])
all.equal(
gammainc(a[3], x),
-(x^a[3] * exp(-x))/a[3] +
(-(x^(a[3] + 1) * exp(-x))/(a[3] + 1) + expint_E1(x)/(a[3] + 1))/a[3])
})
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.