Preface

---

This is a tutorial made solely for the purpose of education and it was designed for students taking Applied Math 0330. It is primarily for students who have very little experience or have never used Mathematica before and would like to learn more of the basics for this computer algebra system. As a friendly reminder, don't forget to clear variables in use and/or the kernel.

Finally, the commands in this tutorial are all written in bold black font, while Mathematica output is in normal font. This means that you can copy and paste all commands into Mathematica, change the parameters and run them. You, as the user, are free to use the scripts for your needs to learn the Mathematica program, and have the right to distribute this tutorial and refer to this tutorial as long as this tutorial is accredited appropriately.

The most important single result in numerical computations is Taylor's theorem, which we now state below.

Theorem: Let f(x) have n+1 continuous derivatives on the closed interval [a.b] for some positive n, and let \( x, x_0 \in [a,b] . \) Then

for and where ξ is a point between x and x₀. ■

Brook Taylor (1685--1731) was educated at St. John's College of Cambridge University (England), from which he graduated in 1709. He wrote two very significant books, ‘Methodus incrementorum directa et inversa’ and ‘Linear Perspective,’ which were published in 1715. He was the one to invent ‘Integration of Parts’ and also a series called the ‘Taylor’s Expansion.’ The Taylor’s Theorem is based on the letter written by him to Machin in which he tells about the origination of this idea. Although it appears that he did not appreciate its importance and he certainly did not bother with a proof of Taylor's theorem, which was discovered and published in 1694 by Johann Bernoulli (1667--1748). Taylor acknowledged that his work was based on that of Newton and Kepler. However, Johann Bernoulli claimed his priority in discovering both of his main results (Taylor's expansion and integration by parts).

Polynomial Approximations

---

The Taylor series method is of general applicability, and it is a standard to which we compare the accuracy of the various other numerical methods for solving an initial value problem for ordinary differential equations: In what follows, it will be assumed that the given slope function f(x,y) has as many derivatives in a neighborhood of the initial point as needed. For simplicity, we consider the set of uniformly distributed grid points starting from x₀: \( x_n = x_0 + n\,h , \quad n=0,1,2,\ldots . \)

Theorem: Suppose that a function y(t) has m+1 continuous derivatives on the interval [a,b] containing grid points \( \{ x_n \}_{n\ge 0} . \) Assume that y(t) has a Taylor series expansion of order m about fixed value \( x_n \in [a,b] : \)

where the increment function is and In particular,

The approximate numerical solution to the initial value problem \( y' = f(x,y ) , \quad y(x_0 ) = y_0 \) over the interval [a,b] is derived by using the Taylor series expansion on each subinterval \( [ x_n , x_{n+1}] . \) The general step for Taylor's method of order m is

The Taylor method of order m has the property that the final global error is of the order \( O\left( h^{m+1} \right) . \) It is possible theoretically to choose m as large as necessary to make this error as small as desired. However, in practice, we usually compute two sets of approximations using step sizes h and h/2 and compare the results.

 

II. Second Order Polynomial Approximation

---

Since the first order Taylor series approximation is identical with Euler’s method, we start with the second order one:

where the increment function Φ₂ is just adding the second order differential deviation to the next term in the same equation used for the first order Taylor expansion. We show how it works in the following examples.

Example: Consider the initial value problem \( y’=1/(2x-3y+5), \quad y(0)=1 \) and let us try to find out y(1).

To achieve this, we pick up a fixed step length h = 0.1, which means that we need to do 10 iterations to reach the end point x = 1.

At each step in the sequence, the expansion requires 3 addition, 1 subtraction, 3 multiplication, and 1 division operations. Since there are 8 steps to this sequence, in order to obtain the desired approximation, the code performs 80 operations, which would be hectic work to do by hand. The approximation is y₁₀ = 1.435290196. This algorithm requires at least 22 operations per step, which means that the entire sequence requires 220 steps. However, this is not significant since Mathematica does all the work.

To check the answer, we write a subroutine to evaluate numerically the initial value problem for arbitrary slope function, initial conditions, and step size:

Now we apply this subroutine to our slope function:

Using Mathematica command NDSolve, we check the obtained value:

Example: Let us take a look at the initial value problem \( y’=(x^2 -y^2) \sin (y), \quad y(0)=1 , \) for which it is impossible to find the actual solution. First we find two derivatives:

We check the formulas with Mathematica: Then we ask Mathematica to do all calculations We compare the second order Taylor series approximation y₁₀ = 3.10788 with those obtained with the aid of NDSolve (one needs to clear variables before invoking function NDSolve): So we see that Taylor's approximation gives two correct significant digits. Now we use another approach involving NestList command. Here is the direct application of polynomial approximation with many loops.

 

 

III. Third Order Polynomial Approximation

---

The third order Taylor approximation is adding a third order differential deviation to the equation for the 2nd order expansion.

where the increment function Φ₃ adds just one term to Φ₂. We will show how it works in two examples that were considered previously.

Example: We reconsider the IVP: \( y’=1/(2x-3y+5), \quad y(0)=1 . \) To use the third order polynomial approximation, we choose a uniform grid with constant step length h = 0.1.

Here is the syntax for the third order Taylor approximation (which takes too long to execute):

The approximation for this value is y₁₀=1.435283295 while the true value is 1.4352866048707.

The total number of numerical operations here is 420. It is obvious at this point why using mathematical programs such as Mathematica is a desired approach for such problems.

Example: We reconsider the IVP: \( y’=(x^2 -y^2) \sin (y), \quad y(0)=1 . \) To use the third order polynomial approximation, we choose a uniform grid with constant step length h = 0.1.

Another example of code (with the same output):

 

IV. Fourth Order Polynomial Approximation

---