plot::ode
-- plot the numerical
solution of an ordinary differential equationplot::ode
([t0, t1, ...], f, Y0, [G])
computes a mesh of numerical sample points Y(t0), Y(t1), ...
representing the solution Y(t) of the first order
differential equation (dynamical system)
dY/dt = f(t,Y), Y(t0) = Y0, t0 and t real, Y and Y0 complex vectors.
The procedure
G: (t, Y) -> [x(t,Y), y(t,Y)] or G: (t, Y) -> [x(t,Y), y(t,Y), z(t,Y)]
maps these solution points (t.i, Y(t.i)) in R x C^n to a mesh of 2D plot points [x.i, y.i] or 3D plot points [x.i, y.i, z.i], respectively,
plot::ode([t0, t1, ...], f, Y0, <method,> <RelativeError = tol,> <Stepsize = h,> [G1 <, Style = style1> <, Color = c1>], [G2 <, Style = style2> <, Color = c2>], ... )
t0, t1, ... |
- | the time mesh: real numerical values. If data are
displayed with Style = Splines , these values must be in ascending
order. |
f |
- | the vector field of the ODE: a procedure. See numeric::odesolve . |
Y0 |
- | the initial condition of the ODE: a list or a
1-dimensional array. See numeric::odesolve . |
G1, G2, ... |
- | ``generators of plot data'': procedures mapping a
solution point (t, Y(t)) to a list [x, y] or
[x, y, z] representing a plot point in 2D or 3D,
respectively. |
method |
- | use a specific numerical scheme (see numeric::odesolve ) |
RelativeError =
tol |
- | sets a numerical discretization tolerance (see
numeric::odesolve ) |
Stepsize = h |
- | sets a constant stepsize (see numeric::odesolve ) |
Style = style |
- | sets the style in which the plot data are displayed. The following styles are available: Points, Lines, Splines. The default style is Lines. |
Color = c |
- | sets the RGB color
c in which the plot data are displayed. The default color
is RGB::Black . |
a graphical object of domain type plot::Group
.
numeric::odesolve
, plot
, plot::Group
, plot::Scene
Y1 := numeric::odesolve(t0..t1, f, Y0 <,
Options>)
,
Y2 := numeric::odesolve(t1..t2, f, Y1 <,
Options>)
etc.
is computed where Options
is some combination of
method
, RelativeError=
tol
, and Stepsize= h
. See
numeric::odesolve
for details on the vector field procedure f
, the initial
condition Y0
, and the options.
G1
,
G2
etc. creates a graphical solution curve from the sample
points Y0
, Y1
etc. Each generator
G
, say, is internally called in the form G(t0,
Y0)
, G(t1, Y1)
, ...
to produce a
sequence of plot points in 2D or 3D.
The solver numeric::odesolve
returns the
solution points Y0
, Y1
etc. as lists or
1-dimensional arrays (the actual type is determined by the initial
value Y0
). Consequently, each generator G
must accept two arguments (t, Y)
: t
is a real
parameter, Y
is a ``vector'' (either a list or a
1-dimensional array).
Each generator must return a list with 2 or 3 elements representing
the (x, y) or (x, y, z) coordinates of the
graphical point associated with a solution point (t, Y)
of
the ODE. All generators must produce graphical data of the same
dimension, i.e., either 2D data as lists with 2 elements, or 3D data as
lists with 3 elements.
Some examples:
G := (t, Y) -> [t, Y[1]]
creates a 2D plot of the
first component of the solution vector along the y-axis,
plotted against the time variable t
along the
x-axis
G := (t, Y) -> [Y[1], Y[2]]
creates a 2D phase plot,
plotting the first component of the solution along the
x-axis and the second component along the y-axis.
The result is a solution curve in phase space (parametrized by the time
t
).
G := (t, Y) -> [Y[1], Y[2], Y[3]]
creates a 3D phase
plot of the first three components of the solution curve.
G
.
Cf. example 2.G1
, G2
etc. can be
specified to generate several curves associated with the same numerical
mesh Y0
, Y1
, ...
. E.g., one can
use a generator twice to plot the data with the options Style = Splines
and Style = Points
. This
allows to display a smooth solution curve together with the sample
points. Cf. examples 1, 2,
and 3.plot::ode
returns a graphical object of type plot::Group
containing all
graphical solution curves specified by the generators G1
,
G2
etc. This object can be combined with other plot
objects to a graphical scene via plot::Scene
. Cf. example 3. A final call to plot
calls the renderer to display the
scene.c
c
must be a list of three numerical real values [r, g, b]
between 0
and 1
representing the red, green
and blue contributions according to the RGB color model. A variety of
colors is provided by MuPAD's RGB
data structure.style
G1
, G2
etc. consists of a sequence of mesh
points in 2D or 3D, respectively.
With Style = Points
, the graphical data are displayed as a
discrete set of points.
With Style = Lines
, the graphical data points are displayed as a
curve consisting of straight line segments between the sample points.
The points themselves are not displayed.
With Style = Splines
, the graphical data points are displayed as
a smooth spline curve connecting the sample points. The points
themselves are not displayed.
The following procedure f
together with the
initial value Y0
represent the initial value problem
dY/dt = f(t, Y) = t*Y - Y^2, Y(t0) = 1. In
MuPAD, the function Y is represented by a list with
one element:
>> f := (t, Y) -> [t*Y[1] - Y[1]^2]: Y0 := [1]:
The solution is to be plotted as a function of the time
t
. We define an appropriate generator for the plot
data:
>> G := (t, Y) -> [t, Y[1]]:
The numerical solution is to consist of sample points
over the time mesh t.i = i, i = 0,1,..,10. The
generator G
is used twice with different Style options. This generates the sample points together
with a smooth spline curve connecting these points:
>> p := plot::ode([i $ i = 0..10], f, Y0, [G, Style = Points, Color = RGB::Blue], [G, Style = Splines, Color = RGB::Red])
plot::Group()
The point size of the point objects inside
p
is increased:
>> p:= plot::modify(p, PointWidth = 50):
The object p
is converted to a graphical
scene with various scene
options:
>> s := plot::Scene(p, Labels = ["t", "y"], Ticks = [Steps = [1.0, 1], Steps = [1.0, 1]], GridLines = Automatic)
plot::Scene()
Finally, the scene is rendered by a call to plot
:
>> plot(s)
>> delete f, Y0, G, p, s:
We consider the nonlinear oscillator y'' + y + y^3 = sin(t), y(0) = 0, y'(0) = 0.5. As a dynamical system for Y = [y, y'], we have to solve the following initial value problem dY/dt = f(t,Y), Y(0) = Y0:
>> f := (t, Y) -> [Y[2], sin(t) - Y[1]^3]: Y0 := [0, 0.5]:
The following generator produces a phase plot in the (x,y) plane, embedded in a 3D plot:
>> G1 := (t, Y) -> [Y[1], Y[2], 0]:
Further, we use the z coordinate of the 3D plot to display the value of the ``energy'' function E = y^2/2 + y'^2/2 over the phase curve:
>> G2 := (t, Y) -> [Y[1], Y[2], (Y[1]^2 + Y[2]^2)/2]:
The phase curve in the (x, y) plane is combined with the graph of the energy function:
>> p := plot::ode([i/5 $ i = 0..100], f, Y0, [G1, Style = Splines, Color = RGB::Red], [G2, Style = Points, Color = RGB::Black], [G2, Style = Lines, Color = RGB::Blue]):
We increase the size of the points used in the representation of the energy:
>> p:= plot::modify(p, PointWidth = 40):
The renderer is called:
>> plot(p, Axes = Box, Labels = ["y", "y'", "E"], CameraPoint = [10, -15, 5]):
>> delete f, Y0, G1, G2, p:
We consider the initial value problem y' = f(t, y) = t*sin(t + y^2), y(0) = 1.
>> f := (t, y) -> t*sin(t + y^2) - y:
The following vector field is tangent to the solution curves:
>> p1 := plot::vectorfield([1, f(t, y)], t = 1..3, y = 0..1, Grid = [10, 5], Color = RGB::Black):
The following object represents the plot of the solution
as a function of t
:
>> p2 := plot::ode( [i/3 $ i=0..12], (t,Y) -> [f(t, Y[1])], [0], [(t, Y) -> [t, Y[1]], Style = Points, Color = RGB::Red], [(t, Y) -> [t, Y[1]], Style = Splines, Color = RGB::Blue]):
We increase the point size:
>> p2:= plot::modify(p2, PointWidth = 40):
Finally, we combine the vector field and the ODE plot to a scene and call the renderer:
>> plot(p1, p2, Scaling = Constrained, Labels = ["t", "y"], GridLines = [Steps = 0.25, Steps = 0.25], Ticks = [Steps = 1.0, Steps = 0.5]):
>> delete f, p1, p2:
The Lorenz ODE is the system
d/dt [x, y, z] = [p*(y - x), -x*z + r*x - y, x*y - b*z]
with fixed parameters p, r, b. As a dynamical system for Y = [x, y, z], we have to solve the ODE dY/dt = f(t,Y) with the following vector field:
>> f := proc(t, Y) local x, y, z; begin [x, y, z] := Y: [p*(y - x), -x*z + r*x - y, x*y - b*z] end_proc:
We consider the following parameters and the following
initial condition Y0
:
>> p := 10: r := 28: b := 1: Y0 := [1, 1, 1]:
The following generator produces a 3D phase plot:
>> Gxyz := (t, Y) -> Y:
The following generator projects the solution curve to the (y, z) plane with x = -20:
>> Gyz := (t, Y) -> [-20, Y[2], Y[3]]:
The following generator projects the solution curve to the (x, z) plane with y = 30:
>> Gxz := (t, Y) -> [Y[1], 30, Y[3]]:
The following generator projects the solution curve to the (x, y) plane with z = 0:
>> Gxy := (t, Y) -> [Y[1], Y[2], 0]:
With these generators, we create a 3D plot object consisting of the phase curve and its projections. The following call is somewhat time consuming, because it calls the numerical integrator to produce the graphical data:
>> obj := plot::ode([i/10 $ i=1..100], f, Y0, [Gxyz, Style = Splines, Color = RGB::Red], [Gyz, Style = Splines, Color = RGB::LightGrey], [Gxz, Style = Splines, Color = RGB::LightGrey], [Gxy, Style = Splines, Color = RGB::LightGrey]):
Finally, the plot is rendered. Again, this call is somewhat time consuming, because internally, several thousand calls to spline functions occur:
>> plot(obj, Axes = Box, Ticks = 4, CameraPoint = [400, -800, 1200]):
>> delete f, p, r, b, Y0, Gxyz, Gyz, Gxz, Gxy, obj: