Monday, 07/28/03

Minimization Using Derivatives

Assume we want to minimize f(x) and it has a unique minimum at p, a ≤ P ≤ b and we start the search at P.

If f ‘(P₀) > 0 then P is at the left of P₀.

P₀

If f ‘(P₀) < 0 then P is at the right of P₀.

Image2

P₀

If the derivative of f(x) is available, then we can use any of the methods we have covered to solve non-linear equation to solve f ‘(x) = 0. (You can use regular-falsi, bisection, Newton Raphson, secant, etc).

Example:

f(x) = sin(x) f ‘(x) = cos(x)

Find minimum in the interval [-¾π,0].

Use bisection for f ‘(x) = 0

A	f ‘(a)	b	f ‘(b)	c	f ‘(c)
-2.356	-0.7071	0	1	-1/781	.3827
-2.356	-.7071	-1.1781	.3827	-1.7672	-.1951
-1.7672	-1.951	-1.1781	.3827	-1.4727	.0980
…	…	…	…	…	…
				1.5708	0

Finding Minimum for Multiple Variables

Steepest Descent of Gradient Method.

Assume we want to minimize f(x) of N variables where X = (x₁, x₂, … , x_N).

The gradient "f(x) is a vector function defined as:

grad (f(x)) = "f(x) = ( df(X)/dx₁, df(X)/dx₂, …, df(X)/dx_N)

From the concept of gradient we know that the gradient vector points to the direction of the greatest of f(X).

Since we are looking for the minimum, we will use -"f(x), that points in the direction of greatest decrease.

Basic gradient method:

- Start at point P₀ and move along the line in the direction -"f(x) on small increment.

- That is going to take you closer to the minimum

P₁ = P₀ - "f(P₀)h h is a small increment

…

P_k = P_k-1 - "f(P_k-1)h

Example

f(x) = x² + 1

f ‘(x) = 2x

"f(x) = 2x_i x_k = x_k-1 – 2x_k-2 h

Let h = 0.1 "f(x) = 2(-1)_i

x₀ = -1

x₁ = -1 - 2(-1)0.1 = -1 + .2 = -0.8

x₂ = -.8 – 2(-.8)(.1) = -.8 + .16 = -.64

x₃ = -.64 – 2(-.64)(.1) = -.64 + .128 = -.512

…

0

The gradient method does not guarantee a global minimum. The minimum obtain is the closest in the path from the starting point.

A method used to find a global minimum is called “simulated annealing” that “shakes” the curve from a high intensity to a lower intensity until the method arrives to the global minimum. “Simulate annealing” simulated cooking down of a metal or crystalline structure.

Tuesday, 07/29/03

Numerical Solution of Differential Equations

Introduction

Example

y’ = ky

dy/dt = ky → òdy/y = òk dt

ln y = kt + k₂

e^{ln y} = e^k₂ e^kt

y = k₃ e^kt

Very often differential equations do not have an analytic solution, so they have to be approximated using numerical method.

Euler’s Method

Let [a,b] be the interval which we want to find the solution of y’ = f(t,y), which y(a) = y_0. We will find a set of points:

(t_0, y₀), (t₁, y₁), …, (t_k,y_k) that are used to approximate y(t) ≈ y(t,k)

Also t_k = a + kh for k = 0,1,2,…M

We want to solve y’ = f(t,y) over [t₀, t₁, …, t_m] with y(t₀) = 0.

Using Taylor’s expansion to approximate y(t):

y(t) = y(t₀) + y’(t₀)(t-t₀) + y”(c₁)(t-t₀)²/2

y”(c₁)(t-t₀)²/2 error, a < c₁ < b

We use this equation to obtain y(t₁):

y(t) = y(t₀) +y’(t₀)(t₁-t₀) + y”(c₁)(t-t₀)²/2

y₁ y₀ y₀’

If the step size is small enough, we can neglect the second order error.

y₁ = y₀ + hy’(t₀)

y₁ = y₀ + hf(t₀,y₀)

In general:

y_k+1 = y_k + h

y_k+1 = y_k + hf(t_k,y_k) for k=0,1,…M-1

Example

y’ = t² – y y(0) = 1 h = .2

