Quadratic Polynomial

buy-out policy for this assignment Under certain conditions, you may buy-out one or two of the
questions that appear in this assignment. The buyout’s cost is 3 grace days per each question, and
questions cannot be prorated (i.e., you cannot buy 2/3 of a question for two grace days.) Check care-
fully your grace days balance before you make a decision here. you must have turned in all the
previous five assignments in order to be eligible for buyout. If you buyout a question, write
its number, and write next to it:   buy-out: three grace days
Overdraft clause: If you do not have sufficiently many grace days, your grade for that question will
be 0. Otherwise:
(a) If you buyout one question, we will grade the other two questions and will assign you a grade
for the assignment based solely on those other questions. Of course, if you do poorly on those
questions, the buy-out may end in a fiasco.
(b) If you buyout two questions we will grade the remaining one. The rule is still as above.

(1) (This problem deals with numerical integration .) The midpoint rule is constructed to be exact
for all constant polynomials. However, it is exact for all linear polynomials, too. Simpson’s rule
was designed to be exact for all quadratic polynomials . However, it is exact for all cubics, too.
The common property of these two rules is that, for both,

with w(t) := (t − t₁)(t − t₂)..., and with t₁, t₂, . . . the nodes of the scheme.
Rules that are even better can be constructed if one insists that w(t) satisfies additional
conditions such as

etc. The integral dt is called the “jth moment of the function w” (wrto the interval [a, b]).
Thus, the name of the game is to choose the nodes in a way that w has as many zero moments as
possible. If all the moments of w up to the jth are zero, your rule scores j + 1 extra points, i.e.,
while you use n points in your rule (whatever n is), the error looks like you are using n + j + 1
points. It is not hard to see that for a given n, the maximum possible value of j is j = n − 1. In
that case, you use n nodes, and the error looks like you are using 2n nodes. Such rules are called
“Gaußian rules”.

In this question you will construct a Gaußian rule for n = 2 and will observe the usefulness of
it. (For n = 1 the Gaußian rule is midpoint.) Your general goal is to find an integ ration scheme of
the form
(1)

(a) Your choice for t₁, t₂ should be as follows. First, find a quadratic polynomial w(t) = t²+bt+c
that satisfies the fol lowing two conditions :

(Hint: plug in the expression you have for w into these integrals, and carry out the integration.
You will get two equations with the two unknowns b , c.) Then find the two roots of w(t). These
are your t₁, t₂. (t₁ and t₂ are not rational numbers. If you wish to use an approximate decimal
representation of them, be sure to use high precision approximation. If you found a closed
form for those roots, e.g., t₁ = keep it in this form).

(b) Now, find w₁ and w₂. (There are several ways how to find them. One way is to recall the
general approach from class: once you are given t₁, t₂ and asked to derive a rule, you, at least
theoretically, interpolate the underlying function by a (linear) polynomial at t₁ and t₂ and
then integrate the polynomial. If you write that polynomial in the Lagrange form, the weights
are then the integrals of the Lagrange polynomials. Check the notes for more details.)

(c) Your scheme is based on two nodes, hence clearly must be exact for all linears. However, if
you did it right, you should be exact for all cubics (like a scheme that uses four points). Check
that. Is it also correct for all quartics (‘quartic’=degree 4)?

(d) Now apply your great scheme to f(x) = e^x, f(x) = 1/(x+1), and f(x) = sin x. Compare your
accuracy with composite midpoint with two subintervals, and with ( simple ) Simpson. (You
may use either matlab or a calculator here ).

(2) In this question, you compare three methods for solving IVP ’s: (a) The Milne-Simpson
predictor-corrector, (b) the Adams-Bashforth-Moulton method, and (c) the Runge-Kutta method
of order 4. You will use each of the methods in order to approximate y(b), with y(t) the solution
of an initial value problem
y′ = f(t, y),
y(a) = y₀.

Note that Runge-Kutta uses four evaluations of f per each step, while the predictor-corrector
methods use two evaluations per each step, hence once you choose some step-size h for Runge-
Kutta, the ‘fair comparison’ should be with predictor-corrector with step size h/2 (why?).

(a) Write, for each of the above methods, a (small) routine that, given the function f, the numbers
a, b, the initial value y0, and a parameter h, applies the corresponding method for finding
(approximately) y(b). Note that in the predictor-corrector methods you need first to find
y₁, y₂, y₃, before you can start the method, so use Runge-Kutta in order to find these three
required values. Use the given h as the step-size for Runge-Kutta, and use step-size h/2 for
each of the predictor-corrector methods. Test your routines by running them on some known
example (e.g., some example from a text book).

(b) You are now given the IVP

Apply each of your routines with h = 1/4, 1/8 . . . until the difference between two consecutive
approximations (i.e., the one for the current h and the one for the previous h) is smaller then
10^-6 (or until you de termine that the method does not work). Repeat now with two other IVP
of your choice (with respect to the same interval and the same initial value). Draw conclusions
(concerning the relative performance of the three methods here).

(c) Use your calculations in (b) is order to determine the order of each method. (After you
completed (b), you know more or less the correct value of y(2), so you can calculate the errors
committed when using h = 1/4, 1/8, . . .. If you divide the error for a certain h by the error for
the previous h you should be able to determine the order of your method).

(3) Your goal in this question is to practice solving IVP that involve several functions, as well as
IVP of higher order (i.e., that involve derivatives of order 2 or higher.) You will do it in two steps:
in the first step, you will observe that the methods that we discussed in class for solving a single
equation apply as well to the system case . In the second step you will realize that a second order
differential equation can be reduced to a system of two first order equations. You can either use
matlab or a calculator for the sake of the computations.

(a) Show that the functions

satisfy the initial value problem:

Then (as a warm up) apply one step of Euler’s method (with h = .1) in order to approximate
the solution at 1.1. Compare with the exact solution.

(b) Now, use the Runge-Kutta (RK4) method for the same purpose (i.e., approximating the solu-
tion at 1.1). First, do that in one step (h = .1), and then try two steps (h = .05). Find the
ratio between the two errors you obtained here. Any conclusions?

(c) After doing (a),(b) you have a basic idea how systems are solved. Suppose now that instead
of a system of two equations you have a single equation of order 2, e.g., we look for a function
y(t) such that
(1)
(Note that now we specify y′(a) as well in the IVP, not only y(a)). Define a vector x(t) of two
functions: x₁(t) := y(t), and x₂(t) := y′(t). Now, you need to convert the second order IVP
for y to a first order IVP on x (Hint: look at any text book, or wait until we cover it in class.)

(d) At this point, you can solve numerically the original IVP for y, by solving the equivalent
IVP on x. Apply one step of predictor-corrector with Euler the predictor and trapezoid the
corrector in order to find y(1.1) (i.e., h = .1.)

(e) Find the solution y on a dense mesh on the interval [1, 2] and plot its graph . (Note: you can
use any method, but you will need to know that your solution is reasonably accurate.)