Unable to perform assignment because the left and right sides have a different number of elements
2 Ansichten (letzte 30 Tage)
Ältere Kommentare anzeigen
% 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 Kommentar
KSSV
am 10 Mär. 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.
Antworten (1)
Torsten
am 10 Mär. 2022
Bearbeitet: Torsten
am 15 Mär. 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 Kommentare
Torsten
am 15 Mär. 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
am 15 Mär. 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)
Siehe auch
Kategorien
Mehr zu Combustion and Turbomachinery finden Sie in Help Center und File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!