f(t,y), so f(t,y) = y’

t₁ = 0 + 0.2 = 0.2

y₁ = y₀ + hf(t₀,y₀) = 1 + .2(0²-1) = .8

t₂ = .2 + .2 = .4

y₂ = .8 + .2(.2² - .8) = .648

t₃ = .4 + .2 = .6

y₃ = .648 + .2(.4² - .648) = .5504

t₄ = .6 + .2 = .8

y₄ = .5504 + .2(.6² - .5504) = .5123

Analytical solution is: y(t) = -e^-t + t² – 2t + 2

y(.4) = -e^-.4 + (.4)² – 2(.4) + 2

= .6897 (Euler: .648)

y(.8) = -e^-.8 + (.8)² – 2(.8) + 2

= .5407 (Euler: .5123)

Heun’s Method

We want to solve: y’(t) = f(t,y(t)) over [a,b] with y(t₀) = y₀.

We can integrate y’(t) over [t₀,t₁]

ò_to^t1y’(t) dt = ò_to^t1f(t,y(t)) dt

y(t₁) – y(t₀) = ò_to^t1f(t,y(t)) dt

We use numerical integration to approximate the integral in the right side. Use trapezoidal rule

y(t₁) – y(t₀) = h/2 (f(t₀,y(t₀)) + f(t₁,y(t₁)))

then we have

y(t₁) = y(t₀) + h/2 (f(t₀,y(t₀)) + f(t₁,y(t₁)))

Observe that we still need to know f(t₁,y(t₁)), that involves y(t₁) that is what we want to solve. To eliminate this circular reference, we use Euler’s approximation to approximate y(t₁).

y(t₁) = y(t₀) + hf(t₀,y(t₀))

So we have:

y(t₁) = y(t₀) + + h/2 (f(t₀,y(t₀)) + f(t₁,y(t₀)+hf(t₀,y(t₀)))

y₁ = y(t₁) y₀ = y(t₀)

y₁ = y₀ + + h/2 (f(t₀,y(t₀)) + f(t₁,y₀+hf(t₀,y₀)))

P₀

In general:

P_k+1 = y_k + hf(t_k,y_k)

Y_k+1 = y_k + h/2 (f(t_k,y_k) + f(t_y,P_k+1))

Euler approximation is used as a prediction and the integral approximation is used as correction.

Example

f(y,t) = y’ = t² – y y(0) = 1 h = .2

k = 0 y₀ = 1 t₀ = 0

k = 1 t₁ = 0 + .2 = .2

P₁ = y₀ + hf(t₀,y₀)

= y₀ + h(t₀² – y₀)

= 1 + .2(0² – 1) = .8

y₁ = y₀ + h/2 ((t₀² – y₀) + (t₁² + P₁))

= 1 + .2/2 ((0² – 1) + ((.2)² + .8))

= .8240

k = 2 t₂ = .4

P₂ = .8240 + .2(.2² – 0.8240) = 0.6672

y₂ = .8240 + .2/2 ((.2² – 0.8240) + .4² - .6672)

= .6944

k = 3 t₃ = .6

P₃ = .6949 + .2(.4² - .6949) = .5879

y₃ = .6949 + .2/2 ((.4² - .6949) + (.6² - .5879))

= .6186

k = 4 t₄ = .8

P₄ = .6186 + .2(.6² - .6186) = .3669

y₄ = .6186 + .2/2 ((.6² - .6186) + (.8² - .5669))

= .6001

k	tk	yk(Euler)	yk(Heun)	yk(exact)
2	.4	.648	.6949	.6897
4	.8	.5123	.6001	.5907

Wednesday, 07/30/03

Numerical Solution of Differential Equations

Taylor Series Method

Using the Taylor expansion to approximate the solution:

y(t₀+h) = y(t₀) + hy’(t₀) + h²y”(t₀)/2 + … + h^Ny^N(t₀)/N! + O(h^N+1)

We want to solve

y’ = f(y,t) y(t₀) = a

From the Taylor expansion we can find the successive points.

y_k+1 = y_k + hy_k’ + h²y_k”/2 + h³y”’/3! … + h^Ny^N/N!

We still need to compute y_k”, y_k”’ … etc using y’ = f(y,t).

The error in the approximation will be O(h^N+1). The higher the order of the Taylor series, the smaller the error.

Example

Solve: y’ = t² – y y(0) = 1 h = .2

Use N = 3

y_k+1 = y_k + hy_k’ + h²y”(t₀)/2 + … + h³y^’’’(t₀)/6

y’ = t² – y

y” = 2t – y’

y”’ = 2 – y”

k=0 t=0 y₀ = 1

y₀’ = 0² – 1 = -1

y₀” = 2(0) – (-1) = 1

y₀”’ = 2 - 1 = 1

k=1 t=.2 y₁ = 1 + .2(.1) + (.2)²(1)/2 + (.2)³(1)/6 = .821333

y₁’ = (.2)² – .821333 = -.781333

y₁” = 2(.2) – (-.781333) = 1.181333

y₁”’ = 2 – 1.181333 = .818667

k=2 t=.4 y₂ = .821333 + .2(-.781333) + (.2)²(1.181333)/2 + (.2)³(.818667)/6 = .689785

y₂’ = (.4)² – .689785 = -.524785

y₂” = 2(.4) – (-.524785) = 1.329785

y₂”’ = 2 – 1.329785 = .670215

k=3 t=.6 y₃ = .689785 + .2(-.524785) + (.2)²(1.329785)/2 + (.2)³(.670215)/6 = .611317

y₃’ = (.6)² – .611317 = -.251317

y₃” = 2(.6) – (-.251317) = 1.451317

y₃”’ = 2 – 1.451317 = .548633

k=4 t=.8 y₄ = .611317 + .2(-.251317) + (.2)²(1.451317)/2 + (.2)³(.5486333)/6=.590812

y(.8) =

exact: .5907

Euler: .5123

Heun: .6001

To reduce approximation error, increase the order N of the expansion or reduce the value of h. I t is more effective to increase N.

Error = O(h^N+1)

Runge-Kutta method of order 4.

- Taylor method gives a good approximation, but it requires the derivative of the function.

- The Runge-Kutta method of order 4 simulates the accuracy of the Taylor series method using an order N=4 but it does not require the computation of derivatives.

y_k+1 = y_k + h(f₁ + f₂ + f₃ + f₄)/6

where:

f₁ = f(t_k,y_k)

f₂ = f(t_k + h/2, y_k + h/2 f₁)

f₃ = f(t_k + h/2, y_k + h/2 f₂)

f₄ = f(t_k + h, y_k + h f₃)

Example

y’ = f(t,y) = t²-y y(0) = 1 h = .2

t₀=0 y₀ = 1

k=0 f₁ = 0² -1 = -1

f₂ = f(0 + .2/2, 1+.2/2(-1))

= f(.1,.9) = .1² - .9 = -.89

f₃ = f(0 + .2/2, 1+.2/2(-.89))

= f(.1,.911) = .1² - .911 = -.901

f₄ = f(0 + .2, 1+.2(-.901))

= f(.2,.8198) = .2² - .8198 =- .7798

t₁=.2 y₁ = (1 + .2( -1 + 2(-.89) + 2 (-.901) + (-.7798))/6

k=1 = .821273

f₁ = .2² -.821273 = -.781273

f₂ = f(.2 + .2/2, 0.821273+.2/2(-.781273))

= f(.3,.743186) = .3² - .743146 = -.653146

f₃ = f(.2 + .2/2, .821273+.2/2(-.653146))

= f(.3,.755958) = .3² - .755958 = -.665958

f₄ = f(.2 + .2, .821273+.2(-.665958))

= f(.4,.688081) = .4² - .688081 = -.528081

t₁=.4 y₂ = (.821273 + .2( -.781273 + 2(-.653146) + 2 (-.665958) + (-.528001))/6

k=1 = .689688

exact: y(.4) = .6847

Taylor: y(.4) = .689783

Thursday, 07/31/03

System of Differential Equations

Assume we have

dx/dt = f(t,x,y) with x(t₀) = x₀

dy/dt = g(t,x,y) y(t₀) = y₀

The solution to this system are functions x(t) and y(t) that when derivated and substituted in the system of equations give equality.

Example of a system of differential equations

x’ = x + 2y x(0) = 6

y’ = 3x + 2y y(0) = 4

solution:

x(t) = 4e^4t + 2e^-t

y(t) = 6e^4t – 2e^-t

You can verify that this is the solution by derivating and substituting.

Numerical Solution

Euler’s Method

We can extend Euler’s method of a single differential equations to a system of equation.

Dx/dt = f(t,x,y)

Dy/dt = g(t,x,y)

Dx = f(t,x,y)dt

Dy = g(t,x,y)dt

Approximation dx_k = x_k+1 - x_k

dy_k = y_k+1 - y_k

dt_k = t_k+1 - t_k

x_k+1 – x_k = f(t_k,x_k,y_k) (t_k+1-t_k)

y_k+1 – y_k = g(t_k,x_k,y_k) (t_k+1-t_k)

So we have that

x_k+1 = x_k + f(t_k,x_k,y_k) h x(t₀) = x₀

y_k+1 = y_k + g(t_k,x_k,y_k) h y(t₀) = y₀

The Euler’s method does not give good accuracy because it approximates the Taylor expansion only to the first derivative.

Runge-kutta for systems of differential equations

Rk can be extended to systems of linear equations.

The formulas are: (Rk order 4)

x_k+1 = x_k + h/6 (f₁ + 2f₂ + 2f₃ + f₄)

y_k+1 = y_k + h/6 (g₁ + 2g₂ + 2g₃ + g₄)

f₁ = f(t_k, x_k, y_k)

f₂ = f(t_k + h/2, x_k + h/2 f₁, y_k + h/2 g₁)

f₃ = f(t_k + h/2, x_k + h/2 f₂, y_k + h/2 g₂)

f₄ = f(t_k + h, x_k + hf₃, y_k + hg₃)

g₁ = g(t_k, x_k, y_k)

g₂ = g(t_k + h/2, x_k + h/2 f₁, y_k + h/2 g₁)

g₃ = g(t_k + h/2, x_k + h/2 f₂, y_k + h/2 g₂)

g₄ = g(t_k + h, x_k + hf₃, y_k + hg₃)

Higher Order Differential Equations

Higher order differential equations involved higher derivatives x”(t), y”(t).

For example:

mx”(t) + cx’(t) + kx(t) = g(t)

The higher order differential equations can be transformed to a system of differential equations of first order.

Use the substitution:

y(t) = x’(t)

The system

mx”(t) + cx’(t) + kx(t) = g(t)

obtain x”(t)

x”(t) = (g(t) – cx’(t) – kx(t)) / m

substituting x”(t)

y’(t) = (g(t) – cy(t) – kx(t)) / m --------- 1

This is the first differential equation of first order in the system.

The second equation is

x’(t) = y(t) ------------- 2

Example

4x” + 3x’ + 5x = 2 with x(0) = 1

x” = (2 – 3x’ – 5x) / 4 x’(0) = 3

Let y = x’ and substituting in x”

y’ = (2 – 3y – 5x)/4 x(0) = 1, y(0) = 3

x’ = y

Friday, 08/01/03

CS 314 Final Review

Curve Fitting

- Least squares line

- Least squares for non-linear equations

- Transformation for data linearization

- Polynomial fitting

- Splines

+ The 5 properties(I … V)

+ Proof

+ How to use the formula to get the spline coefficient

+ End-point constraints

Numerical Differentiation

- Limit of difference quotient

- Central difference formula of order O(n²)

- Central difference formula of order O(n⁴)

Numerical Integration

- Trapezoidal rule

- Simpson’s rule

Numerical Optimization

- Local/global minimum/maximum point

- Properties of minimum/maximum point

- The golden ratio method

- Minimization using derivatives

- Steepest Descent of Gradient method

Solution of Differential Equations

- Euler’s method

- Heun’s method

- Taylor Series method

- Runge-Kutta method of order 4

- System of Differential equations

+ Euler’s method

+ Runge-kutta of order 4

- Higher order differential equation

Final

20% first half, 80% second half

Study:

Class notes, book, homework, do exercises in book

You may bring:

1 page formula

calculator

Final exam: Friday, Aug 8^th

10:20-12:20

wthr 104