Click on an image to go directly to a system or scroll to see all of the systems in this gallery.

OdeFactory Images and Annotations

An OdeFactory Slide Show

Click on a slide to zoom in.

Click "video" to see a video.

View/Sys/Gal: Ode " Iteration" in "IterationExs."
Range: (vMax,vMin) = (5.000,-5.000), (hMin,hMax) = (-5.000,5.000)
VFld: (0)

Simple Iteration

The purpose of this gallery is to show how the IMap view in OdeFactory can be used (a) to find the zeros of a function F(x) and (b) to solve the system of odes

        dx/dt = F(x,y),

        dy/dt = G(x,y).

(a) Newton's method for solving

        F(x) = 0

is the IMap view of

        x <- x-F(x)/(dF(x)/dx).

(b) Euler's method for solving the system of odes

V1:        dx/dt = F(x,y),

        dy/dt = G(x,y)

is the IMap view of

V2:        x <- x+h*F(x,y),

        y <- y+h*G(x,y)

where h is on the order of 0.01.

NOTE: Both vector fields can be viewed in either the Ode view or the IMap views. The Ode view of V1 is generated using the (default) RK4 solver. V2 is different from, but related to V1. For small h, the IMap view of V2 is the Ode view of V1.

An OdeFactory Slide Show

Click on a slide to zoom in.

Click "video" to see a video.

View/Sys/Gal: Ode " Iteration" in "IterationExs."
Range: (vMax,vMin) = (5.000,-5.000), (hMin,hMax) = (-5.000,5.000)
VFld: (0)

Simple Iteration

The purpose of this gallery is to show how the IMap view in OdeFactory can be used (a) to find the zeros of a function F(x) and (b) to solve the system of odes

        dx/dt = F(x,y),

        dy/dt = G(x,y).

(a) Newton's method for solving

        F(x) = 0

is the IMap view of

        x <- x-F(x)/(dF(x)/dx).

(b) Euler's method for solving the system of odes

V1:        dx/dt = F(x,y),

        dy/dt = G(x,y)

is the IMap view of

V2:        x <- x+h*F(x,y),

        y <- y+h*G(x,y)

where h is on the order of 0.01.

NOTE: Both vector fields can be viewed in either the Ode view or the IMap views. The Ode view of V1 is generated using the (default) RK4 solver. V2 is different from, but related to V1. For small h, the IMap view of V2 is the Ode view of V1.

Resized GIF graphic

Resized GIF graphic

View/Sys/Gal: Ode "( 1) solve x^2-x-1=0 w/Newton's method, in IMap View" in "IterationExs."
Range: (vMax,vMin) = (5.000,-5.000), (hMin,hMax) = (-2.250,7.750)
VFld: (x-(x^2-x-1)/(2*x-1))

Solve

        f(x) = x^2-x-1 = 0

using Newton's method.

Enter the iteration

        x <- x-(x^2-x-1)/(2*x-1).

Since f(x) is a quadratic function in x, we know the solutions are just

        x = (1+-sqrt(5))/2

         ～ -.618034, 1.618034

so we can check OdeFactory's results.

To see how Newton's method works, iterations have been started at seeds

        (0,-.2), (0,-.5) and (0,-.9).

If you click on the starting points, you will see that they all converge to -.618034.

Iterations started at seeds

        (0,1.5), (0,2) and (0,2.5)

all converge to 1.618034.

Image 1: IMap view.

Image 2: Ode view.

NOTE: the ode is not defined at x = .5 so, in the Ode view, the solver, RK4, fails (gives squiggles) when solutions are extended in -t.

Resized GIF graphic

View/Sys/Gal: Ode "( 1) solve x^2-x-1=0 with dx/dt = x^2-x-1" in "IterationExs."
Range: (vMax,vMin) = (5.000,-5.000), (hMin,hMax) = (-5.000,5.000)
VFld: (x^2-x-1)

This ode is defined by the equation:

        dx/dt = x^2-x-1 = f(x)

in the (t,x) coordinate system.

The red curve is f(x) rotated ccw 90 degrees. The asymptotes of the solution curve are the zeros of f(x). Clicking on the arrowhead gives

p(98.7576687116700600) = (-0.6180339887498923)

p(-1.2423312883435582) = (0.9886547811993517)

p(-101.2423312883584400) = (1.6180339887498900)

so f(x) = 0 at:

        x = -0.6180339887498923 and

        x = 1.6180339887498900.

The iteration being used is just RK4.

Image 1: Ode view with S on.

Resized GIF graphic

View/Sys/Gal: IMap "( 1b) plot y(x) = sin(x) w/Euler, in IMap1000 view" in "IterationExs."
Range: (vMax,vMin) = (5.500,-4.500), (hMin,hMax) = (-3.000,7.000)
VFld: (x+h,y+h*cos(x)), h=.1

Plot y(x) = sin(x) using Eyler's method.

