Skip to content

Commit 6445c3c

Browse files
committed
Merge branch 'master' into Kuramoto-Sivashinsky
2 parents 7e63fe2 + 83ac84f commit 6445c3c

File tree

109 files changed

+1297
-979
lines changed

Some content is hidden

Large Commits have some content hidden by default. Use the searchbox below for content that may be hidden.

109 files changed

+1297
-979
lines changed

README.md

Lines changed: 4 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -10,7 +10,7 @@ ODE Test Problems can be installed by running the MATLAB script `install` from t
1010

1111
```matlab
1212
% Create a problem
13-
problem = otp.brusselator.presets.Canonical;
13+
problem = otp.lotkavolterra.presets.Canonical;
1414
1515
% Solve the problem
1616
sol = problem.solve('RelTol', 1e-10);
@@ -19,15 +19,15 @@ sol = problem.solve('RelTol', 1e-10);
1919
problem.plot(sol);
2020
2121
% Adjust a parameter
22-
problem.Parameters.a = 2;
22+
problem.Parameters.preyDeathRate = 2;
2323
2424
% Manually use a MATLAB ODE solver to solve the problem
2525
options = odeset('Jacobian', problem.Rhs.Jacobian);
2626
[t, y] = ode15s(problem.Rhs.F, problem.TimeSpan, problem.Y0, options);
2727
2828
% Plot the phase space with a custom title
29-
problem.plotPhaseSpace(t, y, 'Title', 'Reactant X versus Reactant Y');
29+
problem.plotPhaseSpace(t, y, 'Title', 'The Circle of Life');
3030
3131
% Create a movie and write to file
32-
mov = problem.movie(t, y, 'Save', 'brusselator.avi');
32+
mov = problem.movie(t, y, 'Save', 'lotka-volterra.avi');
3333
```

logo.m

Lines changed: 3 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -5,10 +5,8 @@ function logo()
55
t = linspace(sol.x(1), sol.x(end), 450);
66
y = deval(sol, t);
77

8-
fig = p.plot(t, y', 'Title', 'ODE Test Problems', 'Legend', {});
8+
fig = p.plot(t, y, 'Title', 'ODE Test Problems', 'Legend', {}, 'ColorOrder', hsv(3), 'FontName', 'SansSerif');
99
fig.InvertHardcopy = 'off';
10-
ax = fig.CurrentAxes;
11-
ax.FontName = 'SansSerif';
12-
ax.Title.FontSize = 22;
10+
fig.CurrentAxes.Title.FontSize = 22;
1311

14-
end
12+
end

src/+otp/+allencahn/+presets/Canonical.m

Lines changed: 4 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -4,11 +4,11 @@
44
function obj = Canonical(varargin)
55

66
p = inputParser;
7-
addParameter(p, 'Size', 150, @isscalar);
8-
addParameter(p, 'alpha', 0.1, @isnumeric);
9-
addParameter(p, 'beta', 1, @isnumeric);
7+
p.addParameter('Size', 150, @isscalar);
8+
p.addParameter('alpha', 0.1);
9+
p.addParameter('beta', 1);
1010

11-
parse(p, varargin{:});
11+
p.parse(varargin{:});
1212

1313
s = p.Results;
1414

Lines changed: 32 additions & 32 deletions
Original file line numberDiff line numberDiff line change
@@ -1,37 +1,37 @@
11
classdef ArenstorfProblem < otp.Problem
2-
% The classic 3-body problem with parameters from
3-
% Reference: Hairer, Ernst, Syvert Paul Nørsett, and Gerhard Wanner.
4-
% Solving Ordinary Differential Equations I: Nonstiff Problems Springer-Verlag, 1987.
5-
% CH II, p. 129.
6-
methods
7-
function obj = ArenstorfProblem(timeSpan, y0, parameters)
8-
obj@otp.Problem('Arenstorf Orbit', 4, timeSpan, y0, parameters);
2+
% The classic 3-body problem with parameters from
3+
% Reference: Hairer, Ernst, Syvert Paul Nørsett, and Gerhard Wanner.
4+
% Solving Ordinary Differential Equations I: Nonstiff Problems Springer-Verlag, 1987.
5+
% CH II, p. 129.
6+
methods
7+
function obj = ArenstorfProblem(timeSpan, y0, parameters)
8+
obj@otp.Problem('Arenstorf Orbit', 4, timeSpan, y0, parameters);
9+
end
910
end
10-
end
11-
12-
methods (Access=protected)
1311

14-
function onSettingsChanged(obj)
15-
m1 = obj.Parameters.m1;
16-
m2 = obj.Parameters.m2;
17-
18-
obj.Rhs = otp.Rhs( @(t, y) otp.arenstorf.f(t, y, m1, m2), ...
19-
otp.Rhs.FieldNames.Jacobian, @(t, y) otp.arenstorf.jac(t, y, m1, m2));
12+
methods (Access=protected)
13+
14+
function onSettingsChanged(obj)
15+
m1 = obj.Parameters.m1;
16+
m2 = obj.Parameters.m2;
17+
18+
obj.Rhs = otp.Rhs( @(t, y) otp.arenstorf.f(t, y, m1, m2), ...
19+
otp.Rhs.FieldNames.Jacobian, @(t, y) otp.arenstorf.jac(t, y, m1, m2));
20+
end
21+
22+
23+
function validateNewState(obj, newTimeSpan, newY0, newParameters)
24+
validateNewState@otp.Problem(obj, newTimeSpan, newY0, newParameters);
25+
26+
otp.utils.StructParser(newParameters) ...
27+
.checkField('m1', 'scalar', 'real', 'finite', 'positive') ...
28+
.checkField('m2', 'scalar', 'real', 'finite', 'positive');
29+
end
30+
31+
function mov = internalMovie(obj, t, y, varargin)
32+
mov = otp.utils.movie.PhaseSpaceMovie('title', obj.Name, 'xlabel', 'x', 'ylabel', 'y', varargin{:});
33+
mov.record(t, y);
34+
end
35+
2036
end
21-
22-
23-
function validateNewState(obj, newTimeSpan, newY0, newParameters)
24-
validateNewState@otp.Problem(obj, newTimeSpan, newY0, newParameters);
25-
26-
otp.utils.StructParser(newParameters) ...
27-
.checkField('m1', 'scalar', 'real', 'finite', 'positive') ...
28-
.checkField('m2', 'scalar', 'real', 'finite', 'positive');
29-
end
30-
31-
function mov = internalMovie(obj, t, y, varargin)
32-
mov = otp.utils.movie.PhasePlaneMovie(obj.Name, @obj.index2label, varargin{:});
33-
mov.record(t, y(:,1:2));
34-
end
35-
36-
end
3737
end
Lines changed: 16 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -1,9 +1,9 @@
1-
classdef Canonical < otp.bouncingball.presets.Simple
1+
classdef Canonical < otp.bouncingball.BouncingBallProblem
22
% [Name]
3-
% Canonical
3+
% Flat Terrain
44
%
55
% [Description]
6-
% This simply calls the Flat Terrain example.
6+
% The ball has no horizontal velocity, and bounces on a perfectly flat surface.
77
%
88
% [NoVars]
99
% 4
@@ -12,4 +12,17 @@
1212
%
1313
% [Usage]
1414
%
15+
methods
16+
function obj = Canonical
17+
params.g = otp.utils.PhysicalConstants.EarthGravity;
18+
params.ground= @(x) 0;
19+
params.groundSlope= @(x) 0;
20+
21+
y0 = [0; 1; 0; 0];
22+
tspan = [0 10];
23+
24+
obj = obj@otp.bouncingball.BouncingBallProblem(tspan, y0, params);
25+
26+
end
27+
end
1528
end

src/+otp/+bouncingball/+presets/Parabola.m

Lines changed: 3 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -16,11 +16,11 @@
1616
methods
1717
function obj = Parabola
1818
params.g = otp.utils.PhysicalConstants.EarthGravity;
19-
params.groundFunction = @(x) x^2;
20-
params.dgroundFunction = @(x) 2*x;
19+
params.ground = @(x) x^2;
20+
params.groundSlope = @(x) 2*x;
2121

2222
y0 = [0; 1; 1; 0];
23-
tspan = [0, 50];
23+
tspan = [0, 10];
2424
obj = obj@otp.bouncingball.BouncingBallProblem(tspan, y0, params);
2525
end
2626
end

src/+otp/+bouncingball/+presets/RandomTerrain.m

Lines changed: 8 additions & 44 deletions
Original file line numberDiff line numberDiff line change
@@ -13,56 +13,20 @@
1313
function obj = RandomTerrain
1414
params.g = otp.utils.PhysicalConstants.EarthGravity;
1515

16-
n = 1000;
16+
n = 32;
1717

18-
rng(6);
18+
magnitudes = 0.2 * randn(n, 1);
19+
periods = randn(n, 1);
20+
phases = randn(n, 1);
1921

20-
genf = @cos;
21-
dgenf = @sin;
22+
params.ground = @(x) 0.05*x^2 + sum(magnitudes .* cos(periods * x + phases));
23+
params.groundSlope = @(x) 0.1*x - sum(magnitudes .* periods .* sin(periods * x + phases));
2224

23-
magnitudes = 1/32*randn(n, 1);
24-
periods = 0.9*randn(n, 1);
25-
additions = 0.25*randn(n, 1);
26-
27-
groundFunction = @(x) generateTerrain (x, genf, magnitudes, periods, additions) + 0.01*x^2;
28-
dgroundFunction = @(x) dgenerateTerrain(x, dgenf, magnitudes, periods, additions) + 0.02*x;
29-
30-
params.groundFunction = groundFunction;
31-
params.dgroundFunction = dgroundFunction;
32-
33-
y0 = [0; 5.5; 5; 0];
34-
tspan = [0; 100];
25+
y0 = [0; 5; 5; 0];
26+
tspan = [0; 50];
3527

3628
obj = obj@otp.bouncingball.BouncingBallProblem(tspan, y0, params);
3729

3830
end
3931
end
4032
end
41-
42-
function y = generateTerrain(x, genf, magnitudes, periods, additions)
43-
44-
y = 0;
45-
46-
for i = 1:numel(magnitudes)
47-
m = magnitudes(i);
48-
p = periods(i);
49-
a = additions(i);
50-
51-
y = y + m*genf(p*x + a);
52-
end
53-
54-
end
55-
56-
function y = dgenerateTerrain(x, dgenf, magnitudes, periods, additions)
57-
58-
y = 0;
59-
60-
for i = 1:numel(magnitudes)
61-
m = magnitudes(i);
62-
p = periods(i);
63-
a = additions(i);
64-
65-
y = y - m*p*dgenf(p*x + a);
66-
end
67-
68-
end

src/+otp/+bouncingball/+presets/Simple.m

Lines changed: 0 additions & 28 deletions
This file was deleted.

src/+otp/+bouncingball/BouncingBallMovie.m

Lines changed: 22 additions & 33 deletions
Original file line numberDiff line numberDiff line change
@@ -1,48 +1,37 @@
1-
classdef BouncingBallMovie < otp.utils.movie.Movie
1+
classdef BouncingBallMovie < otp.utils.movie.CometMovie
22
properties (SetAccess = immutable, GetAccess = private)
3-
MovieTitle
4-
GroundFunction
5-
end
6-
7-
properties (Access = protected)
8-
AnimatedLine
9-
Head
3+
Ground
104
end
115

126
methods
13-
function obj = BouncingBallMovie(title, groundFunction, varargin)
14-
obj@otp.utils.movie.Movie(varargin{:});
15-
obj.MovieTitle = title;
16-
obj.GroundFunction = groundFunction;
7+
function obj = BouncingBallMovie(ground, varargin)
8+
obj@otp.utils.movie.CometMovie(otp.utils.PhysicalConstants.TwoD, 'xlabel', 'x', 'ylabel', 'y', varargin{:});
9+
obj.Ground = ground;
1710
end
1811
end
1912

2013
methods (Access = protected)
21-
function init(obj, fig, state)
22-
ax = axes(fig);
23-
xlabel(ax, 'x');
24-
ylabel(ax, 'y');
25-
26-
xLim = otp.utils.FancyPlot.axisLimits('x', state.y(:, 1));
27-
hold(ax, 'on');
28-
29-
groundX = linspace(xLim(1), xLim(2), 1024).';
30-
groundY = arrayfun(@(x) obj.GroundFunction(x), groundX);
31-
plot(ax, groundX, groundY, 'k');
14+
function gObjects = initAxes(obj, ax, state)
15+
gObjects = initAxes@otp.utils.movie.CometMovie(obj, ax, state);
3216

33-
otp.utils.FancyPlot.axisLimits('y', [state.y(:, 2); groundY]);
17+
otp.utils.FancyPlot.axisLimits(ax, 'x', state.y(1, :));
18+
groundX = linspace(ax.XLim(1), ax.XLim(2));
19+
groundY = arrayfun(obj.Ground, groundX);
20+
line(ax, groundX, groundY, 'Color', 'k', 'LineWidth', 1.5);
3421

35-
color = otp.utils.FancyPlot.color(1);
36-
obj.AnimatedLine = animatedline(ax, 'Color', otp.utils.FancyPlot.brighten(color, 0.9));
37-
obj.Head = line(ax, 'Color', 'k', 'MarkerFaceColor', color, 'MarkerSize', 7, 'Marker', 'o', ...
38-
'LineStyle', 'none', 'XData', 0, 'YData', 0);
39-
hold(ax, 'off');
22+
otp.utils.FancyPlot.axisLimits(ax, 'y', [state.y(2, :), groundY]);
23+
end
24+
25+
function numComets = getNumComets(~, ~)
26+
numComets = 1;
27+
end
28+
29+
function x = getXPoints(~, ~, state)
30+
x = state.y(1, state.stepRange);
4031
end
4132

42-
function drawFrame(obj, fig, state)
43-
title(fig.CurrentAxes, sprintf('%s at t=%g', obj.MovieTitle, state.tCur));
44-
obj.AnimatedLine.addpoints(state.y(state.stepRange, 1), state.y(state.stepRange, 2));
45-
set(obj.Head, 'XData', state.yCur(1), 'YData', state.yCur(2));
33+
function y = getYPoints(~, ~, state)
34+
y = state.y(2, state.stepRange);
4635
end
4736
end
4837
end

src/+otp/+bouncingball/BouncingBallProblem.m

Lines changed: 8 additions & 8 deletions
Original file line numberDiff line numberDiff line change
@@ -8,12 +8,12 @@
88
methods (Access = protected)
99
function onSettingsChanged(obj)
1010
g = obj.Parameters.g;
11-
gF = obj.Parameters.groundFunction;
12-
dgF = obj.Parameters.dgroundFunction;
11+
ground = obj.Parameters.ground;
12+
slope = obj.Parameters.groundSlope;
1313

14-
obj.Rhs = otp.Rhs(@(t, y) otp.bouncingball.f(t, y, g, gF, dgF), ...
15-
otp.Rhs.FieldNames.Jacobian, @(t, y) otp.bouncingball.jac(t, y, g, gF, dgF), ...
16-
otp.Rhs.FieldNames.Events, @(t, y) otp.bouncingball.events(t, y, g, gF, dgF), ...
14+
obj.Rhs = otp.Rhs(@(t, y) otp.bouncingball.f(t, y, g, ground, slope), ...
15+
otp.Rhs.FieldNames.Jacobian, otp.bouncingball.jac(g, ground, slope), ...
16+
otp.Rhs.FieldNames.Events, @(t, y) otp.bouncingball.events(t, y, g, ground, slope), ...
1717
otp.Rhs.FieldNames.OnEvent, @otp.bouncingball.onevent);
1818
end
1919

@@ -22,8 +22,8 @@ function validateNewState(obj, newTimeSpan, newY0, newParameters)
2222

2323
otp.utils.StructParser(newParameters) ...
2424
.checkField('g', 'scalar', 'finite') ...
25-
.checkField('groundFunction', 'function') ...
26-
.checkField('dgroundFunction', 'function');
25+
.checkField('ground', 'function') ...
26+
.checkField('groundSlope', 'function');
2727
end
2828

2929
function label = internalIndex2label(~, index)
@@ -40,7 +40,7 @@ function validateNewState(obj, newTimeSpan, newY0, newParameters)
4040
end
4141

4242
function mov = internalMovie(obj, t, y, varargin)
43-
mov = otp.bouncingball.BouncingBallMovie(obj.Name, obj.Parameters.groundFunction, varargin{:});
43+
mov = otp.bouncingball.BouncingBallMovie(obj.Parameters.ground, 'Title', obj.Name, varargin{:});
4444
mov.record(t, y);
4545
end
4646

0 commit comments

Comments
 (0)