CODEWITHSUNDEEP
rdoubleintegration
Zied
Zied

Reputation: 11

How to efficiently integrate this function?

I have written this R code to double integrate a function (where a boundary is infinite):

nom <- function(u){
  tau <- 90   
  r <- 0   
  alpha <- .2
  k <- 1
  f1 <- function(c){(tau-u)*(alpha*k^alpha)^2/(k+c)^(alpha+1)*
                      (1/(k+r*tau+u+c)^(alpha+1)+
                       1/(k+r*tau-u+c)^(alpha+1))}
  integrate(f1,lower=0,upper=Inf)$value }
D <- Vectorize(nom)
q1 <- integrate(D,lower=0,upper=tau)

I obtained the following message

Error in integrate(f1, lower = 0, upper = Inf) : 
  non-finite function value

I have remarked that when changing r above the 0 value there will be no problem. But my task consists of choosing r to be exactly zero. Please can any one help to integrate this function with r=0.

Upvotes: 1

Views: 116

Answers (1)

RHertel
RHertel

Reputation: 23798

You're trying to integrate numerically a function which, for r=0, diverges strongly at c=u-k and at c=-u-k. This won't work. I suggest that you try to change the parameters and/or integration range; or that you search in mathematical textbooks for an analytical solution (if it exists).

tau <- 90   
r <- 0   
alpha <- .2
k <- 1
u<-1 # I set this parameter to 1, the value is not decisive for the problem
f1 <- function(c){(tau-u)*(alpha*k^alpha)^2/(k+c)^(alpha+1)*
      (1/(k+r*tau+u+c)^(alpha+1)+
         1/(k+r*tau-u+c)^(alpha+1))}
plot(f1, xlab="c")

enter image description here

PS: According to WolframAlpha your integral does not coverge for r=0. In other words, it seems that there is no solution other than "infinity" in this case. enter image description here

Upvotes: 3

Related Questions