% Modified version of the AFCoKaRE MATLAB template file SolvNonDif.m used
% to model the recycle PFR in Example 31.1 of "A First Course on Kinetics 
% and Reaction Engineering."
%
function Example_31_1
    % Known quantities and constants (in consistent units)
    CSa = 0.04; % g/cc
    CXa = 0.; % g/cc
    VFRa = 60.; % cc/min
    RR = 0.2;
    mSa = VFRa*CSa;
    mXa = VFRa*CXa;
    VFRd = VFRa;
    VFRr = RR*VFRd;
    VFRc = VFRr + VFRd;
    VFRb = VFRc;

    % Function that evaluates the equations
    function f = evalEqns(z)
        mSb = z(1);
        mXb = z(2);
        [mSc,mXc] = Example_31_1_pfr(mSb,mXb,VFRb);
        f = [
            mSa + RR*mSc/(1+RR) - mSb;
            mXa + RR*mXc/(1+RR) - mXb;
        ];
    end % of internal function evalEqns
    
    % guesses for the solution
    z_guess = [
        (1+RR)*mSa*0.5
        RR*mSa*0.2
    ];

    % Solve the set of algebraic equations
    z = fsolve(@evalEqns, z_guess);
    display('The solver found the following values for the unknowns:');
    z
    display('The corresponding values of the functions being solved are as follows:');
    f = evalEqns(z)
    
    % Solve the PFR design equations
    mSb = z(1);
    mXb = z(2);
    [mSc,mXc] = Example_31_1_pfr(mSb,mXb,VFRb);
    % Calculate the outlet cell mass concentration
    mXd = mXc/(1+RR);
    CXd = mXd/VFRd
	
end % of Example_31_1.m