Numerical integration: Solve upper bound given integral Given the answer L to the integral of a function f(x), how can I solve for the upper bound b using numerical integration?
$$\displaystyle \int_0^{b} f(x) dx = L$$
The anti-derivative of f(x) is of course not trivially solvable (unless I've missed something - in my case it is the square root of a quartic function)
 A: You want to solve $F(b) = L$ for $b$ with $$F(b) = \int_0^bf(x)\,dx\,.$$
There are many methods to numerically solve an equation. As you know the derivative of $F$, I would suggest Newton's method.
Given a starting guess $b_0$, you simply iterate
$$b_{n+1} = b_n - \frac{F(b_n)}{F'(b_n)}=b_n - \frac{1}{f(b_n)} \left(\int_0^{b_n}f(x)\,dx -L\right)\,.$$
A: Well there are actually many ways of estimate your integral, provided that it is typical Riemann integral, you can use the Trapezoidal rule, Simpson's rule, Gauss quadrature. Wikipedia has info, https://en.wikipedia.org/wiki/Trapezoidal_rule but it is better if you consults a book on numerical analysis.
For example refer to section 3.2 of Numerical Analysis Walter Gautschi.
https://books.google.co.cr/books?id=-fgjJF9yAIwC&printsec=frontcover#v=onepage&q&f=false
A: Another method that does not require you to recompute all the work if your guess for $b$ was not accurate enough is following. Consider some grid $x_n = n h$ with stepsize $h$. Do the following steps
$$
J_0 = 0\\
J_{n+1} = J_n + \int_{x_n}^{x_{n+1}} f(x) dx.
$$
Here $J(x)$ is simply $\int_0^x f(t) dt$ and $J_n = J(x_n)$.
You may use any single-interval quadrature that you like, no need to split interval $[x_n, x_{n+1}]$ into smaller pieces.
After each step check that $J_{n+1} < L$. If it is not true, then you've localized $b$ in the $[x_n, x_{n+1}]$ interval. Indeed, $J(x_n) < L$ (otherwise you would stop at the previous step) and $J(x_{n+1}) \geq L$. The solution of the
$$
J(b) = L
$$
should be somewhere between the $x_n$ and $x_{n+1}$.
Now you need to solve
$$
\int_{x_n}^b f(x) dx = L - J_n.
$$
The interval $b - x_n$ is small (less than $h$) and you have a decent initial guess, that can be derived from a linear approximation:
$$
\frac{b_0 - x_n}{x_{n+1} - x_n} = \frac{L - J_n}{J_{n+1} - J_{n}}
$$
Let's approximate $\int_{x_n}^b f(x) dx$ with some quadrature formula with weights $w_i$ and abscissae $\xi_i \in [0, 1]$:
$$
\int_{x_n}^b f(x) dx \approx (b-x_n) \sum_{i=0}^M w_i f(x_n + \xi_i (b - x_n))
$$
Rewriting the equation in $F(b) = 0$ form we have:
$$
F(b) = J_n + (b-x_n)\sum_{i=0}^M w_i f(x_n + \xi_i (b - x_n)) - L = 0\\
F'(b) = \sum_{i=0}^M w_i f(x_n + \xi_i (b - x_n)) + (b - x_n)
\sum_{i=0}^M w_i \xi_i f'(x_n + \xi_i (b - x_n))
$$
and the Newton iterations for this equation will be
$$
b_{k+1} = b_{k} - \frac{F(b_k)}{F'(b_k)}.
$$
If exact value for $F'(b)$ is hard to compute, $F'(b) \approx f(b)$ can be used instead.
