Can someone helps me to solve the error in this code?
1 次查看(过去 30 天)
显示 更早的评论
il = 4;
R = 8.314;
F = 96485;
Alpha = 0.5;
n = 2;
Eexc = 53990.065;
sigma_ref = 0.020;
Epro = 18912.42;
delta = 7*2.54*10^-3;
I = 0:50:350;
for i = 1:length(I)
t=0;
for T = [313 333 353]
t=t+1;
% Anode
ia = (10^-7)(exp(-(Eexc/R)((1/T)-(1/298))));
Eta_Anode(i,t) = ((R*T)/(Alpha*n*F))*(log(I(i)/ia));
% Cathode
ic = (10^-3)(exp(-(Eexc/R)((1/T)-(1/353))));
Eta_Cathode(i,t) = -((R*T)/(Alpha*n*F))*(log(I(i)/ic));
%Ohmic
sigma(i,T) = sigma_ref*(exp(-(Epro/R)*((1/T)-(1/298))));
Rohm = delta/sigma(i,T);
Eta_Ohmic(i,t) = I(i)*Rohm;
%Concentration
Eta_Conc(i,t) = ((R*T)/(Alpha*n*F))*(log(il/(il-I(i))));
%Voltage
H_h2(i,T)=0+(28.84*(T-298));
H_o2(i,T)=0+(28.91*(T-298));
H_h2o(i,T)=-285830+(75.37*(T-298));
H(i,T)=H_h2(i,T)+(0.5*H_o2(i,T))-H_h2o(i,T);
S_h2(i,T)=130.68+(28.84*(log(T/298)));
S_o2(i,T)=205+(28.91*(log(T/298)));
S_h2o(i,T)=69.95+(75.37*(log(T/298)));
S(i,T)=S_h2(i,T)+(0.5*S_o2(i,T))-S_h2o(i,T);
G(i,T)=H(i,T)-(T*S(i,T));
E(i,T)=(G(i,T)/(2*F));
ET(i,t) = 1.5184-(1.5421*10^-3*T)+(9.523*10^-5*T*log(T))+(9.84*10^-8*T^2);
EP(i,t)= 1.23+(((R*T)/(n*F))*(log(0.5*0.5^0.5)));
VT(i,t) = ET(i,t) + Eta_Anode(i,t) - Eta_Cathode(i,t) + Eta_Ohmic(i,t) + Eta_Conc(i,t);
VP(i,t) = EP(i,t) + Eta_Anode(i,t) - Eta_Cathode(i,t) + Eta_Ohmic(i,t) + Eta_Conc(i,t);
V(i,t) = E(i,t) + Eta_Anode(i,t) - Eta_Cathode(i,t) + Eta_Ohmic(i,t) + Eta_Conc(i,t);
%Mass Flow Rate
A(i) = I(i)* 50;
M_H2(i) = (A(i)*0.99)/(2*F);
M_O2(i) = (A(i)*0.99)/(4*F);
M_H2O(i) = (A(i)*0.99)/(2*F);
disp (['At Current = ',num2str(I(i)),'A',', Mass Flow Rate of H2 = ', num2str(M_H2(i)),' mol/s',', Mass Flow Rate of O2 = ', num2str(M_O2(i)),' mol/s',', Mass Flow Rate of H2O = ', num2str(M_H2O(i)),' mol/s']);
end
end
figure(1);
plot (I,V(:,1),I,V(:,2),I,V(:,3),I,V(:,4),I,V(:,5),I,V(:,6))
legend('T=25','T=40','T=56','T=64','T=55','T=48')
xlabel('Current Density (A/cm2)');
ylabel('Cell Voltage (Volt)');
title('I-V Curve (General Equation of E reversible)');
figure(2);
plot (I,VP(:,1),I,VP(:,2),I,VP(:,3),I,VP(:,4),I,VP(:,5),I,VP(:,6))
legend('T=25','T=40','T=56','T=64','T=55','T=48')
xlabel('Current Density (A/cm2)');
ylabel('Cell Voltage (Volt)');
title('I-V Curve (Pressure & Temperature variation)');
figure(3);
plot (I,VT(:,1),I,VT(:,2),I,VT(:,3),I,VT(:,4),I,VT(:,5),I,VT(:,6))
legend('T=25','T=40','T=56','T=64','T=55','T=48')
xlabel('Current Density (A/cm2)');
ylabel('Cell Voltage (Volt)');
title('I-V Curve (Temperature variation)');
figure(4);
plot(I,M_H2(:),I,M_H2(:),I,M_H2(:),I,M_H2(:),I,M_H2(:),I,M_H2(:))
xlabel('Current Density (A/cm2)');
ylabel('Molar Flow Rate (mol/s)');
title('Current vs Molar Flow Rate of H2');
figure(5);
plot(I,Eta_Anode(:,1),I,Eta_Anode(:,2),I,Eta_Anode(:,3),I,Eta_Anode(:,4),I,Eta_Anode(:,5),I,Eta_Anode(:,6),I,Eta_Cathode(:,1),I,Eta_Cathode(:,2),I,Eta_Cathode(:,3),I,Eta_Cathode(:,4),I,Eta_Cathode(:,5),I,Eta_Cathode(:,6),I,Eta_Ohmic(:,1),I,Eta_Ohmic(:,2),I,Eta_Ohmic(:,3),I,Eta_Ohmic(:,4),I,Eta_Ohmic(:,5),I,Eta_Ohmic(:,6),I,Eta_Conc(:,1),I,Eta_Conc(:,2),I,Eta_Conc(:,3),I,Eta_Conc(:,4),I,Eta_Conc(:,5),I,Eta_Conc(:,6))
legend('T=25','T=40','T=56','T=64','T=55','T=48')
xlabel('Current Density (A/cm2)');
ylabel('Overpotential (Volt)');
title('Different Overpotential');
figure(6);
plot(I,Eta_Anode(:,1),I,Eta_Anode(:,2),I,Eta_Anode(:,3),I,Eta_Anode(:,4),I,Eta_Anode(:,5),I,Eta_Anode(:,6))
xlabel('Current Density (A/cm2)');
ylabel('Overpotential (Volt)');
title('Activation Anode Overpotential');
figure(7); plot(I,Eta_Cathode(:,1),I,Eta_Cathode(:,2),I,Eta_Cathode(:,3),I,Eta_Cathode(:,4),I,Eta_Cathode(:,5),I,Eta_Cathode(:,6))
xlabel('Current Density (A/cm2)');
ylabel('Overpotential (Volt)');
title('Activation Cathode Overpotential');
figure(8); plot(I,Eta_Ohmic(:,1),I,Eta_Ohmic(:,2),I,Eta_Ohmic(:,3),I,Eta_Ohmic(:,4),I,Eta_Ohmic(:,5),I,Eta_Ohmic(:,6))
xlabel('Current Density (A/cm2)');
ylabel('Overpotential (Volt)');
title('Ohmic Overpotential');
figure(9); plot(I,Eta_Conc(:,1),I,Eta_Conc(:,2),I,Eta_Conc(:,3),I,Eta_Conc(:,4),I,Eta_Conc(:,5),I,Eta_Conc(:,6))
xlabel('Current Density (A/cm2)');
ylabel('Overpotential (Volt)');
title('Concentration Overpotential');
0 个评论
回答(2 个)
Torsten
2024-1-26
编辑:Torsten
2024-1-26
You mean
ia = (10^-7)*(exp(-(Eexc/R)*((1/T)-(1/298))));
ic = (10^-3)*(exp(-(Eexc/R)*((1/T)-(1/353))));
instead of
ia = (10^-7)(exp(-(Eexc/R)((1/T)-(1/298))));
ic = (10^-3)(exp(-(Eexc/R)((1/T)-(1/353))));
?
And the second dimension of V, VP, VT, Eta_anode,... is only 3, not 6, so that you will get a dimension error in the plots.
0 个评论
Dyuman Joshi
2024-1-26
There are missing operators at many places in your code, e.g. -
% vv vv
ia = (10^-7)(exp(-(Eexc/R)((1/T)-(1/298))));
% vv vv
ic = (10^-3)(exp(-(Eexc/R)((1/T)-(1/353))));
Use appropriate operators and this issue will be resolved.
2 个评论
Dyuman Joshi
2024-1-26
Additionally, V is 8x3 where as you trying to access the 4th, 5th and 6th column of it, which don't exist.
Same for VT and VP.
I suggest you go through your code thoroughly once and check each line.
另请参阅
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!