Matlab - 線形微分方程式のオイラー法を記述して解く

Question

見出しからわかるように、私はODEについて何か質問したいと思います。私たちの先生はインターネット上にいくつかの例をアップロードしましたが、それらは単に解決策で構成されています。

ここに、次の方程式を解くことの1つがあります。

m_1 l^2 ϕ''_1 + m_1glϕ_1 + cl(ϕ_1 - ϕ_2) = 0 
m_2 l^2 ϕ''_2 + m_2glϕ_2 + cl(ϕ_2 - ϕ_1) = 0 

解決策がわかるように、D行列が定義されています。したがって、このD行列の最初の行は一般に1と0で構成され、残りはdy/dtを離した後に発生する定数で構成されます。このD行列は何に基づいていますか？私はそれをどのように定義しなければなりませんか？

%m1=100*m2 
% Parametres 
m_2=1;    
m_1=100*m_2;   
l=1;     
c=1; 
g=10; 
% Initial Conditions 
phi_1=1; 
phi_2=-0.5; 
dphi_1=0; 
dphi_2=0; 

% Iteration Parameters 
dt=0.01;    
t_end=10;   

% Time Vector 
t=0:dt:t_end;  
n_max=t_end/dt; 

y(:,1)=[phi_1;phi_2;dphi_1;dphi_2];            
D=[0 0 1 0;0 0 0 1;-g/l-c/(m_1*l) c/(m_1*l) 0 0;c/(m_2*l) -g/l-c/m_2/l 0 0];  

y_exp=y;  
for n = 1:n_max 
    y_exp(:,n+1)=dt*D*y_exp(:,n)+y_exp(:,n);   
end 

y_imp=y; 
for n = 1:n_max 
    y_imp(:,n+1)=inv(eye(4)-D*dt)*y_imp(:,n); 
end 

Answer 1

あなたのシステムは、数値積分の目的のために、オーダー1のシステム及び寸法4を作るため2次元2の一定の係数を有する線形であります変数

y = [ phi_1; phi_2; dphi_1; dphi_2 ] 

マトリックスDは、その注文システム1のヤコビアンです。拡張式

phi_1' = dphi_1 
phi_2' = dphi_2 
dphi_1' = - (m_1*g*l * phi_1 + c*l * (phi_1 - phi_2))/(m_1 * l^2) 
dphi_2' = - (m_2*g*l * phi_2 + c*l * (phi_2 - phi_1))/(m_2 * l^2) 

以下

y' = D * y 

として製剤化マトリックス