Mercurial > octave

annotate scripts/ode/private/integrate_adaptive.m @ 21546:f7f97d7e9294

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
doc: Wrap m-file docstrings to 79 characters + newline (80 total). * isrecording.m, soundsc.m, delaunay3.m, cell2mat.m, cumtrapz.m, del2.m, inputParser.m, interp1.m, interp3.m, narginchk.m, profile.m, validateattributes.m, delaunayn.m, tsearchn.m, voronoin.m, brighten.m, cmunique.m, colorcube.m, imfinfo.m, imshow.m, edit.m, orderfields.m, run.m, warning_ids.m, ode23.m, ode45.m, odeget.m, integrate_adaptive.m, kahan.m, ode_struct_value_check.m, runge_kutta_23.m, fminunc.m, fsolve.m, fzero.m, pkg.m, build.m, specular.m, view.m, bar.m, barh.m, contour3.m, isosurface.m, line.m, pie.m, pie3.m, quiver3.m, scatter.m, scatter3.m, stem3.m, stemleaf.m, surfl.m, tetramesh.m, isfigure.m, mkpp.m, pchip.m, residue.m, splinefit.m, rmpref.m, unique.m, eigs.m, ilu.m, factor.m, factorial.m, gallery.m, hankel.m, histc.m, ols.m, finv.m, fpdf.m, kruskal_wallis_test.m, weekday.m: Wrap m-file docstrings to 79 characters + newline (80 total).
author Rik <rik@octave.org>
date Sun, 27 Mar 2016 15:50:01 -0700
parents 516bb87ea72e
children 0f3e875d9078
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
20621
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
1 ## Copyright (C) 2015 Carlo de Falco
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
2 ## Copyright (C) 2013 Roberto Porcu' <roberto.porcu@polimi.it>
20536
6256f6e366ac Fix copyright text in private ode functions
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20533
diff changeset
3 ##
6256f6e366ac Fix copyright text in private ode functions
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20533
diff changeset
4 ## This file is part of Octave.
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
5 ##
20536
6256f6e366ac Fix copyright text in private ode functions
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20533
diff changeset
6 ## Octave is free software; you can redistribute it and/or modify it
6256f6e366ac Fix copyright text in private ode functions
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20533
diff changeset
7 ## under the terms of the GNU General Public License as published by
6256f6e366ac Fix copyright text in private ode functions
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20533
diff changeset
8 ## the Free Software Foundation; either version 3 of the License, or (at
6256f6e366ac Fix copyright text in private ode functions
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20533
diff changeset
9 ## your option) any later version.
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
10 ##
20536
6256f6e366ac Fix copyright text in private ode functions
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20533
diff changeset
11 ## Octave is distributed in the hope that it will be useful, but
6256f6e366ac Fix copyright text in private ode functions
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20533
diff changeset
12 ## WITHOUT ANY WARRANTY; without even the implied warranty of
6256f6e366ac Fix copyright text in private ode functions
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20533
diff changeset
13 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
6256f6e366ac Fix copyright text in private ode functions
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20533
diff changeset
14 ## General Public License for more details.
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
15 ##
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
16 ## You should have received a copy of the GNU General Public License
20536
6256f6e366ac Fix copyright text in private ode functions
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20533
diff changeset
17 ## along with Octave; see the file COPYING. If not, see
6256f6e366ac Fix copyright text in private ode functions
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20533
diff changeset
18 ## <http://www.gnu.org/licenses/>.
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
19
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
20 ## -*- texinfo -*-
20852
516bb87ea72e 2015 Code Sprint: remove class of function from docstring for all m-files.
Rik <rik@octave.org>
parents: 20847
diff changeset
21 ## @deftypefn {} {@var{solution} =} integrate_adaptive (@var{@@stepper}, @var{order}, @var{@@func}, @var{tspan}, @var{x0}, @var{options})
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
22 ##
20548
25623ef2ff4f doc: Rewrite docstrings for ode* family of functions.
Rik <rik@octave.org>
parents: 20543
diff changeset
23 ## This function file can be called by an ODE solver function in order to
20621
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
24 ## integrate the set of ODEs on the interval @var{[t0, t1]} with an
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
25 ## adaptive timestep.
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
26 ##
20621
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
27 ## The function returns a structure @var{solution} with two fieldss: @var{t}
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
28 ## and @var{y}. @var{t} is a column vector and contains the time stamps.
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
29 ## @var{y} is a matrix in which each column refers to a different unknown
20621
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
30 ## of the problem and the row number is the same as the @var{t} row number.
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
31 ## Thus, each row of the matrix @var{y} contains the values of all unknowns at
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
32 ## the time value contained in the corresponding row in @var{t}.
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
33 ##
20621
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
34 ## The first input argument must be a function handle or inline function
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
35 ## representing the stepper, i.e., the function responsible for step-by-step
20548
25623ef2ff4f doc: Rewrite docstrings for ode* family of functions.
Rik <rik@octave.org>
parents: 20543
diff changeset
36 ## integration. This function discriminates one method from the others.
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
37 ##
20548
25623ef2ff4f doc: Rewrite docstrings for ode* family of functions.
Rik <rik@octave.org>
parents: 20543
diff changeset
38 ## The second input argument is the order of the stepper. It is needed
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
39 ## to compute the adaptive timesteps.
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
40 ##
20621
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
41 ## The third input argument is a function handle or inline function that
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
42 ## defines the ODE:
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
43 ##
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
44 ## @ifhtml
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
45 ## @example
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
46 ## @math{y' = f(t,y)}
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
47 ## @end example
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
48 ## @end ifhtml
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
49 ## @ifnothtml
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
50 ## @math{y' = f(t,y)}.
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
51 ## @end ifnothtml
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
52 ##
20621
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
53 ## The fourth input argument is the time vector which defines the integration
20716
1ecee53513d7 doc: Peridodic grammar check of documentation.
Rik <rik@octave.org>
parents: 20636
diff changeset
54 ## interval, i.e., @var{[tspan(1), tspan(end)]} and all intermediate elements
20621
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
55 ## are taken as times at which the solution is required.
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
56 ##
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
57 ## The fifth argument represents the initial conditions for the ODEs and the
21546
f7f97d7e9294 doc: Wrap m-file docstrings to 79 characters + newline (80 total).
Rik <rik@octave.org>
parents: 20852
diff changeset
58 ## last input argument contains some options that may be needed for the
f7f97d7e9294 doc: Wrap m-file docstrings to 79 characters + newline (80 total).
Rik <rik@octave.org>
parents: 20852
diff changeset
59 ## stepper.
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
60 ##
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
61 ## @end deftypefn
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
62 ##
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
63 ## @seealso{integrate_const, integrate_n_steps}
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
64
20621
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
65 function solution = integrate_adaptive (stepper, order, func, tspan, x0,
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
66 options)
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
67
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
68 fixed_times = numel (tspan) > 2;
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
69
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
70 t_new = t_old = t = tspan(1);
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
71 x_new = x_old = x = x0(:);
20602
756b052037fb avoid stepping beyond end of thspan in ode solvers
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20601
diff changeset
72
20621
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
73 ## Get first initial timestep
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
74 dt = odeget (options, "InitialStep", [], "fast");
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
75 if (isempty (dt))
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
76 dt = starting_stepsize (order, func, t, x, options.AbsTol, options.RelTol,
20634
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
77 options.normcontrol);
20621
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
78 endif
20602
756b052037fb avoid stepping beyond end of thspan in ode solvers
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20601
diff changeset
79
20634
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
80 dir = odeget (options, "direction", [], "fast");
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
81 dt = dir * min (abs (dt), options.MaxStep);
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
82
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
83 options.comp = 0.0;
20602
756b052037fb avoid stepping beyond end of thspan in ode solvers
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20601
diff changeset
84
20552
eb9e2d187ed2 maint: Use Octave coding conventions in scripts/ode/private dir.
Rik <rik@octave.org>
parents: 20548
diff changeset
85 ## Factor multiplying the stepsize guess
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
86 facmin = 0.8;
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
87 facmax = 1.5;
20552
eb9e2d187ed2 maint: Use Octave coding conventions in scripts/ode/private dir.
Rik <rik@octave.org>
parents: 20548
diff changeset
88 fac = 0.38^(1/(order+1)); # formula taken from Hairer
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
89
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
90 ## Initialize the OutputFcn
20634
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
91 if (options.haveoutputfunction)
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
92 if (options.haveoutputselection)
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
93 solution.retout = x(options.OutputSel,end);
20602
756b052037fb avoid stepping beyond end of thspan in ode solvers
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20601
diff changeset
94 else
20634
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
95 solution.retout = x;
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
96 endif
20634
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
97 feval (options.OutputFcn, tspan, solution.retout,
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
98 "init", options.funarguments{:});
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
99 endif
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
100
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
101 ## Initialize the EventFcn
20634
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
102 if (options.haveeventfunction)
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
103 ode_event_handler (options.Events, tspan(end), x,
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
104 "init", options.funarguments{:});
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
105 endif
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
106
20634
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
107 if (options.havenonnegative)
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
108 nn = options.NonNegative;
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
109 endif
20602
756b052037fb avoid stepping beyond end of thspan in ode solvers
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20601
diff changeset
110
20634
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
111 solution.cntloop = 2;
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
112 solution.cntcycles = 1;
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
113 solution.cntsave = 2;
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
114 solution.unhandledtermination = true;
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
115 ireject = 0;
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
116
20602
756b052037fb avoid stepping beyond end of thspan in ode solvers
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20601
diff changeset
117 k_vals = [];
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
118 iout = istep = 1;
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
119 while (dir * t_old < dir * tspan(end))
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
120
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
121 ## Compute integration step from t_old to t_new = t_old + dt
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
122 [t_new, options.comp] = kahan (t_old, options.comp, dt);
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
123 [t_new, x_new, x_est, new_k_vals] = ...
20621
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
124 stepper (func, t_old, x_old, dt, options, k_vals, t_new);
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
125
20636
6e81f4b37e13 Performance improvements for ODE functions.
Rik <rik@octave.org>
parents: 20634
diff changeset
126 solution.cntcycles += 1;
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
127
20634
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
128 if (options.havenonnegative)
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
129 x_new(nn, end) = abs (x_new(nn, end));
20604
a260a6acb70f fix test failures introduced by a22d8a2eb0e5
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20602
diff changeset
130 x_est(nn, end) = abs (x_est(nn, end));
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
131 endif
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
132
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
133 err = AbsRel_Norm (x_new, x_old, options.AbsTol, options.RelTol,
20634
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
134 options.normcontrol, x_est);
20602
756b052037fb avoid stepping beyond end of thspan in ode solvers
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20601
diff changeset
135
20621
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
136 ## Accept solution only if err <= 1.0
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
137 if (err <= 1)
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
138
20636
6e81f4b37e13 Performance improvements for ODE functions.
Rik <rik@octave.org>
parents: 20634
diff changeset
139 solution.cntloop += 1;
6e81f4b37e13 Performance improvements for ODE functions.
Rik <rik@octave.org>
parents: 20634
diff changeset
140 ireject = 0; # Clear reject counter
20602
756b052037fb avoid stepping beyond end of thspan in ode solvers
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20601
diff changeset
141
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
142 ## if output time steps are fixed
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
143 if (fixed_times)
20543
3339c9bdfe6a Activate FSAL property in dorpri timestepper
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20536
diff changeset
144
20604
a260a6acb70f fix test failures introduced by a22d8a2eb0e5
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20602
diff changeset
145 t_caught = find ((dir * tspan(iout:end) > dir * t_old)
a260a6acb70f fix test failures introduced by a22d8a2eb0e5
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20602
diff changeset
146 & (dir * tspan(iout:end) <= dir * t_new));
20601
43822bda4f65 fix indexing bug introduced with a22d8a2eb0e5
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20600
diff changeset
147 t_caught = t_caught + iout - 1;
20602
756b052037fb avoid stepping beyond end of thspan in ode solvers
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20601
diff changeset
148
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
149 if (! isempty (t_caught))
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
150 t(t_caught) = tspan(t_caught);
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
151 iout = max (t_caught);
20601
43822bda4f65 fix indexing bug introduced with a22d8a2eb0e5
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20600
diff changeset
152 x(:, t_caught) = ...
20634
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
153 runge_kutta_interpolate (order, [t_old t_new], [x_old x_new],
20847
ddc18b909ec7 codesprint: ode: first and second order interpolator for dense_output
jcorno <jacopo.corno@gmail.com>
parents: 20716
diff changeset
154 tspan(t_caught), new_k_vals, dt, func,
ddc18b909ec7 codesprint: ode: first and second order interpolator for dense_output
jcorno <jacopo.corno@gmail.com>
parents: 20716
diff changeset
155 options.funarguments);
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
156
20636
6e81f4b37e13 Performance improvements for ODE functions.
Rik <rik@octave.org>
parents: 20634
diff changeset
157 istep += 1;
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
158
20621
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
159 ## Call Events function only if a valid result has been found.
20634
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
160 ## Stop integration if eventbreak is true.
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
161 if (options.haveeventfunction)
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
162 break_loop = false;
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
163 for idenseout = 1:numel (t_caught)
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
164 id = t_caught(idenseout);
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
165 td = t(id);
20634
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
166 solution.event = ...
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
167 ode_event_handler (options.Events, t(id), x(:, id), [],
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
168 options.funarguments{:});
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
169 if (! isempty (solution.event{1}) && solution.event{1} == 1)
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
170 t(id) = solution.event{3}(end);
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
171 t = t(1:id);
20634
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
172 x(:, id) = solution.event{4}(end, :).';
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
173 x = x(:,1:id);
20634
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
174 solution.unhandledtermination = false;
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
175 break_loop = true;
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
176 break;
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
177 endif
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
178 endfor
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
179 if (break_loop)
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
180 break;
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
181 endif
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
182 endif
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
183
20621
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
184 ## Call OutputFcn only if a valid result has been found.
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
185 ## Stop integration if function returns false.
20634
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
186 if (options.haveoutputfunction)
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
187 cnt = options.Refine + 1;
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
188 approxtime = linspace (t_old, t_new, cnt);
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
189 approxvals = interp1 ([t_old, t(t_caught), t_new],
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
190 [x_old, x(:, t_caught), x_new] .',
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
191 approxtime, 'linear') .';
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
192 if (options.haveoutputselection)
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
193 approxvals = approxvals(options.OutputSel, :);
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
194 endif
20634
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
195 for ii = 1:numel (approxtime)
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
196 pltret = feval (options.OutputFcn, approxtime(ii),
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
197 approxvals(:, ii), [],
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
198 options.funarguments{:});
20604
a260a6acb70f fix test failures introduced by a22d8a2eb0e5
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20602
diff changeset
199 endfor
20634
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
200 if (pltret) # Leave main loop
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
201 solution.unhandledtermination = false;
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
202 break;
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
203 endif
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
204 endif
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
205
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
206 endif
20602
756b052037fb avoid stepping beyond end of thspan in ode solvers
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20601
diff changeset
207
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
208 else
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
209
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
210 t(++istep) = t_new;
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
211 x(:, istep) = x_new;
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
212 iout = istep;
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
213
20621
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
214 ## Call Events function only if a valid result has been found.
20634
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
215 ## Stop integration if eventbreak is true.
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
216 if (options.haveeventfunction)
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
217 solution.event = ...
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
218 ode_event_handler (options.Events, t(istep), x(:, istep), [],
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
219 options.funarguments{:});
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
220 if (! isempty (solution.event{1}) && solution.event{1} == 1)
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
221 t(istep) = solution.event{3}(end);
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
222 x(:, istep) = solution.event{4}(end, :).';
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
223 solution.unhandledtermination = false;
20621
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
224 break;
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
225 endif
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
226 endif
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
227
20621
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
228 ## Call OutputFcn only if a valid result has been found.
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
229 ## Stop integration if function returns false.
20634
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
230 if (options.haveoutputfunction)
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
231 cnt = options.Refine + 1;
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
232 approxtime = linspace (t_old, t_new, cnt);
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
233 approxvals = interp1 ([t_old, t_new],
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
234 [x_old, x_new] .',
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
235 approxtime, 'linear') .';
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
236 if (options.haveoutputselection)
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
237 approxvals = approxvals(options.OutputSel, :);
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
238 endif
20634
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
239 for ii = 1:numel (approxtime)
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
240 pltret = feval (options.OutputFcn, approxtime(ii),
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
241 approxvals(:, ii), [], options.funarguments{:});
20604
a260a6acb70f fix test failures introduced by a22d8a2eb0e5
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20602
diff changeset
242 endfor
20634
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
243 if (pltret) # Leave main loop
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
244 solution.unhandledtermination = false;
20552
eb9e2d187ed2 maint: Use Octave coding conventions in scripts/ode/private dir.
Rik <rik@octave.org>
parents: 20548
diff changeset
245 break;
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
246 endif
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
247 endif
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
248
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
249 endif
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
250
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
251 ## move to next time-step
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
252 t_old = t_new;
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
253 x_old = x_new;
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
254 k_vals = new_k_vals;
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
255
20634
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
256 solution.cntloop += 1;
20602
756b052037fb avoid stepping beyond end of thspan in ode solvers
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20601
diff changeset
257
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
258 else
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
259
20636
6e81f4b37e13 Performance improvements for ODE functions.
Rik <rik@octave.org>
parents: 20634
diff changeset
260 ireject += 1;
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
261
20621
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
262 ## Stop solving because, in the last 5,000 steps, no successful valid
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
263 ## value has been found
20621
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
264 if (ireject >= 5_000)
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
265 error (["integrate_adaptive: Solving was not successful. ", ...
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
266 " The iterative integration loop exited at time", ...
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
267 " t = %f before the endpoint at tend = %f was reached. ", ...
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
268 " This happened because the iterative integration loop", ...
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
269 " did not find a valid solution at this time stamp. ", ...
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
270 " Try to reduce the value of 'InitialStep' and/or", ...
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
271 " 'MaxStep' with the command 'odeset'.\n"],
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
272 t_old, tspan(end));
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
273 endif
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
274
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
275 endif
20602
756b052037fb avoid stepping beyond end of thspan in ode solvers
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20601
diff changeset
276
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
277 ## Compute next timestep, formula taken from Hairer
20621
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
278 err += eps; # avoid divisions by zero
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
279 dt *= min (facmax, max (facmin, fac * (1 / err)^(1 / (order + 1))));
20602
756b052037fb avoid stepping beyond end of thspan in ode solvers
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20601
diff changeset
280 dt = dir * min (abs (dt), options.MaxStep);
756b052037fb avoid stepping beyond end of thspan in ode solvers
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20601
diff changeset
281
756b052037fb avoid stepping beyond end of thspan in ode solvers
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20601
diff changeset
282 ## make sure we don't go past tpan(end)
756b052037fb avoid stepping beyond end of thspan in ode solvers
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20601
diff changeset
283 dt = dir * min (abs (dt), abs (tspan(end) - t_old));
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
284
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
285 endwhile
20602
756b052037fb avoid stepping beyond end of thspan in ode solvers
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20601
diff changeset
286
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
287 ## Check if integration of the ode has been successful
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
288 if (dir * t(end) < dir * tspan(end))
20634
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
289 if (solution.unhandledtermination == true)
20601
43822bda4f65 fix indexing bug introduced with a22d8a2eb0e5
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20600
diff changeset
290 error ("integrate_adaptive:unexpected_termination",
20621
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
291 [" Solving was not successful. ", ...
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
292 " The iterative integration loop exited at time", ...
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
293 " t = %f before the endpoint at tend = %f was reached. ", ...
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
294 " This may happen if the stepsize becomes too small. ", ...
20604
a260a6acb70f fix test failures introduced by a22d8a2eb0e5
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20602
diff changeset
295 " Try to reduce the value of 'InitialStep'", ...
20634
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
296 " and/or 'MaxStep' with the command 'odeset'."],
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
297 t(end), tspan(end));
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
298 else
20601
43822bda4f65 fix indexing bug introduced with a22d8a2eb0e5
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20600
diff changeset
299 warning ("integrate_adaptive:unexpected_termination",
20621
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
300 ["Solver was stopped by a call of 'break'", ...
20634
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
301 " in the main iteration loop at time", ...
20621
b92f8e148936 maint: Continued clean-up of functions in ode/private dir.
Rik <rik@octave.org>
parents: 20604
diff changeset
302 " t = %f before the endpoint at tend = %f was reached. ", ...
20634
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
303 " This may happen because the @odeplot function", ...
80e630b37ba1 maint: Remove unnecessary 'v' prefix before variables in ODE m-files.
Rik <rik@octave.org>
parents: 20621
diff changeset
304 " returned 'true' or the @event function returned 'true'."],
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
305 t(end), tspan(end));
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
306 endif
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
307 endif
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
308
20601
43822bda4f65 fix indexing bug introduced with a22d8a2eb0e5
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20600
diff changeset
309 ## Set up return structure
43822bda4f65 fix indexing bug introduced with a22d8a2eb0e5
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20600
diff changeset
310 solution.t = t(:);
20600
a22d8a2eb0e5 fix adaptive strategy in ode solvers.
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20565
diff changeset
311 solution.x = x.';
20602
756b052037fb avoid stepping beyond end of thspan in ode solvers
Carlo de Falco <carlo.defalco@polimi.it>
parents: 20601
diff changeset
312
20533
fcb792acab9b Moving ode45, odeset, odeget, and levenshtein from odepkg to core.
jcorno <jacopo.corno@gmail.com>
parents:
diff changeset
313 endfunction
20552
eb9e2d187ed2 maint: Use Octave coding conventions in scripts/ode/private dir.
Rik <rik@octave.org>
parents: 20548
diff changeset
314