Differential Equation Problem using ode45

22 views (last 30 days)
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

Accepted Answer

Stephan
Stephan on 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

More Answers (0)

Categories

Find more on Numerical Integration and Differential Equations in Help Center and File Exchange

Community Treasure Hunt

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

Start Hunting!