## Taylor method with automatic differentiation

I teach numerical analysis course and in the one of the chapters, I teach my students about the numerical solution of ordinary differential equations (ODEs). One of numerical solution methods I teach is the so-called Taylor method. This method works by applying the Taylor series to the solution of the differential equation. Suppose we have one-dimensional differential equation as follow:

and .

Let’s assume that we only consider an autonomous ODE first. To find a solution at the next step, , we need to apply Taylor expansion, as follow:

.

Therefore, depending on how accurate we want our solution is, all we need is to compute the derivatives of . Unfortunately, computing all the derivatives is not an easy task for a computer. We have to derive all the derivatives by ourselves and then input them to the computer.

In this note, I will introduce an automatic differentiation as a tool to compute the derivatives of so that we can compute as accurately as possible. Let us consider the series expansion of as above, but I will use constant to replace the th derivative of ,

.

Let us derive the above series with respect to .

.

Notice that the derivative at the left hand side of the above equation is actually equal to since we know that satisfies the ordinary differential equation. Then we have,

.

Since is a series of then is also a series of . Then both sides of the above equation are series of , thus we can equate each of the coefficient of to get , and so on. The most difficult task in this problem is to find the series of since the function can be anything.

Let’s take an example. Suppose we want to integrate the following differential equation,

and .

Let us choose and compute . As before, we assume that is a series of as follows,

.

Differentiating both sides of the above equation, we get

.

The right hand side of the above is equal to

.

Therefore, we get,

, ,, and so on.

The real question is how to find , and so on. If the function only involves operations such as addition, subtraction and multiplication, it would be easy to find the ‘s. But if the function has operations such as division, exponential, logarithm and trigonometric functions then we have to think a bit harder on how to find ‘s. We go back to the previous example. Since the operations in this example are only multiplication and subtraction, we can easily compute , and so on. Using algebraic modification, we get the following,

, thus

, thus

, thus

and so on.

The more ‘s we compute, the more accurate our solution is. Using Taylor method this way, we can easily apply this method in a program. It should be clear why this method is called automatic differentiation because we don’t need to manually derive all the derivatives.

Next, I tried to numerically integrate the previous example. Below is the graphic of the error of Taylor method with automatic differentiation versus time. As we can see from the figure, that the error of our numerical integrator is so small and it is almost exact. The number of I used in this graph is 10.

I use python to numerically integrate this example. You can find the code to produce the above figure in the end of this note. One of the drawbacks of this method is when you change your ODE, you have to change your code as well. In other numerical integrators such as Runge Kutta, we don’t need to change the whole code, we just need to change the ODE and we should get the solution. Together with my student, I am trying to tackle this problem and I hope it will be posted in here.

import math as m import numpy as np import pylab as plt #preliminary constants order = 10 h = 0.1 t = 0.0 t_end = 5 y0 = 0.5 #defining solution sol = [] sol.append(y0) time = [] time.append(t) #start of the program while t<t_end: a = [] a.append(sol[-1]) for k in range(1,order): if k == 1: a.append(a[k-1]*(1-a[k-1])) else: sum = 0.0 for j in range(0,k-1): sum = sum + a[j]*(-a[k-1-j]) a.append((sum+a[k-1]*(1.0-a[0]))/k) sum = 0.0 for i in range(0,order): sum = sum + a[i]*h**i sol.append(sum) t = t + h time.append(t) #end of the program #real solution computation real_sol = [] for t in time: temp = m.exp(t)/(m.exp(t)+1.0) real_sol.append(temp) #plotting the error versus time sol = np.array(sol) real_sol = np.array(real_sol) time = np.array(time) error = sol-real_sol plt.plot(time[:],error[:]) plt.show()

## Mixture problem in Calculus

I keep forgetting the way to solve this problem. So, when I in future have no clue I know where to look.

**The problem**

A 200-galon tank is half full of distilled water. At time , a solution containing 0.5 pound/galon of salt enters the tank at the rate of 5 gal/minutes, and the well-stirred mixture is withdrawn at the rate of 3 gal/minutes.

a. Give a mathematical model representation of the above problem.

b. At what time, will the tank be full?

c. At the time the tank is full, how many pounds of salt will it contain?

**The mathematical model**

The model describing the above process is based on the following formula:

