Unable to perform assignment because the left and right sides have a different number of elements
2 vues (au cours des 30 derniers jours)
Afficher commentaires plus anciens
% Input Data Parameters
h=1; % Step size [s]
tf=50; % Solve Space of t=[0,tf] [s]
Ta=297;
P=[1197, 1.8, 15.8];
I = 2;
C=4;
% Simulation time interval
t=0:h:tf;
% Coefficient of irreversible and reversible terms
A=[-1298*10^(-3),4038*10^(-3),-4674*10^(-3),2500*10^(-3),-618.6*10^(-3),128.5*10^(-3)];
B=[8.814*10^(-3),23.01*10^(-3),-19.09*10^(-3),4.099*10^(-3),1.011*10^(-3),-0.1937*10^(-3)];
% Intial conditions
Ts=zeros(1,length(t));
Ts(1)=Ta;
SOC=zeros(1,length(t));
Relec = zeros(1,length(t));
dUOC = zeros(1,length(t));
Qgen = zeros(1,length(t));
SOC(1)=1;
% ODE equation to solve
f=@(t,Ts) (Ta-Ts)/(P(1).*(P(2)+P(3)))+(Qgen.*P(3))/(P(1).*(P(2)+P(3)));
% RK4 Loop
for i=1:(length(t)-1)
SOC(i+1) = SOC(i)-(1/C)*I*((t(i+1)-t(i))/3600);
Relec = A(1).*(SOC.^5)+A(2).*(SOC.^4)+A(3).*(SOC.^3)+A(4).*(SOC.^2)+A(5).*(SOC)+A(6);
dUOC = B(1).*(SOC.^5)+B(2).*(SOC.^4)+B(3).*(SOC.^3)+B(4).*(SOC.^2)+B(5).*(SOC)+B(6);
Qgen =(I^2) .* Relec-I.*dUOC;
k1=f(t(i),Ts(i));
k2=f(t(i)+0.5*h,Ts(i)+0.5*k1*h);
k3=f(t(i)+0.5*h,Ts(i)+0.5*k2*h);
k4=f(t(i)+h,Ts(i)+k3*h);
Ts(i+1)=Ts(i)+(h/6)*(k1+2*k2+2*k3+k4);
1 commentaire
KSSV
le 10 Mar 2022
You need to rethink on your code. This function:
f=@(t,Ts) (Ta-Ts)/(P(1).*(P(2)+P(3)))+(Qgen.*P(3))/(P(1).*(P(2)+P(3)));
Gives you a array as output. It should give you a scalar as output. Check this function.
Réponses (1)
Torsten
le 10 Mar 2022
Modifié(e) : Torsten
le 15 Mar 2022
% Input Data Parameters
h=1; % Step size [s]
tf=50; % Solve Space of t=[0,tf] [s]
Ta=297;
P=[1197, 1.8, 15.8];
I = 2;
C=4;
% Simulation time interval
t=0:h:tf;
% Coefficient of irreversible and reversible terms
A=[-1298*10^(-3),4038*10^(-3),-4674*10^(-3),2500*10^(-3),-618.6*10^(-3),128.5*10^(-3)];
B=[8.814*10^(-3),23.01*10^(-3),-19.09*10^(-3),4.099*10^(-3),1.011*10^(-3),-0.1937*10^(-3)];
% Intial conditions
Ts=zeros(1,length(t));
Ts(1)=Ta;
soc = zeros(1,length(t));
Relec = zeros(1,length(t)-1);
dUOC = zeros(1,length(t)-1);
Qgen = zeros(1,length(t)-1);
soc(1)=1;
% ODE equation to solve
f=@(t,T,Q) (Ta-T)/(P(1).*(P(2)+P(3)))+(Q.*P(3))/(P(1).*(P(2)+P(3)));
for i=1:(length(t)-1)
SOC = soc(i);
Relec(i) = A(1).*(SOC.^5)+A(2).*(SOC.^4)+A(3).*(SOC.^3)+A(4).*(SOC.^2)+A(5).*(SOC)+A(6);
dUOC(i) = B(1).*(SOC.^5)+B(2).*(SOC.^4)+B(3).*(SOC.^3)+B(4).*(SOC.^2)+B(5).*(SOC)+B(6);
Qgen(i) =(I^2) .* Relec(i)-I.*dUOC(i);
k1=f(t(i),Ts(i),Qgen(i));
k2=f(t(i)+0.5*h,Ts(i)+0.5*k1*h,Qgen(i));
k3=f(t(i)+0.5*h,Ts(i)+0.5*k2*h,Qgen(i));
k4=f(t(i)+h,Ts(i)+k3*h,Qgen(i));
Ts(i+1)=Ts(i)+(h/6)*(k1+2*k2+2*k3+k4);
soc(i+1) = soc(i)-(1/C)*I*((t(i+1)-t(i))/3600);
end
plot(t,Ts)
end
Usually, you would also have to modify Qgen for each call of the function f depending on time.
For simplicity, I assumed Qgen to be constant on t(i) <= t <= t(i+1)
3 commentaires
Torsten
le 15 Mar 2022
Qgen(i) =(I^2) .* Relec(i)-I.*dUOC(i);
If Qgen is heat generation, you can set this term (which is heat generation at time t(i), I guess) to a value of your choice.
Jan
le 15 Mar 2022
As usual I mention, that -4674*10^(-3) is an expensive power operation, while -4674e-3 is a cheap constant. I've learned programming on a ZX81 with 1kB RAM and avoiding EXP is my destiny.
A(1).*(SOC.^5)+A(2).*(SOC.^4)+A(3).*(SOC.^3)+A(4).*(SOC.^2)+A(5).*(SOC)+A(6)
% Horner scheme:
SOC * (SOC * (SOC * (SOC * (A(1) * (SOC) + A(2)) + A(3)) + A(4)) + A(5)) + A(6)
Voir également
Catégories
En savoir plus sur Combustion and Turbomachinery 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!