FazBrowse GitHub Viewer
|
Trending
|
URL:
|
Home
Tools:
[Download Repo ZIP]
[View Raw Code]
[Original HTTPS Page]
phase_plot/phase_plot.m at master · auralius/phase_plot · GitHub
auralius
/
phase_plot
Public
Notifications
You must be signed in to change notification settings
Fork
9
Star
18
Code
Issues
0
Pull requests
0
Actions
Projects
Security and quality
0
Insights
Additional navigation options
Code
Issues
Pull requests
Actions
Projects
Security and quality
Insights
Expand file tree
Breadcrumbs
phase_plot
/
phase_plot.m
Copy path
More file actions
More file actions
Latest commit
History
History
History
115 lines (98 loc) · 2.96 KB
Breadcrumbs
phase_plot
/
phase_plot.m
Copy path
File metadata and controls
115 lines (98 loc) · 2.96 KB
Raw
Copy raw file
Download raw file
Open symbols panel
Edit and raw actions
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
function
phase_plot(
sys
,
intial_values
,
range
,
simtime
,
scale
)
%
Phase portrait plot for SECOND and THIRD order ODE
%
sys is the system transfer function (in s-domain)
%
%
intial_values is ithe initial states of th system (vector of nx1)
%
where n is the order of the system
%
%
range is the minimum and the maximum boundary for the states
%
e.g: 3rd order system, with states: x1, x2, and x3
%
[x1_min x1_max; x2_min x2_max; x3_min x3_max]
%
%
simtime is the simulation time
%
%
scale is used to adjust the dimension of the arrows
%
this corresponds to the AutoScale property of the quiver function
%
%
References:
%
http://matlab.cheme.cmu.edu/2011/08/09/phase-portraits-of-a-system-of-odes/
if
nargin
<
5
scale
=
0.5
;
end
%%
Solve the ODE
[
num
,
den
]
=
tfdata(
sys
);
n
=
length(den{
:
})
-
1
;
fprintf(
'
System order:
%i\n
'
,
n
);
if
(
n
<
2
||
n
>
3
)
disp(
'
Cannot continue...
'
);
return
end
[
~
,
x
]
=
ode45(@
f
,
0
:
0.001
:
simtime
,
intial_values
);
%
ode45 at 1 khz
%%
Vector field
%
2nd order ODE
if
n
==
2
x1
=
linspace(range(
1
,
1
), range(
1
,
2
),
10
);
x2
=
linspace(range(
2
,
1
), range(
2
,
2
),
10
);
[
X1
,
X2
]
=
meshgrid(
x1
,
x2
);
u
=
zeros(size(
X1
));
v
=
zeros(size(
X2
));
t
=
0
;
for
i
=
1
:
numel(
X1
)
X_DOT
=
f(
t
,[X1(
i
); X2(
i
)]);
Vmod
=
sqrt(X_DOT(
1
)^
2
+
X_DOT(
2
)^
2
);
u(
i
)
=
X_DOT(
1
)/
Vmod
;
v(
i
)
=
X_DOT(
2
)/
Vmod
;
end
%
Drawing
h
=
quiver(
X1
,
X2
,
u
,
v
,
'
r
'
);
h.AutoScaleFactor
=
scale
;
axis
tight equal
;
hold
on
;
plot(x(
:
,
1
), x(
:
,
2
),
'
b
'
,
'
LineWidth
'
,
3
);
axis
tight equal
;
xlim(range(
1
,
:
));
ylim(range(
2
,
:
));
xlabel(
'
$x_1$
'
,
'
interpreter
'
,
'
latex
'
)
ylabel(
'
$x_2$
'
,
'
interpreter
'
,
'
latex
'
)
hold
off
;
%
3rd order ODE
elseif
n
==
3
x1
=
linspace(range(
1
,
1
), range(
1
,
2
),
10
);
x2
=
linspace(range(
2
,
1
), range(
2
,
2
),
10
);
x3
=
linspace(range(
3
,
1
), range(
3
,
2
),
10
);
[
X1
,
X2
,
X3
]
=
meshgrid(
x1
,
x2
,
x3
);
u
=
zeros(size(
X1
));
v
=
zeros(size(
X2
));
w
=
zeros(size(
X3
));
t
=
0
;
for
i
=
1
:
numel(
X1
)
X_DOT
=
f(
t
,[X1(
i
); X2(
i
); X3(
i
)]);
Vmod
=
sqrt(X_DOT(
1
)^
2
+
X_DOT(
2
)^
2
+
X_DOT(
3
)^
2
);
u(
i
)
=
X_DOT(
1
)/
Vmod
;
v(
i
)
=
X_DOT(
2
)/
Vmod
;
w(
i
)
=
X_DOT(
3
)/
Vmod
;
end
%
Drawing
h
=
quiver3(
X1
,
X2
,
X3
,
u
,
v
,
w
,
'
r
'
);
h.AutoScaleFactor
=
scale
;
hold
on
;
plot3(x(
:
,
1
), x(
:
,
2
), x(
:
,
3
),
'
b
'
,
'
LineWidth
'
,
3
);
axis
tight equal
;
xlim(range(
1
,
:
));
ylim(range(
2
,
:
));
zlim(range(
3
,
:
));
xlabel(
'
$x_1$
'
,
'
interpreter
'
,
'
latex
'
)
ylabel(
'
$x_2$
'
,
'
interpreter
'
,
'
latex
'
)
zlabel(
'
$x_3$
'
,
'
interpreter
'
,
'
latex
'
)
view([
1
1
1
]);
hold
off
;
end
fracparts
=
regexp(evalc(
'
sys
'
),
'
([^
\n
]*)
\n
[ ]*-[-]+[ ]*
\n
([^
\n
]*)
'
,
'
tokens
'
);
TFlatex
=
[
'
$$G(s)=
\f
rac{
'
fracparts{
1
}{
1
}
'
}{
'
fracparts{
1
}{
2
}
'
}$$
'
];
title(
TFlatex
,
'
interpreter
'
,
'
latex
'
);
%%
Callback function for the ODE solver
function
x_dot
=
f(
t
,
x
)
A
=
tf2ss(num{
:
}, den{
:
});
x_dot
=
A
*
x
;
end
end
Back
|
FazBrowse Home
|
New Git URL