rate of change of the amount of salt = (rate of salt arrived) – (rate of salt departed)

Suppose is the volume of the liquid in the tank at time . The volume of liquid inside the tank will increase as the inflow (5 gal/minutes) is bigger than the outflow (3 galon/minutes). Thus, we have

100 galon + (5 gal/min – 3 gal/min) min,

gal,

as we know that initially the tank is half full (100 gal).

Suppose is the amount of salt inside the tank at time . The rate in of salt is as follows,

rate in (pound/min) (0.5 pound/gal) (5 gal/min) 2.5 pound/min,

while the rate out is

rate out (pound/min) outflow (gal/min)

(pound/min)

(pound/min).

Therefore, the mathematica model is given by,

. (1)

**The solution**

First we determine time at which the tank will be full. We use the following equation,

,

We solve for when the volume is 200 galon, then we obtain .

To solve the differential equation (1) we need the initial condition

as follow,

.

We are going to use Maple to solve the problem

` > model := diff(y(t),t) = 2.5 - (3*y(t))/(100+2*t); `

` > init := y(0)=0; `

` > dsolve({model,init}); `

We obtain the solution that is given by

When , the amount of salt inside the tank will be

.

## A geometric proof of trigonometric derivatives

As we already know, the first derivatives of and are and respectively. The question is how do we prove it?

Since high school, the only proof that I knew is to use two of the trigonometric identities which are:

- , and
- .

Hence,

.

The derivative of can also be derived using the second identity.

However, there is another proof, which is (I think) more elegant. It is using a geometric interpretation of and . Let us draw a unit circle in .

Thus we have and . We are looking for and, however we are not going to discuss the first derivative of as it can be done in the same way as .

We are now going to estimate , see the second figure. The following relationship applies:

. or .

Thus,

,

where is the chord PQ and is the arc PQ. As goes to zero, the ratio of the arc PQ and the chord PQ tends to 1 and the chord PQ tends to become a tangent line of the unit circle at , therefore . Assuming that the limit of multiplication is the same as the multiplication of limit, we have the following:

.

As mentioned before, the first derivative of can be done the same way.

Also, I should thank Pak Yudi Soeharyadi as I got the idea of this geometric proof from him.

## Fundamental Theorem of Calculus is fundamental !

Most of my students do not seem to understand the two concepts of integrals. I have been there as well, so I did not blame them for not understanding the concept of integrals. Even I sometimes forget. From that reason I try to post what the fundamental theorem of calculus means exactly.

It is all started with the definition first definition of integrals, which is formulated by the famous mathematician Bernhard Riemann. This is called the definite integral, or .

We can think that as the area below the curve , see the picture below.

Thus, is the area S. As a first approximation, we shall divide the interval into an equal number of segments, say . Then we can approximate the area of S as follows:

.

The right hand side of the equation above is called the Riemann sum. Our final step to find the value of definite integral is to take limit the Riemann sum as ,

.

Now, the second concept about integrals will be explained. This concept is called the indefinite integral or is also usually called the anti-derivative. Suppose are functions of such that (or we can say that the derivative of is or, the anti-derivative of is ), then

,

where is a constant real number.

Hence, these two are two different concepts that are not related one to another. However, there is a great theorem that relates those two. The theorem perfectly and comfortably links the limit of Riemann sum and the anti-derivative. Actually there are two great theorems.

**First Fundamental Theorem of Calculus** Suppose is a continuous function defined on closed interval [a,b] and define as

,

then is continuous and differentiable, moreover, .

**Second Fundamental Theorem of Calculus** Suppose is a function defined on closed interval [a,b] and is the anti-derivative of . If is integrable on [a,b], then

.

Therefore, if we want to compute, for examples, , we don’t have to find the limit of Riemann sum of as goes to zero.

## Result of mid semester exam 2 of Calculus II (UTS 2)

As I promised yesterday, you can now see the result of mid semester exam 2 of calculus II.

Please go to this page.

Next wednesday is gonna be our last tutorial. I hope you can do well in the final examination.

## Result of mid calculus II mid semester exam (UTS 1)

To all of my students in my tutorial class, subject Calculus IIa (T-19), Bu Mul and I finally finished the marking of your mid semester exams.

Klik here to go to my teaching page and go to Calculus IIa to find your exam results.

## Your says