1########################################################################
2##
3## Copyright (C) 2006-2021 The Octave Project Developers
4##
5## See the file COPYRIGHT.md in the top-level directory of this
6## distribution or <https://octave.org/copyright/>.
7##
8## This file is part of Octave.
9##
10## Octave is free software: you can redistribute it and/or modify it
11## under the terms of the GNU General Public License as published by
12## the Free Software Foundation, either version 3 of the License, or
13## (at your option) any later version.
14##
15## Octave is distributed in the hope that it will be useful, but
16## WITHOUT ANY WARRANTY; without even the implied warranty of
17## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
18## GNU General Public License for more details.
19##
20## You should have received a copy of the GNU General Public License
21## along with Octave; see the file COPYING.  If not, see
22## <https://www.gnu.org/licenses/>.
23##
24########################################################################
25
26## -*- texinfo -*-
27## @deftypefn {} {@var{stop_solve} =} odeplot (@var{t}, @var{y}, @var{flag})
28##
29## Open a new figure window and plot the solution of an ode problem at each
30## time step during the integration.
31##
32## The types and values of the input parameters @var{t} and @var{y} depend on
33## the input @var{flag} that is of type string.  Valid values of @var{flag}
34## are:
35##
36## @table @option
37## @item @qcode{"init"}
38## The input @var{t} must be a column vector of length 2 with the first and
39## last time step (@code{[@var{tfirst} @var{tlast}]}.  The input @var{y}
40## contains the initial conditions for the ode problem (@var{y0}).
41##
42## @item @qcode{""}
43## The input @var{t} must be a scalar double specifying the time for which
44## the solution in input @var{y} was calculated.
45##
46## @item @qcode{"done"}
47## The inputs should be empty, but are ignored if they are present.
48## @end table
49##
50## @code{odeplot} always returns false, i.e., don't stop the ode solver.
51##
52## Example: solve an anonymous implementation of the
53## @nospell{@qcode{"Van der Pol"}} equation and display the results while
54## solving.
55##
56## @example
57## @group
58## fvdp = @@(t,y) [y(2); (1 - y(1)^2) * y(2) - y(1)];
59##
60## opt = odeset ("OutputFcn", @@odeplot, "RelTol", 1e-6);
61## sol = ode45 (fvdp, [0 20], [2 0], opt);
62## @end group
63## @end example
64##
65## Background Information:
66## This function is called by an ode solver function if it was specified in
67## the @qcode{"OutputFcn"} property of an options structure created with
68## @code{odeset}.  The ode solver will initially call the function with the
69## syntax @code{odeplot ([@var{tfirst}, @var{tlast}], @var{y0}, "init")}.  The
70## function initializes internal variables, creates a new figure window, and
71## sets the x limits of the plot.  Subsequently, at each time step during the
72## integration the ode solver calls @code{odeplot (@var{t}, @var{y}, [])}.
73## At the end of the solution the ode solver calls
74## @code{odeplot ([], [], "done")} so that odeplot can perform any clean-up
75## actions required.
76## @seealso{odeset, odeget, ode23, ode45}
77## @end deftypefn
78
79function stop_solve = odeplot (t, y, flag)
80
81  ## No input argument checking is done for better performance
82  persistent hlines num_lines told yold;
83
84  ## odeplot never stops the integration
85  stop_solve = false;
86
87  if (isempty (flag))
88    ## Default case, plot and return a value
89    told = [told; t(:)];
90    yold = [yold, y];
91    for i = 1:num_lines
92      set (hlines(i), "xdata", told, "ydata", yold(i,:));
93    endfor
94    drawnow;
95
96    retval = false;
97
98  elseif (strcmp (flag, "init"))
99    ## t is either the time slot [tstart tstop] or [t0, t1, ..., tn]
100    ## y is the initial value vector for the ode solution
101    told = t(1);
102    yold = y(:);
103    figure ();
104    hlines = plot (told, yold, "o-");
105    xlim ([t(1), t(end)]);  # Fix limits which also speeds up plotting
106    num_lines = numel (hlines);
107
108  elseif (strcmp (flag, "done"))
109    ## Cleanup after ode solver has finished.
110    hlines = num_lines = told = yold = [];
111
112  endif
113
114endfunction
115
116
117%!demo
118%! ## Solve an anonymous implementation of the Van der Pol equation
119%! ## and display the results while solving
120%! fvdp = @(t,y) [y(2); (1 - y(1)^2) * y(2) - y(1)];
121%! opt = odeset ("OutputFcn", @odeplot, "RelTol", 1e-6);
122%! sol = ode45 (fvdp, [0 20], [2 0], opt);
123