|
%{
|
|
File: test_phasePortrait.m
|
|
Author: Jonathan C. Teixeira
|
|
Email: jonta.teixeira@gmail.com
|
|
Description: Generate phase portrait of ode system
|
|
%}
|
|
|
|
% as usual
|
|
clc; clear; close all;
|
|
|
|
ode = 'nonlinear'; % ode selections: 'linear' or 'nonlinear'
|
|
|
|
switch lower(ode)
|
|
case 'linear'
|
|
% homogeneous linear system
|
|
% name of variables
|
|
varname = {'x', 'y'};
|
|
portrait_lim = {[-5,-5],[5,5]};
|
|
num_curves = 100;
|
|
[x0, y0] = deal(0,0);
|
|
|
|
a = [-1, +1, -1, +1];
|
|
b = [-5, -5, +5, +5];
|
|
|
|
for i=1:numel(a)
|
|
% define ode system
|
|
dy = @(x) [a(i)*x(1); b(i)*(x(1) + x(2))];
|
|
|
|
% plot phase portrait
|
|
phasePortrait(dy,portrait_lim{1},portrait_lim{2},num_curves)
|
|
hold on
|
|
plot(x0,y0,'sr');
|
|
xlabel(varname{1}); ylabel(varname{2}); grid on
|
|
print('-dpng', '-r300', ['phase_portrait_',ode,'_',num2str(i)]);
|
|
end
|
|
|
|
case 'nonlinear'
|
|
%% nonlinear system
|
|
% parameter definition
|
|
tau = 60;
|
|
k0 = 1.7e5;
|
|
G = 90.3;
|
|
E_R = 6307;
|
|
Ca0 = 1.5;
|
|
T0 = linspace(250,350,12);
|
|
|
|
% name of variables
|
|
varname = {'Concentration', 'Temperature'};
|
|
portrait_lim = {[0,200],[3,550]};
|
|
num_curves = 300;
|
|
|
|
for i=1:numel(T0)
|
|
% define ode system
|
|
dy = @(x) [(Ca0 - x(1))/tau - k0*exp(-E_R/x(2))*x(1); ...
|
|
(T0(i) - x(2))/tau + G*k0*exp(-E_R/x(2))*x(1)];
|
|
|
|
% critical points
|
|
opts = optimoptions('fsolve','Display','off','MaxFunctionEvaluations',3000, 'functiontolerance',1e-6);
|
|
x = [fsolve(dy,[0;0], opts), ...
|
|
fsolve(dy,[0;500], opts), ...
|
|
fsolve(dy,[1;380], opts)]
|
|
j = [sqrt(sum(dy(x(:,1)).^2))<1e-3; ...
|
|
sqrt(sum(dy(x(:,2)).^2))<1e-3; ...
|
|
sqrt(sum(dy(x(:,3)).^2))<1e-3]
|
|
disp(i);
|
|
y0 = x(2,j); x0 = x(1,j);
|
|
|
|
|
|
% plot phase portrait
|
|
phasePortrait(dy,portrait_lim{1},portrait_lim{2},num_curves, false)
|
|
hold on
|
|
plot(x0,y0,'sr','markerfacecolor','r');
|
|
xlabel(varname{1}); ylabel(varname{2}); grid on;
|
|
title(['T0 = ',num2str(T0(i))])
|
|
print('-dpng', '-r300', ['phase_portrait_',ode,'_',num2str(i)]);
|
|
end
|
|
otherwise
|
|
error('** select valid model')
|
|
end
|