1

Prove that the Simpson Rule for integrating the function $f$ on the interval $[a,b]$: $$I_2(f)=\frac{b-a}{6}(f(a)+4f(\frac{a+b}{2})+f(b))$$ is actually the Gauss quadrature for the weight $g=1$ and $3$ nodes.

I can not get why when we have $3$ nodes in the Gauss quadrature the nodes have to be the ones from the formula above! Can anyone give me some hints on how to prove this?

user53969
  • 145

1 Answers1

2

Start from the definition of a Gauss quadrature:

$$\int_{-1}^1 f(x)\ dx \approx \sum_{i=1}^3 g_i f(x_i).$$

We perform a change of variables to change the range of integration, so

$$\int_a^b f(x)\ dx = \frac{b-a}{2}\int_{-1}^1 f\left(\frac{b-a}{2}z+\frac{b+a}{2}\right)\ dz.$$

Apply gaussian quadrature to the new integral on the right.

$$\frac{b-a}{2}\int_{-1}^1 f\left(\frac{b-a}{2}z+\frac{b+a}{2}\right)\ dz \approx \frac{b-a}{2}\sum_{i=1}^3 g_i f\left(\frac{b-a}{2}z_i+\frac{b+a}{2}\right).$$

With a three-point quadrature, our points are $x_1 = a, x_2 = \frac{a+b}{2}, x_3 = b$. These correspond to $z_1 = -1, z_2 = 0, z_3 = 1$.

Plugging these in, we find $$\begin{align*} \int_a^b f(x)\ dx \approx \frac{b-a}{2}& \left[g_1 f\left(\frac{b-a}{2}\cdot -1 + \frac{a+b}{2}\right)\right. \\ &+\left.g_2 f\left(\frac{b-a}{2}\cdot 0 + \frac{a+b}{2}\right)\right. \\ &+ \left.g_3 f\left(\frac{b-a}{2}\cdot 1 + \frac{a+b}{2}\right)\right]. \end{align*}$$

The right-hand side simplifies to $$\frac{b-a}{2}\left[ g_1 f(a) + g_2 f\left(\frac{a+b}{2}\right) + g_3 f(b)\right].$$

Now you can see the function arguments are in the right form, and now all you have to do is compute the weights, which is straightforward since $g(x) = 1$. You can use $g_i = \int_\alpha^\beta g(x) \prod_{\substack{1 \le j \le 3 \\ j \neq i}} \frac{x-x_j}{x_i-x_j}\ dx$. But note that we've already done the change of variables to shift $[a,b]$ to $[-1,1]$ in the quadrature, so $\alpha = -1$ and $\beta = 1$. It is very straightforward to compute the integrals.

For $g_1$, we have

$$g_1 = \int_{-1}^1 \frac{x-0}{-1-0}\cdot \frac{x-1}{-1 - 1}\ dx = \frac12\left[\left.\frac{x^3}{3}\right|_{-1}^1 - \left.\frac{x^2}{2}\right|_{-1}^1\right] = \frac12\left[\frac13+\frac13-\left(\frac12-\frac12\right)\right] = \frac13.$$

You should find that $g_2$ and $g_3$ are similar. The result follows.

Emily
  • 35,688
  • 6
  • 93
  • 141
  • you've found the coefficients for quadrature formula using points $a, \frac{a+b}{2}$ and $b$ as nodes. However, for the quadrature to be Gaussian, these points need to be the roots of orthogonal polynomial. See here, for example

    http://en.wikipedia.org/wiki/Gaussian_quadrature, the part about fundamental theorem. Am I not right?

    – tenpercent Dec 04 '13 at 17:32
  • I should say, you're right, but in a limited case. The polynomial $P_n$ is normalized so that $P_1 = 1$. But we need not choose this normalization. In any case, the roots for $n=3$ are symmetric, so we can just pick a normalization that re-scales them to where we want, namely $\pm 1, 0$, or $a,\frac12(a+b),b$, or whatever. In the end, the scaling will wash out. – Emily Dec 04 '13 at 17:41
  • Well, multiplying a polynomial by a constant doesn't seem to change either scalar product of polynomials being equal to 0 or roots of polynomials. To adjust outer nodes to $\pm 1$ you would have to widen the segment $[-1, 1]$. – tenpercent Dec 04 '13 at 17:46
  • Yes, that's true. I would go from $\int_a^b f(x)\ dx$ to $\int_{-1}^1 f(z)\ dz/dx\ dx$ to $\int_{\alpha}^\beta f(\zeta)\ d\zeta/dz\ dz/dx\ dx$ but in the end I still get something of the desired form. I'm just skipping all the work that's ultimately unnecessary. – Emily Dec 04 '13 at 17:51
  • I'm just trying to say that orthogonal polynomials in $L_2([a, b])$ are not orthogonal in $L_2([a',b'])$. – tenpercent Dec 04 '13 at 17:59
  • They are up to scalar change of variables. I'm not shifting domains, I'm changing variables. – Emily Dec 04 '13 at 18:07
  • For Legendre polynomials, you have $\int_{-1}^1 P_m(x)P_n(x) = \frac{2}{2n+1} \delta_{mn}$. Let $x = 2y$. $\int_{-1/2}^{1/2} P_m(y)P_n(y) dy = \frac{1}{2n+1}\delta_{mn}$. Doesn't change orthogonality, just changes constants and domains. – Emily Dec 04 '13 at 18:09
  • I might have multiplied instead of divided by 2 there, or vice versa. Whatever -- the point still stands that you can choose your Gauss quadrature to be exact for some re-scaled family of Legendre polynomials. You don't have to choose $P_n = 1$ and live with those weights. – Emily Dec 04 '13 at 18:13
  • Well. I've understood what I've been doing wrong. We don't have the segment fixed in the problem statement. So the quadrature is actually Gaussian, but the segment, on which we build orthogonal polynomials, is not $[a,b]$. My bad. – tenpercent Dec 04 '13 at 18:26