Problem facing in solving xdot=Ax system when A is unstable.
11 vues (au cours des 30 derniers jours)
Afficher commentaires plus anciens
When solving a problem using ode45 for the system
A is special unstable matrix. I encountered a situation where the solution x turned out to be on the order of $10^29$. Since I'm interested in the error between two states, theoretically the error should be zero. However, up to 4 decimal places, the values of the states seems to be the same, i.e., for example,
and
and but e=x1−x2 is not equal to zero. I understand that x1 and x2 have mismatched in decimal points, and the values are amplified by the exponent. How can I restrict these numerical instabilities? Please help.
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/1628198/image.png)
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/1628213/image.png)
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/1628218/image.png)
3 commentaires
Réponse acceptée
Sam Chak
le 26 Fév 2024
Modifié(e) : Sam Chak
le 26 Fév 2024
This is an educated guess. If the unstable state matrix is
and the initial values are
, then the error between the two states is perfectly zero.
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/1628328/image.png)
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/1628333/image.png)
%% Analytical approach
syms x(t) y(t)
eqns = [diff(x,t) == y,
diff(y,t) == x];
cond = [x(0)==1, y(0)==1];
S = dsolve(eqns, cond)
%% Numerical approach
[t, x] = ode45(@odefun, [0 10], [1 1]);
plot(t, x), grid on,
xlabel('t'), ylabel('\bf{x}')
legend({'$x(t)$', '$\dot{x}(t)$'}, 'interpreter', 'latex', 'fontsize', 16, 'location', 'NW')
%% error between two states
error = x(:,1) - x(:,2);
norm(error)
function dxdt = odefun(t, x)
A = [0 1;
1 0];
dxdt = A*x;
end
3 commentaires
Sam Chak
le 27 Fév 2024
Without the analytical solution, how can you be certain that states
and
are precisely identical, as depicted in the example I provided in my previous answer? Furthermore, considering the random initial values, if
differs from
, it is impossible for
to be equivalent to
.
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/1628863/image.png)
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/1628868/image.png)
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/1628873/image.png)
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/1628878/image.png)
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/1628883/image.png)
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/1628888/image.png)
N = 3;
A = [0 1 0;0 0 1;3 5 6];
% [v, d] = eig(A)
B = [0; 0; 1];
D = [1 -1 0;-1 2 -1;0 -1 1];
K = [107, 71, 20];
A1 = kron(eye(N,N),A);
A2 = kron(D,B*K);
Aa = A1 - A2
eig(Aa);
x_0 = rand(3*N, 1);
%x_0=rand(6,1);
tspan = 0:1e-6:1.;
% opts = odeset('AbsTol',1e-8,'RelTol',1e-4, 'MaxStep',0.01);
% e1=eig(A-1*B*K)
% e2=eig(A-3*B*K)
[t, x] = ode45(@(t,x) odefcn(t,x,D,A,B,K,N), tspan, x_0);
plot(t, x(:,1), t, x(:,4)), grid on
legend('x_1', 'x_4', 'fontsize', 16, 'location', 'NW')
e112 = x(:,1) - x(:,4);
norm(e112)
%% Error between x1 and x4
plot(t, e112), grid on, title('Error')
function xdot=odefcn(t,x,D,A,B,K,N)
A1 = kron(eye(N,N),A);
A2 = kron(D,B*K);
xdot= (A1 - A2)*x;
end
Plus de réponses (0)
Voir également
Catégories
En savoir plus sur Matrix Indexing dans Help Center et File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!