Math428: Problem Set #1

1.

Consider a refinement study where the exact solution, I, is not known. Suppose one suspects that an approximation I_n that uses npanels has an error bound of the form

$$\vert I - I_n\vert \leq \frac{C}{n^p}$$

for some integer p. Prove that one can perform a refinement study without knowing I by using successive terms. That is, prove that

$$\vert I_n - I_{n/2}\vert \leq \frac{D}{n^p}.$$

How are C and D related to one another?

This problem highlights a basic theme in bounding errors: You begin with the term to be bounded and use the given information.

$$\begin{align*}\vert I_n - I_{n/2}\vert & = \vert I_n \textcolor{blue}{- I + I} -... ...2^p)}& \text{\textcolor{blue}{So we see that $D = C(1+2^p)$ .}} \\ \end{align*}$$

2.

Use midpoint rule estimate

$$\int_{\frac{1}{10}}^1 \frac{\cos\left(\frac{\ln(x)}{x}\right)}{x} dx$$

with 100, 200, 400 and 800 panels. Verify through this refinement that you are achieving the correct convergence rate. (This problem is related to the first ``Hundred Digit Challenge'' problems, but $\epsilon$ is fixed at 0.1.) Estimate the number of panels required to compute the first eight digits accurately.

See the solution to problem 3.

3.

Repeat problem 2 using Simpson's Rule.

This is just a question of computing things. To compute Simpson's rule, I took the intermediate step of using the trapezoid rule because it is conceptually simpler to do this. All but the last line of the table was posted in the lecture notes, so I complete the table here.

N	MN	M2N-MN	TN	T2N-TN	SN	S2N-SN
100	2.8664e-01	8.8986e-03	3.1909e-01	-1.6223e-03	2.9746e-01	5.2458e-04
200	2.9554e-01	1.8820e-03	3.0286e-01	-3.6624e-03	2.9798e-01	3.3855e-05
400	2.9742e-01	4.4824e-04	2.9920e-01	-8.9022e-04	2.9801e-01	2.0871e-06
800	2.9787e-01	1.1069e-04	2.9832e-01	-2.2099e-04	2.9802e-01	1.2988e-07

Here, we see that Simpson's rule picks up more than a digit of accuracy at every refinement while the other require more than two refinements to capture one more digit of accuracy.

4.

Prove that the Trapezoid Rule

$$I \approx I_n = \frac{1}{2} f(x_1) h + \sum_{i=2}^{n} f(x_i) h + \frac{1}{2} f(x_{n+1}) h$$

has an error bound

$$\vert I - I_n\vert \leq \frac{M(b-a)^3}{12 n^2}$$

where $\vert f''(x)\vert \leq M$ for all x on the domain of integration.

There are many ways to do this problem. In class, we discussed either using Taylor series or polynomial test functions.

