I came to an interesting question that I cannot easily answer.
There is a uniform grid $x = {x_1, \dots, x_{i-1}, x_{i}, \dots, x_n}$. Values of two functions that are not known analytically are given on that grid, $f = f_1, \dots, f_n$ and $g = g_1, \dots, g_n$. The task is to compute
$$ I_i = \int_{x_{i-1}}^{x_i}\,f \cdot g\, \mathrm{d}x$$
on every segment. For that purpose we have to approximate locally $f$ and $g$. The approximation has to be polynomial. There are two strategies.
- First approximate $f$ by a polynomial $F_n$ of the order $n$, approximate $g$ by a polynomial $G_m$ of the order $m$ ($n$ may be $= m$, but it is not required). For example,
$$F_2 = \varphi_0 + \varphi_1 (x - x_{i-1}) + \varphi_2 (x - x_{i-1})^2$$
$$G_1 = \gamma_0 + \gamma_1 (x - x_{i-1})$$
where $\varphi$'s and $\gamma$'s are, for example, determined by appropriate Lagrange interpolation formula depending on the local stencils for $f$ and for $g$.
Then multiply analytically $F_n$ and $G_m$ (which is a polynomial of the order $n+m$), then integrate that polynomial analytically and evaluate it between the given boundaries.
- Alternatively, we can first multiply the values of $f$ and $g$ in every grid point, $h_i = f_i \cdot g_i$, then approximate the function $h$ with polynomial $H_o$, for example:
$$H_1 = \chi_0 + \chi_1 (x - x_{i-1})$$
then integrate this polynomial analytically and, finally, evaluate it again between the boundaries. If this $o = 1$, the resulting analytical expression for the integral is exactly the trapezoidal rule.
If we consider only the edge points of the segment $(i-1, i)$, in the first approach we use all 4 values ($f_{i-1}, f_i, g_{i-1}, g_i$) to define the second-order polynomial ($F_1 \cdot G_1$). In the second case by multiplying $f_i \cdot g_i$ and $f_{i-1} \cdot g_{i-1}$ we apparently reduce the number of data points from 4 to 2 and we use these 2 to define the first-order polynomial, i.e. the linear function $H_1$.
Of course, in principle, $m$, $n$ or $o$, may be larger than 1 leading to higher order approximating polynomials.
My questions:
a) Which of the two solutions would be your preference and why?
b) One would expect the solution 1) to have higher order of accuracy as the unknown function under the integral ($f\cdot g$) is approximated by a higher order polynomial. Would you agree with that?
c) Which solution would you expect to be more numerically stable?
d) Are you aware of any reference book dealing with this problem?
Many thanks.
Edit: I've corrected the title to reflect the question more accurately.