First, compute dy/dx = cos(x) then write the 1D ode as the 2D system in dx/dt and dy/dt as:

        dx/dt = 1,

        dy/dt = cos(x).

Next use the Euler iterative method to solve the 2D system:

        x <- x+h,

        y <- y+h*cos(x).

Parameters are:

        h=.1

Use ICs: (x,y) = (0,0).

You need to get into the IMap view to see the plot of:

        y(x) = sin(x).

Image 1: IMap view, the plot of y(x) = sin(x) for x >= 0. IMaps are not invertable so they only give images for t >= 0.

Of course if you view the system:

        dx/dt = 1,

        dy/dt = cos(x)

in the Ode view, using the default RK4 method, you get a better plot.

Resized GIF graphic

View/Sys/Gal: Ode "( 1c) dx/dt = 1, dy/dt = cos(x) w/RK4 in Ode view" in "IterationExs."
Range: (vMax,vMin) = (5.500,-4.500), (hMin,hMax) = (-3.000,7.000)
VFld: (1,cos(x))

This system of odes is defined by the equations:

        dx/dt = 1,

        dy/dt = cos(x).

The solution to dx/dt = 1 is x = t and the solution to dy/dt = cos(t) is y = sin(t)+c so

        y(x) = sin(x)+c.

ICs are (x(0),y(0)) = (0,0) so c = 0.

Image 1: Ode view, trajectory fot ICs (0,0) is the curve y(x) = sin(x).

Resized GIF graphic

View/Sys/Gal: Ode "( 2) ode dx/dt = sin(x) w/RK4" in "IterationExs."
Range: (vMax,vMin) = (5.000,-5.000), (hMin,hMax) = (-4.000,6.000)
VFld: (sin(t))

This ode is defined by the equation:

        dx/dt = sin(t).

The solution is

        x(t) = -cos(t)+c.

Using ICs (t,x) = (0,0) gives c = 1 so the solution curve is the cos(t) function flipped over about the t axis and then slid up by 1.

The default algorithm used to generate the solution curve is RK4.

Resized GIF graphic

View/Sys/Gal: IMap "( 3) ode dx/dt = sin(x) w/Euler, in IMap1000 view" in "IterationExs."
Range: (vMax,vMin) = (5.000,-5.000), (hMin,hMax) = (-4.000,6.000)
VFld: (x+h,y+h*sin(x)), h=.01

This is Euler's method applied to

        dy/dx = sin(x)

which has been rewritten as

        dx/dt = 1,

        dy/dt = sin(x).

The iteration is:

        x <- x+h,

        y <- y+h*sin(x)

Parameters are:

        h=.01.

The seed is (0,0).

Try increasing h in the IMap view.

Resized GIF graphic

View/Sys/Gal: EMap "( 4) the EMapCT3 view of ( 3) with h = .13" in "IterationExs."
Range: (vMax,vMin) = (10.000,-10.000), (hMin,hMax) = (-9.436,10.564)
VFld: (x+h,y+h*sin(x)), h = .13;

This iteration is defined by:

        x <- x+h,

        y <- y+h*sin(x).

Parameters are:

        h = .13;

The EMap view gives an interesting image.

Click "Adj Ctrl Params ..." and check out the various color tables. Also try adjusting h.

Resized GIF graphic

View/Sys/Gal: Ode "( 5) ode (y,-x) w/RK4" in "IterationExs."
Range: (vMax,vMin) = (5.000,-5.000), (hMin,hMax) = (-5.000,5.000)
VFld: (y,-x)

This system of odes is defined by the equations:

        dx/dt = y,

        dy/dt = -x.

