Differential Equation Problem using ode45

10 vues (au cours des 30 derniers jours)
Luqman Hardiyan
Luqman Hardiyan le 25 Nov 2019
I was solving the molar flow rate or dF/dW using this equation, it says that there are error on line 32.
Why there's an error on line 32 saying not enough input arguments? where did I go wrong? thank you
function PBR_Isothermal
clc; clear;
%A = Dammitol
%B = Oxygen
%C = Valualdehyde
%D = Carbon Dioxide
%E = Water
%I = Nitrogen
%Feed based on inlet as 100 kmol/h
Fa0 = 0.1*100;
Fb0 = 0.07*100;
Fc0 = 0;
Fd0 = 0;
Fe0 = 0.02*100;
Fi0 = 0.81*100;
[Fa W] = ode45(FunA,[0 1000],Fa0);
% [Fb W] = ode45(FunA,[0 1000],Fb0);
% [Fc W] = ode45(FunA,[0 1000],Fc0);
% [Fd W] = ode45(FunA,[0 1000],Fd0);
% [Fe W] = ode45(FunA,[0 1000],Fe0);
end
function A = FunA(F,W)
%Total Pressure
Ptot = 114.7;
%Defining Constant
FT = F(1)+F(2)+F(3)+F(4)+F(5)+F(6);
RT = (1/353-1/373)/1.987;
%Partial Pressure of Each
PDT = (F(1)/FT)*Ptot;
PO2 = (F(2)/FT)*Ptot;
PVA = (F(3)/FT)*Ptot;
PCO = (F(4)/FT)*Ptot;
PWA = (F(5)/FT)*Ptot;
PN2 = (F(6)/FT)*Ptot;
%Constants
k1 = 1.771*10^-3;
k2 = 23295;
k3 = 0.5;
k4 = 1.0;
k5 = 0.8184;
k6 = 0.5;
k7 = 0.2314;
k8 = 1.0;
k9 = 1.25;
k10 = 2.0;
k11 = 2.795*10^-4;
k12 = 33000;
k13 = 0.5;
k14 = 2.0;
%Rate of Reaction
r1 = (k1*exp(-k2*RT)*PO2^k3*PDT^k4)/(1+k5*PO2^k6+k7*PDT^k8+k9*PVA^k10);
r2 = (k11*exp(-k12*RT)*PO2^k13*PVA^k14)/(1+k5*PO2^k6+k7*PDT^k8+k9*PVA^k10);
%Molar Flow Rate of Each Component
F(1) = -r1;
F(2) = -1/2*r1;
F(3) = r1;
F(4) = 2*r2;
F(5) = r1+2*r2;
F(6) = Fi0;
end

Réponse acceptée

Stephan
Stephan le 25 Nov 2019
[W,Fa] = PBR_Isothermal;
% plot results
subplot(3,2,1)
plot(W,Fa(:,1))
subplot(3,2,2)
plot(W,Fa(:,2))
subplot(3,2,3)
plot(W,Fa(:,3))
subplot(3,2,4)
plot(W,Fa(:,4))
subplot(3,2,5)
plot(W,Fa(:,5))
subplot(3,2,6)
plot(W,Fa(:,6))
function [W,Fa] = PBR_Isothermal
clc; clear;
%A = Dammitol
%B = Oxygen
%C = Valualdehyde
%D = Carbon Dioxide
%E = Water
%I = Nitrogen
%Feed based on inlet as 100 kmol/h
Fa0 = 0.1*100;
Fb0 = 0.07*100;
Fc0 = 0;
Fd0 = 0;
Fe0 = 0.02*100;
Fi0 = 0.81*100;
[W, Fa] = ode45(@FunA,[0 1000],[Fa0 Fb0 Fc0 Fd0 Fe0 Fi0]);
% [Fb W] = ode45(FunA,[0 1000],Fb0);
% [Fc W] = ode45(FunA,[0 1000],Fc0);
% [Fd W] = ode45(FunA,[0 1000],Fd0);
% [Fe W] = ode45(FunA,[0 1000],Fe0);
end
function A = FunA(W,F)
%Total Pressure
Ptot = 114.7;
%Defining Constant
FT = F(1)+F(2)+F(3)+F(4)+F(5)+F(6);
RT = (1/353-1/373)/1.987;
%Partial Pressure of Each
PDT = (F(1)/FT)*Ptot;
PO2 = (F(2)/FT)*Ptot;
PVA = (F(3)/FT)*Ptot;
PCO = (F(4)/FT)*Ptot;
PWA = (F(5)/FT)*Ptot;
PN2 = (F(6)/FT)*Ptot;
%Constants
k1 = 1.771*10^-3;
k2 = 23295;
k3 = 0.5;
k4 = 1.0;
k5 = 0.8184;
k6 = 0.5;
k7 = 0.2314;
k8 = 1.0;
k9 = 1.25;
k10 = 2.0;
k11 = 2.795*10^-4;
k12 = 33000;
k13 = 0.5;
k14 = 2.0;
%Rate of Reaction
r1 = (k1*exp(-k2*RT)*PO2^k3*PDT^k4)/(1+k5*PO2^k6+k7*PDT^k8+k9*PVA^k10);
r2 = (k11*exp(-k12*RT)*PO2^k13*PVA^k14)/(1+k5*PO2^k6+k7*PDT^k8+k9*PVA^k10);
%Molar Flow Rate of Each Component
A(1) = -r1;
A(2) = -1/2*r1;
A(3) = r1;
A(4) = 2*r2;
A(5) = r1+2*r2;
Fi0 = 0.81*100;
A(6) = Fi0;
A=A';
end

Plus de réponses (0)

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by