• Taylor series.

  We expand the general integral in a Taylor series so that we acquire a trapezoid rule plus error terms. There are several ways to set up this integral, and I present one of them here.
  

  $$\int_{x_i}^{x_{i+1}} f(x) dx = \int_{x_i}^{x_{i+1}} \left[ \frac{1}{2}f(x)+ \frac{1}{2}f(x) \right] dx$$
  
  
  The f(x) in the first term on the left will be expanded about x_iand the second term will be expanded about x_i+1.
  $$\begin{align*}\frac{1}{2}f(x) & = \frac{1}{2} \left[ f(x_i) + f'(x_i)(x-x_i) + \... ...})(x-x_{i+1})^2 + \frac{1}{3!} f'''(\xi_1)(x-x_{i+1})^3 \right] \\ \end{align*}$$
  Now, here are some integrals we will need.

  
  $$\begin{align*}\int_{x_i}^{x_{i+1}} (x-x_i) dx, & = \frac{1}{2} h^2, & \int_{x_i... ...h^4, & \int_{x_i}^{x_{i+1}} (x-x_{i+1})^3 dx, & = -\frac{1}{4} h^4. \end{align*}$$
  

  Thus,
  $$\begin{multline}\int_{x_i}^{x_{i+1}} \frac{1}{2}f(x)+ \frac{1}{2}f(x) dx = \fra... ...{x_{i+1}} [f'''(\xi_1)(x-x_i)^3 + f'''(\xi_2)(x-x_{i+1})^3]dx. \end{multline}$$
  

  On the first line, we see that we have successfully acquired the trapezoid rule. The remaining terms describe the error in this quadrature.

  Now, we need to expand in a Taylor series again.
  $$\begin{align}f'(x_{i+1}) & = f'(x_i) + h f''(x_i) + \frac{1}{2} h^2 f'''(\xi_3), \\ f''(x_{i+1}) & = f''(x_i) + h f'''(\xi_4). \end{align}$$
  

  Substituting (2-3) into (1), we have
  $$\begin{multline*}\int_{x_i}^{x_{i+1}} \frac{1}{2}f(x)+ \frac{1}{2}f(x) dx = \fr... ...x_{i+1}} [f'''(\xi_1)(x-x_i)^3 + f'''(\xi_2)(x-x_{i+1})^3]dx. \end{multline*}$$
  And so,
  $$\begin{multline*}\int_{x_i}^{x_{i+1}} \frac{1}{2}f(x)+ \frac{1}{2}f(x) dx = \fr... ...x_{i+1}} [f'''(\xi_1)(x-x_i)^3 + f'''(\xi_2)(x-x_{i+1})^3]dx. \end{multline*}$$
  We can use the triangle inequality to bound the remaining integrals, and we see that all the extra terms are at a higher order O(h⁴).

  Thus,
  

  $$\int_{x_i}^{x_{i+1}} f(x) dx = \frac{1}{2} h \left[f(x_i) + f(x_{i+1}) \right] + -\frac{1}{12} h^3 f''(x_i) + O(h^4).$$
  
  
  Therefore,
  
  $$\left\vert\int_{x_i}^{x_{i+1}} f(x) dx - \frac{1}{2} h \left... ...right\vert \leq \frac{1}{12} h^3 \vert f''(x_i)\vert + O(h^4),$$
  
  
  establishing the result that we need.