The trajectories are all circles centered at the origin. With ICs (x(0),y(0) = (0,1), the trajectory is the unit circle:

p(t) = (x(t),y(t)) = (sin(t),cos(t)).

RK4 uses a default iteration step size of h = 0.01 but h can vary. When the trajectory bends rapidly, h is decreased and where it bends slowly h is increased.

The RK4 computed value of

        p(100) = (-0.506366,0.862319)

and the true value is

        p(100) = (-0.506365,0.862319)

It takes 10,000 iterations to get to t = 100.

If you center at (0,1) and zoom in several times, you can see the error in the RK4 method.

Resized GIF graphic

View/Sys/Gal: IMap "( 6) ode (y,-x) w/Euler, h=.001, in IMap view" in "IterationExs."
Range: (vMax,vMin) = (5.000,-5.000), (hMin,hMax) = (-5.000,5.000)
VFld: (x+h*y,y+h*(-x)), h=.001

In this example an approximate solution to the system of odes

        dx/dt = y, dy/dt = -x

is computed using the Euler iteration method:

        x <- x+ Δt*dx/dt,

        y <- y+ Δt*dy/dt.

To see how it works, let Δt = h be the Euler iteration time step and use the IMap view to see the approximate solution to the system of odes defined by:

        dx/dt = x+h*y,

        dy/dt = y+h*(-x).

With ICs

        (x(0),y(0)) = (0,1), and h = .001,

the Euler method gives

        p(10000) = (-.546745,-.843279).

10000 is the number of iterations needed to get to the point p(10000).

The actual time for the 10000 steps is

t = (10000steps)*(.001time/step) = 10.

As in example ( 5), solving

        dx/dt = y, dy/dt = -x with

        ICs (x(0),y(0)) = (0,1) gives

        (x(t),y(t)) = (sin(t),cos(t)).

The trajectory is the unit circle and the exact solution is:

(x(10),y(10)) = (-.5440221,-839071).

As h increases, the computation is faster but the approximate Euler solution gets worse.

Resized GIF graphic

View/Sys/Gal: EMap "( 6b) ode (y,-x) w/Euler, h=.21, in EMapCT3 view" in "IterationExs."
Range: (vMax,vMin) = (5.000,-5.000), (hMin,hMax) = (-5.000,5.000)
VFld: (x+h*y,y+h*(-x)), h = .210;

This iteration is defined by:

        x <- x+h*y,

        y <- y+h*(-x).

Parameters are:

        h = .210;

EMap CT: 3

Resized GIF graphic

View/Sys/Gal: IMap "( 7) ode (y,-x) w/Euler, h=.01, IMap view" in "IterationExs."
Range: (vMax,vMin) = (5.000,-5.000), (hMin,hMax) = (-5.000,5.000)
VFld: (x+h*y,y-h*x), h=.01

The IMap view of this system shows Euler's method applied to

        dx/dt = y,

        dy/dt = -x.

The system of odes is defined by the equations:

        dx/dt = x+h*y,

        dy/dt = y-h*x

Here

        h = .01.

ICs are (0,1) and the trajectory should be the unit circle. It is easy to see the error in the Euler method with this h value.

Also see the EMap view for h = .25.

Resized GIF graphic

View/Sys/Gal: IMap "( 8) ode (y,-x) w/Euler, h=.1, in IMap1000 view" in "IterationExs."
Range: (vMax,vMin) = (5.000,-5.000), (hMin,hMax) = (-5.000,5.000)
VFld: (x+h*F,y+h*G), h=.1, F=y; G=-x

This system of odes is defined by the equations:

        dx/dt = x+h*F,

        dy/dt = y+h*G.

Parameters are:

        h=.1

Functions are:

        F=y; G=-x

ICs are: (x,y) = (0,1).

Click "Flow" in the IMap view.

In the EMap view, with h = .34, try the various color tables.

Resized GIF graphic

Resized GIF graphic

View/Sys/Gal: EMap "( 9) ode (y,x^2-y) w/RK4" in "IterationExs."
Range: (vMax,vMin) = (5.000,-5.000), (hMin,hMax) = (-5.000,5.000)
VFld: (y,x^2-y)

The Ode view shows a trajectory starting at (x,y) = (0,1) generated using the RK4 algorithm for the system:

        dx/dt = y,

        dy/dt = x^2-y.

The IMap view shows the orbit starting at seed (x,y) = (0,1) generated using the simple iteration algorithm:

        x <- y,

        y <- x^2-y.

Note: this is not the solution curve for the ode generated using the Euler method - see example (10).

For this vector field, the IMap view is not very interesting but the EMap view is somewhat interesting.

If you click the border button, "Bdr," then click the arrowhead, you will find that the solution curve for the ode crosses the upper border at (3.909638,5).

Resized GIF graphic

View/Sys/Gal: IMap "(10) ode (y,x^2-y) w/Euler, h=.1, in IMap1000 view" in "IterationExs."
Range: (vMax,vMin) = (5.000,-5.000), (hMin,hMax) = (-5.000,5.000)
VFld: (x+h*y,y+h*(x^2-y)), h=.1

In general, applying Euler's method to approximate a solution to the ode

        dx/dt = F(x,y),

        dy/dt = G(x,y)

involves the iteration

        x <- x+h*F(x,y),

        y <- y+h*G*x,y)

using a small value of h = Δt.

If F = y and G = x^2-y consider the system of odes defined by the equations:

        dx/dt = x+h*y,

        dy/dt = y+h*(x^2-y)

with

        h=0.1.

The IMap view of this system shows the Euler method applied to the ode system

        dx/dt = y,

        dy/dt = x^2-y

If h is 0.01 we get reasonable agreement with RK4 for the +t part of the trajectory. To get the -t part of the trajectory, just change the sign on h = -0.1.

In the EMap view, the system gives a fractal. Try CT9 with h = -.1.

Resized GIF graphic

View/Sys/Gal: EMap "(10b) ode (y,x^2-y) w/Euler, h=-.04, in EMapCT3 view" in "IterationExs."
Range: (vMax,vMin) = (7.066,-2.934), (hMin,hMax) = (-4.759,5.241)
VFld: (x+h*y,y+h*(x^2-y)), h = -.040;

This iteration is defined by:

        x <- x+h*y,

        y <- y+h*(x^2-y).

Parameters are:

        h = -.040;

EMap CT: 3