> Dear all, > I have some issues with integrate in R thus I would like to request > your help. I am trying to calculate the integral of f(x)*g(x). > The f(x) is a step function while g(x) is a polynomial. > If f(x) (step function) changes its value only few times (5 or 6 'steps') > everything is calulated ok(verified results in scrap paper) but if f(x) > takes like 800 different values I receive the error > > "Error in integrate number of subdivisions reached" > > I did some checks on the internet and I found that I can increase the > number of subdivisions (as this is a parameter in integrate(). > Thus I raised it from 100 to 1000 (and then to 10000). > > A. Does this makes the error produced higher or does it only stress > the computer? > B. When the number was raised to 10.000 I started getting the error message > "roundoff error was detected" > > What do you think I should do to solve that? > I would like to thank u in advance for your help > > Best Regards > Alex
There's obviously a more numerically stable approach. If g(x) is a polynomial you do know its polynomial antiderivative. Take that and sum up all intervals where the step function is constant. Example: g(x) = 1 constant, integrate x^2 over 0..1 in 10000 subdivisions. The antiderivative is 1/3 * x^3, thus g <- function(x) 1 x <- seq(0, 1, len=10001) sum((x[2:10001]^3 - x[1:10000]^3)*g(2:10001))/3 #=> 0.3333333 The antiderivative of a polynomial a_1 x^n + ... + a_{n+1} given as vector c(a1,...) can also be determined automatically, no manual work is necessary. Hans Werner ______________________________________________ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.