• Polynomial test function.

  We anticipate a second order method, so we set up a second order polynomial test function.
  

  $$p(x) = c_0 + c_1 x + \frac{c_2}{2} x^2.$$
  
  

  We seek the difference between the exact integral,
  

  $$I = \int_0^L p(x) dx = c_0 L + c_1 \frac{L^2}{2} + c_2 \frac{L^3}{6},$$
  
  
  and the trapezoid quadrature,
  
  $$T_N = \frac{h}{2} \left[p(0) + p(L)\right] + h \sum_{i=2}^{N} p(x_i).$$
  
  

  Changing the indices on the sum, and substituting $h = \frac{L}{N}$and x_i= h(i-1),
  

  $$T_N = \frac{L}{2N} \left[p(0) + p(L)\right] + \frac{L}{N} \sum_{i=1}^{N-1} p\left(\frac{iL}{N}\right).$$
  
  

  Now, we evaluate the polynomial and see what happens.
  $$\begin{align*}T_N & = \frac{L}{2N} \left[c_0 + c_0 + c_1 L + \frac{1}{2} c_2 L^... ...ac{L^2}{2} + c_2 \left[ \frac{L^3}{6} + \frac{L^3}{12 N^2} \right] \end{align*}$$
  

  Thus, we see that
  

  $$I - T_N = -\frac{c_2 L^3}{12 N^2}$$
  
  
  which is the desired result because we interpret L as b-a and c₂ as M, the upper bound on the second derivative. Cool, eh?

5.

[related to KC 6.1:5] Derive a procedure for determining a polynomial of degree at most 2 for which p(0), p(1) and $p'(\xi)$is prescribed, $\xi$ being any preassigned point.

For simplicity, I introduce the following notation:

$$p(0) = p_0, \ \ \ \ p(1) = p_1, \ \ \ \ p'(\xi) = p_\xi'.$$

Since we have three pieces of information, we anticipate that we are going to be able to fit this with a second order polynomial. (In retrospect, we know that there is something fishy about $p_\xi'$. Remember, we lose some information when we differentiate.) Anyway, here we go!

p(x) = c₀ + c₁ x + c₂ x².

$$\begin{align*}c_0 & = p_0, \\ c_0 + c_1 + c_2 & = p_1, \\ c_1 + 2 c_2 \xi & = p_\xi'. \end{align*}$$
This can be solved by brute force and substitution, or we can examine the corresponding linear system.

$$\left[ \begin{array}{ccc} 1 & 0 & 0 \\ 1 & 1 & 1 \\ 0 & 1 &... ...ft[ \begin{array}{c} p_0 \\ p_1 \\ p_\xi' \end{array}\right]$$

If we study the matrix, we see that its determinant is $2\xi - 1$. Thus, the system is not invertible if $\xi = \frac{1}{2}$. If $\xi \neq \frac{1}{2}$, the coefficients are uniquely determined from the data about p.

In the degenerate case, $\xi = \frac{1}{2}$, we must conclude that we must assume that some of the coefficients of the second order polynomial are linearly dependent, and we can introduce an extra order to the polynomial. Thus,

p(x) = c₀ + c₁ x + c₂ x² + c₃ x³.

And,
$$\begin{align*}c_0 & = p_0, \\ c_0 + c_1 + c_2 + c_3& = p_1, \\ c_1 + c_2 + \frac{3}{4} c_3 & = p_\xi'. \end{align*}$$
By inspection, we see that c₁ and c₂ only appear together as a sum. Thus, we can define a new variable, d = c₁ + c₂, and we must accept that we have no specific information about c₁ and c₂. If we know either c₁ or c₂, the other is then determined by d. The new linear system is

$$\left[ \begin{array}{ccc} 1 & 0 & 0 \\ 1 & 1 & 0 \\ 0 & 1 &... ...t[ \begin{array}{c} p_0 \\ p_1 \\ p_\xi' \end{array}\right].$$

The determinant of the corresponding matrix is $\frac{3}{4}$, so we know that c₀, d and c₃ can always be determined uniquely. Again, there is some freedom is choosing the pair c₁ and c₂ as long as the sum is d.

6.

[KC 7.2: 1] Derive the Newton-Cotes formula for

$$\int_0^1 f(x) dx$$

based on the nodes $0, \frac{1}{3}, \frac{2}{3}, 1$.

Here, since we must integrate the polynomial explicitly, I would recommend using a divided difference table as shown in the class notes.

0	f(0)
		$3 \left[f\left(\frac{1}{3}\right)-f(0)\right]$
$\frac{1}{3}$	$f\left(\frac{1}{3}\right)$		$\frac{9}{2} \left[f\left(\frac{2}{3}\right)-2 f\left(\frac{1}{3}\right)+ f(0)\right]$
		$3 \left[f\left(\frac{2}{3}\right)-f\left(\frac{1}{3}\right)\right]$		$\frac{9}{2}\left[f(1)-3 f\left(\frac{2}{3}\right)+ 3 f\left(\frac{1}{3}\right)-f(0)\right]$
$\frac{2}{3}$	$f\left(\frac{2}{3}\right)$		$\frac{9}{2}\left[f(1)- 2 f\left(\frac{2}{3}\right)+f\left(\frac{1}{3}\right)\right]$
		$3 \left[f(1)-f\left(\frac{2}{3}\right)\right]$
1	f(1)

Thus, the correct interpolating polynomial is
$$\begin{multline*}p(x) = f(0) + 3 \left[f\left(\frac{1}{3}\right)-f(0)\right] x +... ...right] x \left(x-\frac{1}{3}\right) \left(x-\frac{2}{3}\right). \end{multline*}$$

To find the quadrature coefficients, we integrate over the panel.
$$\begin{multline*}\int_0^1 p(x) dx = f(0) + \frac{3}{2} \left[f\left(\frac{1}{3}\... ...eft(\frac{2}{3}\right)+ 3 f\left(\frac{1}{3}\right)-f(0)\right]. \end{multline*}$$
Now, we collect terms for each evaluation of f.

$$\int_0^1 p(x) dx = \frac{1}{8} f(0) + \frac{3}{8} f\left(\fra... ...t) + \frac{3}{8} f\left(\frac{2}{3}\right) + \frac{1}{8} f(1).$$

This is ``Simpson's $\frac{3}{8}$ Rule'' which is discussed a little bit in class. See? Building your own quadrature is not so bad, eh?