MATLAB: Non-linear solver

non-linear solver

I have a equation i=a*exp(-c*F*E/(R*T))+b*exp(-d*F*E/(R*T)) *bold that I'm solving for a,b,c and d. The constants are F=96485, T=298, R=8.314. I'm trying to find a,b,c and d with given E and i value. How do I set it up in MATLAB? I've tried using fsolve but there was an error.

Best Answer

See if this does what you want:

F=96485; 
T=298; 
R=8.314;
% % % i=a*exp(-c*F*E/(R*T))+b*exp(-d*F*E/(R*T));                % Original Equation
% % % a = p(1),  b = p(2),  c = p(3),  d = p(4)                 % Define Parameters
fcn = @(p,E,i) p(1).*exp(-p(3)*F*E/(R*T)) + p(2).*exp(-p(4)*F*E/(R*T)) - i;
E = 42;                                                         % Use Correct Value

i = 24;                                                         % Use Correct Value
p0 = rand(4,1)*1E-6;
[Params,fval] = fsolve(@(p)fcn(p,E,i), p0)

If you have vectors of values for ‘E’ and ‘i’, use a nested loop to iterate through them:

E = [...];
i = [...];
for k1 = 1:numel(E)
    for k2 = 1:numel(i);
        Params{k1,k2} = fsolve(@(p)fcn(p,E(k1),i(k2)), p0);
    end
end

Related Solutions

MATLAB: How to half the value of a constant inside an iteration

maxiter = 100;
iter = 1;
x = 1:1:99;
n = 8*(1/2).^x; % n = 8 , 4, 2, 0.5,...
while iter< maxiter
    gamma = 10 * 10^(-n(iter)); % 10e-8, 10e-4,...
    f = gamma * p;
    iter = iter + 1;
end

edit: changed x to 1:1:99 instead of 1:1:100

MATLAB: What could be the issue if a global variable gets an error “Undefined function or variable ‘F_FEED’ ” when called from a function. The global variable is defined in the file from which the function is called.

Maybe you mean

