現在この質問をフォロー中です
- フォローしているコンテンツ フィードに更新が表示されます。
- コミュニケーション基本設定に応じて電子メールを受け取ることができます。
State space modeling of LTI
3 ビュー (過去 30 日間)
古いコメントを表示
I have a LTI system is given in the attached picture. For which I have A=[0 1;1 2] B=[1;1], C=[3,4], D=0.1 and L=[0 1]. I need to find z(t) over an interval [0,20] and initial condition is assumed as x0=[1;-1] with an input w=0.5sin(2pi)t and v=t, how can I make a code for this??
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/282863/image.jpeg)
採用された回答
Ameer Hamza
2020 年 4 月 8 日
編集済み: Ameer Hamza
2020 年 4 月 8 日
Such LTI system can be solved using lsim
C=[3,4];
D=0.1;
L=[0 1];
sys = ss(A,B,C,D);
t = linspace(0,20,1000)';
x0 = [1; -1];
[t,x] = ode45(@odeFun, t, x0);
v = t;
y = C*x' + D*v';
z = L*x';
plot(t,z)
function dxdt = odeFun(t, x)
A=[0 1;1 2];
B=[1;1];
w = 0.5*sin(2*pi*t);
dxdt = A*x+B*w;
end
However, the system is unstable and the output diverges to infinity.
19 件のコメント
Muhammad Atif
2020 年 4 月 8 日
thanks sir, but how can I use initial condition here and what's about disturbance 'v'??
Ameer Hamza
2020 年 4 月 8 日
編集済み: Ameer Hamza
2020 年 4 月 8 日
Ok. I just noticed that your LTI system is not of a standard form. So regular tools for LTI systems will not work here. I gave a solution by solving the differential equations. Please check the updated answer. It also considers disturbance v and the initial conditions.
Muhammad Atif
2020 年 4 月 11 日
Hi Sir, I modified the code that you created above, In this I added the input w(t) and noise signal v(t). The plot of (z=L*v' ) should start when input starts to take effects. I attached two images for references. This plot I obtained from this code, while I wanna start it from t=5. How can I change this?![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/283588/image.jpeg)
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/283588/image.jpeg)
%%function to get dxdt
function dxdt = odeFun(t, x)
ms = 973;
ks = 42720;
cs = 3000;
ku = 101115;
mu = 114;
G = 512*10^(-6)
q = 0.1;
v1 = 12.5;
A = [0 0 1 -1 ; 0 0 0 1 ; -ks/ms 0 -cs/ms cs/ms; ks/mu -ku/mu cs/mu -cs/mu];
B = [0 0 -2*pi*q*(G*v1)^(1/2) 0]';
a = 0.2 ;
l = 5 ;
%%input w(t)=1 for 5<=t<=10 and w(t)=-1 for 15<=t<=20
w_int1 = t>=5 & t<=10;
w_int2 = t>=15 & t<=20;
w = zeros(1,numel(t));%add this
rng(0);
w(w_int1)=1;
w(w_int2)=-1;
dxdt = A*x+B*w';
end
C = [0 0 0 1];
D = 0.1;
L = [0 0 1 0];
t = linspace(0,30,1000)';
x0 = [0.1; 0; -0.1; -0.2];
[t,x] = ode45(@odeFun, t, x0);
T = linspace(0,20,1000)';
%%noise v(t) distributed over an interval [0,20]
v_int = t>=0 & t<=20;
v = zeros(1,numel(t));
rng(0);
v(v_int)= 0.2*rand(1,length(t(v_int)));
y = C*x' + D*v;
z = L*x';
plot(t,z)
Ameer Hamza
2020 年 4 月 11 日
Start from initial condition of zero to see the effect. Alsi it is better to use if-block to calculate w based on t
C = [0 0 0 1];
D = 0.1;
L = [0 0 1 0];
t = linspace(0,30,1000)';
x0 = [0; 0; -0; -0];
[t,x] = ode45(@odeFun, t, x0);
T = linspace(0,20,1000)';
%%noise v(t) distributed over an interval [0,20]
% v_int = t>=0 & t<=20;
% v = zeros(1,numel(t));
% rng(0);
% v(v_int)= 0.2*rand(1,length(t(v_int)));
% y = C*x' + D*v;
% z = L*x';
plot(t,x)
%%function to get dxdt
function dxdt = odeFun(t, x)
ms = 973;
ks = 42720;
cs = 3000;
ku = 101115;
mu = 114;
G = 512*10^(-6);
q = 0.1;
v1 = 12.5;
A = [0 0 1 -1 ; 0 0 0 1 ; -ks/ms 0 -cs/ms cs/ms; ks/mu -ku/mu cs/mu -cs/mu];
B = [0 0 -2*pi*q*(G*v1)^(1/2) 0]';
a = 0.2 ;
l = 5 ;
%%input w(t)=1 for 5<=t<=10 and w(t)=-1 for 15<=t<=20
w = 0;
if t>=5 && t<=10
w = 1;
elseif t>=15 && t<=20
w = -1;
end
dxdt = A*x+B*w';
end
Muhammad Atif
2020 年 4 月 13 日
Hi sir, As in the perviously question we found z(t) and y(t); Now in this question I need to use the y(t) as the in put of this state-space equation for the given Af = [-2.8682 0.5240 0.1093 0.4576; 0.6676 -2.3337 -0.1285 -0.5333; 4.2978 1.3131 -13.1636 3.3939; -30.1511 55.6790 48.9362 -218.6780];
Bf = [0.0167; -0.0400; -0.1676; -3.6279];
Cf = [1.6042 1.1359 -0.7995 -0.1854];
and xf=[0;0;0;0] is the intial state for this system
and the state space equation is given in the picture
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/284078/image.jpeg)
Ameer Hamza
2020 年 4 月 13 日
You can use lsim for this system
Af = [-2.8682 0.5240 0.1093 0.4576; 0.6676 -2.3337 -0.1285 -0.5333; 4.2978 1.3131 -13.1636 3.3939; -30.1511 55.6790 48.9362 -218.6780];
Bf = [0.0167; -0.0400; -0.1676; -3.6279];
Cf = [1.6042 1.1359 -0.7995 -0.1854];
Df = 0;
sys = ss(Af,Bf,Cf,Df);
lsim(sys,y_hat,t); % t should be a column vector and y_hat should have same number of rows as t and 2 columns
Muhammad Atif
2020 年 4 月 13 日
Af = [-2.8682 0.5240 0.1093 0.4576; 0.6676 -2.3337 -0.1285 -0.5333; 4.2978 1.3131 -13.1636 3.3939; -30.1511 55.6790 48.9362 -218.6780];
Bf = [0; -0; -0; 0];
Cf = [1.6042 1.1359 -0.7995 -0.1854];
Df = 0;
sys = ss(Af,Bf,Cf,Df);
y_hat = zeros(size(t));
lsim(sys,y_hat,t)
why this system gives zero output when I choose y_hat=0??
Ameer Hamza
2020 年 4 月 13 日
Because the input to the system is zero so the states cannot change. You need to give some non-zero input to change the state.
Muhammad Atif
2020 年 4 月 18 日
Sorry to bother you again. Now I have no idea how to use the the state of xf0 for the other system. It's a sort of switching condition in which two subsystems. When 't' reaches the defined interval the system one will run, otherwise the 2nd system will run. But how can I use the state xf0 of the first system for system 2 when system 1 run for the last time??
clc
clear all
Af1 = [-2.8682 0.5240 0.1093 0.4576; 0.6676 -2.3337 -0.1285 -0.5333; 4.2978 1.3131 -13.1636 3.3939; -30.1511 55.6790 48.9362 -218.6780];
Bf1 = [0.0167; -0.0400; -0.1676; -3.6279];
Bf2 = [0;0;0;0];
Af2 = [-2.8682 0.5240 0.1093 0.4576; 0.6676 -2.3337 -0.1285 -0.5333; 4.2978 1.3131 -13.1636 3.3939; -30.1511 55.6790 48.9362 -218.6780];
Cf1 = [1.6042 1.1359 -0.7995 -0.1854];
Cf2 = [0.9791 0.7264 -0.6448 -0.0632];
Df1=0;
Df2=0;
C = [0 0 0 1];
D = 0.1;
L = [0 0 1 0];
t = linspace(0,30,1000);
x0 = [0; 0; -0.00; -0.00];
[t,x] = ode45(@odeFun12, t, x0);
%%noise v(t) distributed over an interval [0,20]
v = 0.05*rand(size(t));
v = v.*(t<20);
y = C*x' + D*v';
z = L*x';
%plot(t,y)
%grid on
%%system 1
if t<=1 | (t>=2 & t<=7) | (t>=12 & t<=17) | (t>=22 & t<=27)
sys = ss(Af1,Bf1,Cf1,Df1);
xf0 = [0; 0; 0.00; 0.00];
lsim(sys,y,t,xf0);
%%system 2 I need to use xf0 of the 1st system when program switches to system 2
else
sys1 = ss(Af2,Bf2,Cf2,Df2);
y1=zeros(size(t));
lsim(sys1,y1,t,xf0);
end
Muhammad Atif
2020 年 4 月 18 日
function dxdt = odeFun12(t, x)
ms = 973;
ks = 42720;
cs = 3000;
ku = 101115;
mu = 114;
G = 512*10^(-6);
q = 0.1;
v1 = 12.5;
A = [0 0 1 -1; 0 0 0 1; -ks/ms 0 -cs/ms cs/ms; ks/mu -ku/mu cs/mu -cs/mu];
B = [0 -2*pi*q*(G*v1)^(1/2) 0 0]';
a = 0.2 ;
l = 25;
%%input w = a*pi*v/l.*sin(2*pi*v/l.*t) for 5<=t<=l/v and zero otherwise, while 'v' is uniformely distributed noise between[0,0.1] for an interval 0 to 20
w = zeros(size(t));
mask = t>=5 & t<=5+l/v1;
w = a*pi*v1/l.*sin(2*pi*v1*t/l).*mask;
dxdt = A*x+B*w';
end
Ameer Hamza
2020 年 4 月 18 日
I guess something like this will work
Cf1 = [1.6042 1.1359 -0.7995 -0.1854];
Cf2 = [0.9791 0.7264 -0.6448 -0.0632];
Df1=0;
Df2=0;
t = linspace(0,30,1000);
x0 = [0; 0; -0.00; -0.00];
[t,x] = ode45(@odeFun12, t, x0);
%%noise v(t) distributed over an interval [0,20]
mask = t<=1 | (t>=2 & t<=7) | (t>=12 & t<=17) | (t>=22 & t<=27);
y(mask) = Cf1*x(mask,:).';
y(~mask) = Cf2*x(~mask,:).';
function dxdt = odeFun12(t, x)
Af1 = [-2.8682 0.5240 0.1093 0.4576; 0.6676 -2.3337 -0.1285 -0.5333; 4.2978 1.3131 -13.1636 3.3939; -30.1511 55.6790 48.9362 -218.6780];
Bf1 = [0.0167; -0.0400; -0.1676; -3.6279];
Bf2 = [0;0;0;0];
Af2 = [-2.8682 0.5240 0.1093 0.4576; 0.6676 -2.3337 -0.1285 -0.5333; 4.2978 1.3131 -13.1636 3.3939; -30.1511 55.6790 48.9362 -218.6780];
ms = 973;
ks = 42720;
cs = 3000;
ku = 101115;
mu = 114;
G = 512*10^(-6);
q = 0.1;
v1 = 12.5;
a = 0.2 ;
l = 25;
%%input w = a*pi*v/l.*sin(2*pi*v/l.*t) for 5<=t<=l/v and zero otherwise, while 'v' is uniformely distributed noise between[0,0.1] for an interval 0 to 20
w = a*pi*v1/l.*sin(2*pi*v1*t/l);
if t<=1 || (t>=2 && t<=7) || (t>=12 && t<=17) || (t>=22 && t<=27)
dxdt = Af1*x+Bf1*w';
else
dxdt = Af2*x+Bf2*w';
end
end
Muhammad Atif
2020 年 4 月 19 日
I think I didn't explain well and you couldn't get me. y(t) is already calculated from the ode fuction. In the last question I need to use y(t) as an input for the next two systems. In which y(t) for the 2nd system is zero, while the first system will use the y(t) calculated before. lsim(sys,y,t,xf0),as you see 'y' is used as an input for this system. y1=zeros(size(t)),lsim(sys1,y1,t,xf0); in this it is zero, but this system will use the last state of first system when it will run
Ameer Hamza
2020 年 4 月 19 日
I think something like this will work. You can extend it to the full time duration.
clc
clear all
Af1 = [-2.8682 0.5240 0.1093 0.4576; 0.6676 -2.3337 -0.1285 -0.5333; 4.2978 1.3131 -13.1636 3.3939; -30.1511 55.6790 48.9362 -218.6780];
Bf1 = [0.0167; -0.0400; -0.1676; -3.6279];
Bf2 = [0;0;0;0];
Af2 = [-2.8682 0.5240 0.1093 0.4576; 0.6676 -2.3337 -0.1285 -0.5333; 4.2978 1.3131 -13.1636 3.3939; -30.1511 55.6790 48.9362 -218.6780];
Cf1 = [1.6042 1.1359 -0.7995 -0.1854];
Cf2 = [0.9791 0.7264 -0.6448 -0.0632];
Df1=0;
Df2=0;
C = [0 0 0 1];
D = 0.1;
L = [0 0 1 0];
t = linspace(0,30,1000);
x0 = [0; 0; -0.00; -0.00];
[t,x] = ode45(@odeFun12, t, x0);
%%noise v(t) distributed over an interval [0,20]
v = 0.05*rand(size(t));
v = v.*(t<20);
y = C*x' + D*v';
z = L*x';
%plot(t,y)
%grid on
%%system 1
% if t<=1 | (t>=2 & t<=7) | (t>=12 & t<=17) | (t>=22 & t<=27)
sys1 = ss(Af1,Bf1,Cf1,Df1);
sys2 = ss(Af2,Bf2,Cf2,Df2);
xf0 = [0; 0; 0.00; 0.00];
t1 = t(t<=1);
y1 = y(t<=1);
[~,Y1,X1] = lsim(sys1,y1,t1,xf0);
t2 = t(1<t & t<=2);
y2 = 0*y(1<t & t<=2);
[~,Y2,X2] = lsim(sys2,y2,t2,X1(end,:));
t3 = t(2<t & t<=7);
y3 = y(2<t & t<=7);
[~,Y3,X3] = lsim(sys1,y3,t3,X2(end,:));
t4 = t(7<t & t<=12);
y4 = 0*y(7<t & t<=12);
[~,Y4,X4] = lsim(sys2,y4,t4,X3(end,:));
Y = [Y1;Y2;Y3;Y4];
X = [X1;X2;X3;X4];
function dxdt = odeFun12(t, x)
ms = 973;
ks = 42720;
cs = 3000;
ku = 101115;
mu = 114;
G = 512*10^(-6);
q = 0.1;
v1 = 12.5;
A = [0 0 1 -1; 0 0 0 1; -ks/ms 0 -cs/ms cs/ms; ks/mu -ku/mu cs/mu -cs/mu];
B = [0 -2*pi*q*(G*v1)^(1/2) 0 0]';
a = 0.2 ;
l = 25;
%%input w = a*pi*v/l.*sin(2*pi*v/l.*t) for 5<=t<=l/v and zero otherwise, while 'v' is uniformely distributed noise between[0,0.1] for an interval 0 to 20
w = zeros(size(t));
mask = t>=5 & t<=5+l/v1;
w = a*pi*v1/l.*sin(2*pi*v1*t/l).*mask;
dxdt = A*x+B*w';
end
Muhammad Atif
2020 年 4 月 20 日
That's it, finally I got the results. Thank you very much for your help
その他の回答 (0 件)
参考
タグ
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!エラーが発生しました
ページに変更が加えられたため、アクションを完了できません。ページを再度読み込み、更新された状態を確認してください。
Web サイトの選択
Web サイトを選択すると、翻訳されたコンテンツにアクセスし、地域のイベントやサービスを確認できます。現在の位置情報に基づき、次のサイトの選択を推奨します:
また、以下のリストから Web サイトを選択することもできます。
最適なサイトパフォーマンスの取得方法
中国のサイト (中国語または英語) を選択することで、最適なサイトパフォーマンスが得られます。その他の国の MathWorks のサイトは、お客様の地域からのアクセスが最適化されていません。
南北アメリカ
- América Latina (Español)
- Canada (English)
- United States (English)
ヨーロッパ
- Belgium (English)
- Denmark (English)
- Deutschland (Deutsch)
- España (Español)
- Finland (English)
- France (Français)
- Ireland (English)
- Italia (Italiano)
- Luxembourg (English)
- Netherlands (English)
- Norway (English)
- Österreich (Deutsch)
- Portugal (English)
- Sweden (English)
- Switzerland
- United Kingdom (English)
アジア太平洋地域
- Australia (English)
- India (English)
- New Zealand (English)
- 中国
- 日本Japanese (日本語)
- 한국Korean (한국어)