Построение неполного эллиптического интеграла 1-го рода

Question

Я хотел установить небольшой фрейм данных, чтобы построить для себя несколько точек неполного эллиптического интеграла 1-го рода для разных значений амплитуды phi и модуля k. Функция для интеграции - 1/sqrt(1 - (k*sin(x))^2) между 0 и phi. Вот код, который я себе представил:

v.phi <- seq(0, 2*pi, 1)
n.phi <- length(v.phi)
v.k <- seq(-1, +1, 0.5)
n.k <- length(v.k)
k <- rep(v.k, each = n.phi, times = 1)
phi <- rep(v.phi, each = 1, times = n.k)
df <- data.frame(k, phi)
func <- function(x, k) 1/sqrt(1 - (k*sin(x))^2)
df$area <- integrate(func,lower=0, upper=df$phi, k=df$k)

Но это приводит к ошибкам, и я, очевидно, ошибаюсь при построении новой переменной df $ area ... Может кто-нибудь правильно меня описать?

Answer 1

Вы можете использовать mapply:

df$area <- mapply(function(phi,k){
  integrate(func, lower=0, upper=phi, k=k)$value
}, df$phi, df$k)

Однако это приводит к ошибке, поскольку существуют некоторые значения k, равные 1 или -1, в то время как допустимые значения -1 < k < 1. Вы не можете оценить этот интеграл для k = +/- 1.

Обратите внимание, что есть лучший способ оценить этот интеграл: неполная эллиптическая функция первого рода реализована в пакете gsl:

> integrate(func, lower=0, upper=6, k=0.5)$value
[1] 6.458877
> gsl::ellint_F(6, 0.5)
[1] 6.458877

Как я уже сказал, эта функция не определена для k=-1 или k=1:

> gsl::ellint_F(6, 1)
[1] NaN
> gsl::ellint_F(6, -1)
[1] NaN
> integrate(func, lower=0, upper=6, k=1)
Error in integrate(func, lower = 0, upper = 6, k = 1) : 
  non-finite function value