Numerical Solution to Inverse Integral (Pseudo Random Number Generation) If I have an arbitrary positive monotonically decreasing function $f(x), x \in [0,\infty]$, is there an 'efficient' method for finding $y$ in: 
$r = \int\limits_0^y f(x) dx $
for a known $r \in [0, \int\limits_0^\infty f(x) dx]$. By efficient I guess I mean more efficient than doing numerical integration until one lands in within a given distance from $r$. 
I particularly care about the cases where the integral of $f$ has no closed form.

Best answer so far (rewrite of below), doubt it gets much better than this:
Start with $s_0 = r$, $y_0 = 0$, then
$y_{i+1} = y_i + \frac{s_i}{f(y_i)} $
$s_{i+1} = s_i - \int\limits_{y_i}^{y_{i+1}}f(x)dx$
 A: I guess the only non-trivial thing about the problem is that:
$$
 x f(0) \geq \int_0^{x} f(t) dt \geq x f(x).
$$
So you start by computing the integral 
$$
 r_1 = \int_0^{y_1} f(t) dt,\quad y_1 = \frac{r}{f(0)}.
$$
Then replace $r$ by $r - r_1$. I think under reasonable assumptions this should converge pretty quickly (and always lower bound). It should be noted that I only use one of the inqualities, one can probably optimize it by using the other one.
A: There's always the option of expanding $f$ as a series, integrating that, and performing Lagrange inversion to obtain a series for the inverse function, from which you can determine suitable function approximations.
But for the purpose of generating non-uniform pseudorandom numbers, transformation isn't always the best choice; you might have better luck with the rejection approach. I suppose it depends on what distribution you're dealing with.
A: Luc Devroye has written some timeago a superb  book "Non-uniform random variate generation".
The whole book is available for free on his webpage.
A: Firstly this problem is not well-defined as stated. Let $f(x)=e^{-100x}$ and let $r=f(0)=1$. Since, $\int_0^{\infty}f(x)dx=\frac{1}{100}<1$ there is no $y$ that solves your problem. I am suspecting that you have $f$ so you should convince yourself that the problem actually has a solution. 
I always solve this sort of problem by root finding:
Let $H(y)=\int_0^y f(x)dx-r$ then the problem becomes finding a root of $H$. Now there are lots of software packages to solve this, these will give you much better accuracy than trial and error.  
In Matlab, the function fzero(fun,x0) will solve this problem. 
A: Newton's method? The derivative should be fairly straightforward to compute...
