fplot error odes help

2 visualizaciones (últimos 30 días)
Easton Olsen
Easton Olsen el 24 de Abr. de 2019
I cant seem to figure out what is wrong with my odes can someone please help
%% Units: g, mol, dm^3, K, Kpa, s
clear
clc
%Specifying operating conditions
R = 8.314472; %Ideal Gas Constant
m_a = 60.052; % Molar mass acetic acid
m_b = 74.1; % Molar mass Butanol
phi = 0.2; %Porosity of Catalyst
rho_c = 19300; %Density of Catalyst
P_0 = input('Please specify an operating pressure in Kpa.');
while P_0 <= 0
P_0 = input('Please specify an operating pressure greater than zero.');
end
T = input('Please specify an operating temperature in K.');
while T <= 0
T = input('Please specify an operating temperature greater than zero.');
end
Dp = input('Please specify a pipe diameter in dm.');
while Dp <= 0
Dp = input('Please specify a pipe diameter greater than zero.');
end
F_t_0 = input('Please specify incoming molar flow rate in mol/s.');
while F_t_0 <= 0
F_t_0 = input('Please specify an incoming molar flow rate greater than zero.');
end
x_b = input('Please specify mole fraction of Butanol in feed.');
while x_b <= 0 || x_b > 1
x_b = input('Please specify a mole fraction between zero and one.');
end
x_a = 1 - x_b; %Mole fraction of acetic acid
%Calculating parameters for pressure drop equation
Ap = pi*(Dp/2)^2; %Cross Sectional Area of Pipe
Cto = P_0/(R*T); %Initial concentration
rho_a = x_a*Cto*m_a; %Density of acetic acid
rho_b = x_b*Cto*m_b; %Density of butanol
rho = x_a*rho_a + x_b*rho_b; %Approximate total density
mu_a = 0.001155; %Viscosity of acetic acid
mu_b = 0.002593; %Viscosity of butanol
mu = x_a*mu_a + x_b*mu_b; %Approximate total viscosity
m_flow = x_a*F_t_0*m_a + x_b*F_t_0*m_b; %Incoming mass flow rate
G = m_flow/Ap; %Mass velocity
beta_1 = (G*(1-phi)/(rho*Dp*phi^3));
beta_2 = ((150*(1-phi)*mu)/Dp) + 1.75*G;
beta = beta_1*beta_2;
rho_bulk_c = rho_c*(1 - phi); %Bulk density of Catalyst
alpha = ((2*beta)/(Ap*rho_bulk_c*P_0));
%Calculating parameters for rate law
k_f = 0.0090953*exp(-45.49/R*T);
k_r = 0.008643*exp(-23.90/R*T);
%ODE solver
syms Fa(W) Fb(W) Fc(W) Fd(W) p(W)
ode1 = diff(Fa) == (Cto^2*p^2/(Fa+Fb+Fc+Fd)^2)*(k_f*Fa*Fb - k_r*Fc*Fd);
ode2 = diff(Fb) == (Cto^2*p^2/(Fa+Fb+Fc+Fd)^2)*(k_f*Fa*Fb - k_r*Fc*Fd);
ode3 = diff(Fc) == -1*(Cto^2*p^2/(Fa+Fb+Fc+Fd)^2)*(k_f*Fa*Fb - k_r*Fc*Fd);
ode4 = diff(Fd) == -1*(Cto^2*p^2/(Fa+Fb+Fc+Fd)^2)*(k_f*Fa*Fb - k_r*Fc*Fd);
ode5 = diff(p) == (-1*alpha/2*p)*((Fa+Fb+Fc+Fc)/F_t_0);
odes = [ode1; ode2; ode3; ode4; ode5];
[FaSol(W), FbSol(W), FcSol(W), FdSol(W), pSol(W)] = dsolve(odes);
%Plotting results
fplot(FaSol)
hold on
fplot(FbSol)
fplot(FcSol)
fplot(FdSol)
fplot(pSol)
grid on
Error using fplot>singleFplot (line 229)
Input must be a function or functions of a single variable.
Error in fplot>@(f)singleFplot(cax,{f},limits,extraOpts,args) (line 193)
hObj = cellfun(@(f) singleFplot(cax,{f},limits,extraOpts,args),fn{1},'UniformOutput',false);
Error in fplot>vectorizeFplot (line 193)
hObj = cellfun(@(f) singleFplot(cax,{f},limits,extraOpts,args),fn{1},'UniformOutput',false);
Error in fplot (line 163)
hObj = vectorizeFplot(cax,fn,limits,extraOpts,args);
Error in reactions (line 60)
fplot(FaSol)
>>

Respuestas (0)

Categorías

Más información sobre Ordinary Differential Equations en Help Center y File Exchange.

Etiquetas

Community Treasure Hunt

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

Start Hunting!

Translated by