function main
global F_FEED T0 P0 WEIGHT_cat R Ua F_G0 F_H2O0 F_N20 F_A0 F_AC0 F_AA0 F_FA0 F_T0;
F_FEED = 100;           %TOTAL FLOWRATE OF THE FEED (mol/s)
T0 = 650;          %INLET TEMPERATURE OF THE REACTOR (K)
P0 = 303975;    %INLET PRESSURE OF THE REACTOR (Pa)
WEIGHT_cat = 100;       %CATALYST WEIGHT (kg)
R = ToBeSpecified;       
Ua = ToBeSpecified;
Feed = zeros(8,1);
Feed(1) = 0.014*F_FEED;           %INLET FLOWRATE OF G (mol/s)
Feed(2) = 0.634*F_FEED;           %INLET FLOWRATE OF H2O (mol/s)
Feed(3) = 0.352*F_FEED;           %INLET FLOWRATE OF N2 (mol/s)
Feed(8) = T0; 
F_G0 = Feed(1);                 
F_H2O0 = Feed(2);              
F_N20 = Feed(3);                
F_A0 = Feed(4);               
F_AC0 = Feed(5);                
F_AA0 = Feed(6);                
F_FA0= Feed(7);                
F_T0 = F_FEED;
[W,F] = ode45(@Reactors , [0 WEIGHT_cat], Feed);
end
function dFdW = Reactors(t,F)
global F_FEED T0 P0 WEIGHT_cat R Ua F_G0 F_H2O0 F_N20 F_A0 F_AC0 F_AA0 F_FA0 F_T0;
dFdW = zeros(8,1);
%dFdW = [r_G]
%       [r_H2O]
%       [r_N2]
%       [r_O2]
%       [r_A]
%       [r_AC]
%       [r_AA]
%       [r_FA]
%       [T]     Temperature
%reaction
Hf_G = -577.9;                    %Heat of formation: G (kJ/mol)
Hf_H2O = -241.818;                %Heat of formation: H2O (kJ/mol)
Hf_N2 = 0;                        %Heat of formation: N2 (kJ/mol)
Hf_A = -65;                       %Heat of formation: A (kJ/mol)
Hf_AC = -370.06;                  %Heat of formation: AC (kJ/mol)
Hf_AA = -166.1;                   %Heat of formation: AA (kJ/mol)
Hf_FA = -108.6;                   %Heat of formation: FA (kJ/mol)
Cp_G = -35.3267443 + 0.645106658*F(8) - 0.00061741*F(8)^2 + 0.000000239229*F(8)^3;               %provide integrated form
Cp_H2O = 32.95266198 - 0.00391977*F(8) + 0.0000222314*F(8)^2 - 0.0000000104953*F(8)^3;
Cp_N2 = 28.883 -0.00157*F(8)+ 0.00000808*F(8)^2 -0.000000002871*F(8)^3;
Cp_A = -14.9617445 + 0.373240705*F(8) - 0.00034559*F(8)^2 + 0.000000129748*F(8)^3;
Cp_AC = -11.7479865 + 0.402868987*F(8) - 0.00031782*F(8)^2 + 0.000000109142*F(8)^3;
Cp_AA = -5.82005966 + 0.23787576*F(8) - 0.00017362*F(8)^2 + 0.0000000543415*F(8)^3;
Cp_FA = 19.10519375 + 0.048680586*F(8) + 0.00000878405*F(8)^2 - 0.0000000150776*F(8)^3;
integraldeltaCp_rxn_1 = 2*(32.95266198*(F(8)-298) - 0.5*0.00391977*(F(8)^2-298^2) + (1/3)*0.0000222314*(F(8)^3-298^3) - (1/4)*0.0000000104953*(F(8)^4-298^4)) + (-14.9617445*(F(8)-298) + 0.5*0.373240705*(F(8)^2-298^2) - (1/3)*0.00034559*(F(8)^3-298^3) + (1/4)*0.000000129748*(F(8)^4-298^4)) - ((-35.3267443*(F(8)-298) + 0.5*0.645106658*(F(8)^2-298^2) - (1/3)*0.00061741*(F(8)^3-298^3) + (1/4)*0.000000239229*(F(8)^4-298^4))); 
integraldeltaCp_rxn_2 = (32.95266198*(F(8)-298) - 0.5*0.00391977*(F(8)^2-298^2) + (1/3)*0.0000222314*(F(8)^3-298^3) - (1/4)*0.0000000104953*(F(8)^4-298^4) + (-11.7479865*(F(8)-298) + 0.5*0.402868987*(F(8)^2-298^2) - (1/3)*0.00031782*(F(8)^3-298^3) + (1/4)*0.000000109142*(F(8)^4-298^4)) - (32.95266198*(F(8)-298) - 0.5*0.00391977*(F(8)^2-298^2) + (1/3)*0.0000222314*(F(8)^3-298^3) - (1/4)*0.0000000104953*(F(8)^4-298^4)));
integraldeltaCp_rxn_3 = (32.95266198*(F(8)-298) - 0.5*0.00391977*(F(8)^2-298^2) + (1/3)*0.0000222314*(F(8)^3-298^3) - (1/4)*0.0000000104953*(F(8)^4-298^4) + (19.10519375*(F(8)-298) + 0.5*0.048680586*(F(8)^2-298^2) + (1/3)*0.00000878405*(F(8)^3-298^3) - (1/4)*0.0000000150776*(F(8)^4-298^4)) + (-5.82005966*(F(8)-298) + 0.5*0.23787576*(F(8)^2-298^2) - (1/3)*0.00017362*(F(8)^3-298^3) + (1/4)*0.0000000543415*(F(8)^4-298^4)) - (-35.3267443*(F(8)-298) + 0.5*0.645106658*(F(8)^2-298^2) - (1/3)*0.00061741*(F(8)^3-298^3) + (1/4)*0.000000239229*(F(8)^4-298^4)));
integraldeltaCp_rxn_4 = (32.95266198*(F(8)-298) - 0.5*0.00391977*(F(8)^2-298^2) + (1/3)*0.0000222314*(F(8)^3-298^3) - (1/4)*0.0000000104953*(F(8)^4-298^4) + (19.10519375*(F(8)-298) + 0.5*0.048680586*(F(8)^2-298^2) + (1/3)*0.00000878405*(F(8)^3-298^3) - (1/4)*0.0000000150776*(F(8)^4-298^4)) + (-5.82005966*(F(8)-298) + 0.5*0.23787576*(F(8)^2-298^2) - (1/3)*0.00017362*(F(8)^3-298^3) + (1/4)*0.0000000543415*(F(8)^4-298^4)) - (-35.3267443*(F(8)-298) + 0.5*0.645106658*(F(8)^2-298^2) - (1/3)*0.00061741*(F(8)^3-298^3) + (1/4)*0.000000239229*(F(8)^4-298^4)));
Hrxn_1 =  2*Hf_H2O + Hf_A - Hf_G + (integraldeltaCp_rxn_1);               %Heat of reaction: reaction 1
Hrxn_2 =  2*Hf_H2O + Hf_AC - Hf_G + (integraldeltaCp_rxn_2);              %Heat of reaction: reaction 2
Hrxn_3 =  Hf_H2O + Hf_AA + Hf_FA - Hf_G + (integraldeltaCp_rxn_3);        %Heat of reaction: reaction 3
Hrxn_4 =  Hf_H2O + Hf_AA + Hf_FA - Hf_G + (integraldeltaCp_rxn_4);        %Heat of reaction: reaction 4
%INITIAL CONCENTRATIONS  
VOLUMETRICFLOWRATEin = F_FEED*R*T0;                         %Initial volumetric flowrate into the reactor (m^3/s) (assuming ideal gas in)   NEGLECTING P/P0 to start...
C_G0 =      F_G0/VOLUMETRICFLOWRATEin;                      %Inital concentration of G
C_H2O0 =    F_H2O0/VOLUMETRICFLOWRATEin;                    %Inital concentration of H2O
C_N20 =     F_N20/VOLUMETRICFLOWRATEin;                     %Inital concentration of N2
C_A0 =  0;                                                  %Outlet concentration of A

