solve multipe eq in a loop
显示 更早的评论
hello
would you plz fix this code
i dont know how to solve multiple variables in a loop
clc
clear all
dt=1 ;V=1000*1000*100 ; pro=0.15 ; Co=10^-6 ; alpha =5.615; Bo=1.35;
beta =1.127 ;Ax=1000*100 ; mio=1 ;Dx=1000;Kx=20*10^-3;
cte=(V*pro*Co)/(alpha*Bo*dt);
T=(beta * Ax* Kx)/(mio *Bo*Dx);
p1(1:360)=5000 ; p2(1)=5000; p3(1)=5000;p4(1)=5000 ; p5(1)=5000; p6(1)=5000;
p7(1)=5000 ; p8(1)=5000; p9(1)=5000;
syms p2(i+1) p3(i+1) p4(i+1) p5(i+1) p6(i+1) p7(i+1) p8(i+1) p9(i+1)
for i=1:1:359
eqn = [
p3(i+1)*T-(cte +2*T)*p2(i+1)+T*p1(i+1)==-(cte *p2(i));
p4(i+1)*T-(cte +2*T)*p3(i+1)+T*p2(i+1)==-(cte *p2(i)-250);
p5(i+1)*T-(cte +2*T)*p4(i+1)+T*p3(i+1)==-(cte *p3(i));
p6(i+1)*T-(cte +2*T)*p5(i+1)+T*p4(i+1)==-(cte *p4(i));
p7(i+1)*T-(cte +2*T)*p6(i+1)+T*p5(i+1)==-(cte *p5(i));
p8(i+1)*T-(cte +2*T)*p7(i+1)+T*p6(i+1)==-(cte *p6(i)-250);
p9(i+1)*T-(cte +2*T)*p8(i+1)+T*p7(i+1)==-(cte *p7(i));
-(cte +2*T)*p9(i+1)+T*p8(i+1)==-(cte *p8(i))
];
S = solve(eqn, [p2(i+1), p3(i+1), p4(i+1), p5(i+1), p6(i+1), p7(i+1),
p8(i+1),p9(i+1)]);
p2(i+1) = S.p2(i+1);
p3(i+1) = S.p3(i+1);
p4(i+1) = S.p4(i+1);
p5(i+1) = S.p5(i+1);
p6(i+1) = S.p6(i+1);
p7(i+1) = S.p7(i+1);
p8(i+1) = S.p8(i+1);
p9(i+1) = S.p9(i+1);
end
time_d = {'90';'180';'270';'360'};
p1= [p1(90);p1(180);p1(270);p1(360)];
p2 = [p2(90);p2(180);p2(270);p2(360)];
p3 = [p3(90);p3(180);p3(270);p3(360)];
p4 = [p4(90);p4(180);p4(270);p4(360)];
p5 = [p5(90);p5(180);p5(270);p5(360)];
p6= [p6(90);p6(180);p6(270);p6(360)];
p7= [p7(90);p7(180);p7(270);p7(360)];
p8= [p8(90);p8(180);p8(270);p8(360)];
p9= [p9(90);p9(180);p9(270);p9(360)];
T1 = table(time_d,p1,p2,p3,p4,p5,p6,p7,p8,p9)
1 个评论
KALYAN ACHARJYA
2021-1-1
Please read the following thread
https://in.mathworks.com/matlabcentral/answers/304528-tutorial-why-variables-should-not-be-named-dynamically-eval
回答（1 个）
Alan Stevens
2021-1-1
Your equations are linear in p (p2 to p9, that is, since p1 is constant for all times), so the following is a simple way to solve for them:
dt=1 ;V=1000*1000*100 ; pro=0.15 ; Co=10^-6 ; alpha =5.615; Bo=1.35;
beta =1.127 ;Ax=1000*100 ; mio=1 ;Dx=1000;Kx=20*10^-3;
cte=(V*pro*Co)/(alpha*Bo*dt);
T=(beta * Ax* Kx)/(mio *Bo*Dx);
N = 360;
p = zeros(8,N);
p1=5000 ; p(:,1) = 5000;
k = -(cte + 2*T);
% Equations are linear in p, so solve using M*p(i) = -cte*p(i-1) - K
% or p(i) = M\(-cte*p(i-1)-K), where M and K are given below.
M = [k, T, 0, 0, 0, 0, 0, 0;
T, k, T, 0, 0, 0, 0, 0;
0, T, k, T, 0, 0, 0, 0;
0, 0, T, k, T, 0, 0, 0;
0, 0, 0, T, k, T, 0, 0;
0, 0, 0, 0, T, k, T, 0;
0, 0, 0, 0, 0, T, k, T;
0, 0, 0, 0, 0, 0, T, k];
K = [T*p1; 250; 0; 0; 0; 250; 0; 0];
for i = 2:N
p(:,i) = M\(-cte*p(:,i-1) - K);
end
t = 1:N;
plot(t,p),grid
xlabel('t'),ylabel('p2 to p9')
legend('p2','p3','p4','p5','p6','p7','p8','p9')
Note that the eight rows of p represent p2 to p9.
类别
在 帮助中心 和 File Exchange 中查找有关 General Applications 的更多信息
2021-1-1
2021-1-1
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
Translated by 
