Differential Equation Problem using ode45

조회 수: 20 (최근 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 CenterFile Exchange에서 Numerical Integration and Differential Equations에 대해 자세히 알아보기

Community Treasure Hunt

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

Start Hunting!

Translated by