C_AC0 = 0;                                                  %Outlet concentration of AC

C_AA0 = 0;                                                  %Outlet concentration of AA

C_FA0 = 0;                                                  %Outlet concentration of FA

C_T0 =      F_FEED/VOLUMETRICFLOWRATEin;                    %Initial total concentration
%FINAL Concentrations
F_T = F(1) + F(2) + F(3) + F(4) + F(5) + F(6) + F(7);       
C_G =   C_T0*(F(1)/F_T)*(T0/F(8));    %Outlet concentration of G              
C_H2O = C_T0*(F(2)/F_T)*(T0/F(8));      %Outlet concentration of H2O
C_N2 =  C_T0*(F(3)/F_T)*(T0/F(8));      %Outlet concentration of N2              NEGLECTING P/P0 to start...
C_A =   C_T0*(F(4)/F_T)*(T0/F(8));    %Outlet concentration of A
C_AC =  C_T0*(F(5)/F_T)*(T0/F(8));      %Outlet concentration of AC
C_AA =  C_T0*(F(6)/F_T)*(T0/F(8));      %Outlet concentration of AA
C_FA =  C_T0*(F(7)/F_T)*(T0/F(8));      %Outlet concentration of FA
%k-VALUES 
k1 = 23.3244*exp(-45007.1/(R*F(8)));           %Reaction rate constant k1 (m^3/kgcat*s)
k2 = 2528.8*exp(-79244.89/(R*F(8)));           %Reaction rate constant k2 (m^3/kgcat*s)
k3 = 2150.595*exp(-73702.78/(R*F(8)));         %Reaction rate constant k3 (m^3/kgcat*s)
k4 = 0.146768*exp(-24596.97/(R*F(8)));         %Reaction rate constant k4 (m^3/kgcat*s)
%reaction rates
r_1 = -k1*C_G;                                     %Rate of reaction 1
r_2 = -k2*C_G;                                     %Rate of reaction 2
r_3 = -k3*C_G;                                     %Rate of reaction 3
r_4 = -k4*C_G;                                     %Rate of reaction 4
r_G = -(k1+k2+k3+k4)*C_G;                  %Rate of formation of G
r_N2 = 0;                                  %Rate of formation of N2
r_H2O = k1*(C_G)^2+(k3+k4)*G_G;            %Rate of formation of H2O
r_A = k1*C_G;                              %Rate of formation of A
r_AC = k2*C_G;                             %Rate of formation of AC
r_AA = (k3+k4)*C_G;                        %Rate of formation of AA
r_FA = (k3+k4)*C_G;                        %Rate of formation of FA
%Differential Equations to Solve
dFdW(1) = (r_G);
dFdW(2) = (r_H2O);
dFdW(3) = (r_N2);
dFdW(4) = (r_A);
dFdW(5) = (r_AC);
dFdW(6) = (r_AA);
dFdW(7) = (r_FA);
dFdW(8) = (Ua(Ta-T0)+(-r_1)*(-Hrxn_1)+(-r_2)*(-Hrxn_2)+(-r_3)*(-Hrxn_3)+(-r_4)*(-Hrxn_4))/(F(1)*Cp_G + F(2)*Cp_H2O + F(3)*Cp_N2 + F(4)*Cp_A + F(5)*Cp_AC + F(6)*Cp_AA + F(7)*Cp_FA); %dTdW from pg 545 in Fogler;   
end

Best Answer

Related Solutions

MATLAB: How to half the value of a constant inside an iteration

MATLAB: What could be the issue if a global variable gets an error “Undefined function or variable ‘F_FEED’ ” when called from a function. The global variable is defined in the file from which the function is called.

Related Question