Differential Equation Problem using ode45

6 ビュー (過去 30 日間)
Luqman Hardiyan
Luqman Hardiyan 2019 年 11 月 25 日
回答済み: Stephan 2019 年 11 月 25 日
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

採用された回答

Stephan
Stephan 2019 年 11 月 25 日
[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

その他の回答 (0 件)

カテゴリ

Help Center および File ExchangeOrdinary Differential Equations についてさらに検索

Community Treasure Hunt

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

Start Hunting!

Translated by