Index exceeds matrix dimensions.

    clearvars
    close all
    clc

    global eta eps rext rint delta rhoc cpi PinSweep
    global DH0 Trif mui ki MW S R Ac Dp p ls Tw

    % Gas Constant [J/(kmol*K)]

    R = 8314.45;

    % Catalyst Density [kgcat/m^3]

    rhoc = 2400;

    % Membrane Thickness [m]

    delta = 5e-6;

    % Reaction efficiency [-]

    eta = 0.28;

    %External Radius [m]

    rext = 0.08; 

    % Internal Radius [m]

    rint = 0.045;

    % Cross Sectional Area [m^2]

    Ac = pi*rint^2;

    % Reactor Length [m]

    L = 10;

    % Void Fraction [-]

    eps = 0.6;

    % Stoichiometric Matrix [-]

    nu = [-1    -1      1      1      0      0      0];
          %CO   %H2O   %CO2    %H2  %Inert %H2perm %sweep

    % Stoichiometric Matrix at Permeate side [-]

    nup = [0      0      0      -1      0     -1      0];
          %CO   %H2O   %CO2    %H2  %Inert %H2perm %sweep

    % Specific Heats matrix [kJ/kmol*K]

    cpi = 1e-3*[2.9108e+04   3.3363e+04   2.9370e+04   2.7617e+04   2.9105e+04    2.7617e+04   3.3363e+04     %C1
                8.7730e+03   2.6790e+04   3.4540e+04   9.5600e+04   8.6149e+03    9.5600e+04   2.6790e+04     %C2
                3.0851e+03   2.6105e+03   1.4280e+03   2.4660e+03   1.7016e+03    2.4660e+03   2.6105e+03     %C3
                8.4553e+03   8.8960e+03   2.6400e+04   3.7600e+03   1.0347e+02    3.7600e+03   8.8960e+03     %C4
                1.5382e+03   1.1690e+03   5.8800e+02   5.6760e+02   9.0979e+02    5.6760e+02   1.1690e+03];   %C5
                   %CO          %H2O         %CO2         %H2          %Inert       %H2perm       %sweep  

    % Viscosity [Pa*s]

    mui = [ 1.1127e-06    1.7096e-08      2.148e-06      1.797e-07   6.5592e-07       1.797e-07     1.7096e-08    %C1
            0.5338        1.1146          0.46           0.685       0.6081           0.685         1.1146        %C2
           94.7           0             290             -0.59        54.714          -0.59          0             %C3
            0             0               0            140           0              140             0        ];   %C4

              %CO            %H2O         %CO2         %H2           %Inert         %H2perm         %sweep 

    % Thermal conductivity [W/m*K]

    ki = [  5.9882e-04     6.2041e-06     3.69         2.653e-03  3.3143e-04      2.653e-03     6.2041e-06     %C1
            0.6863         1.3973        -0.3838       0.7452     0.7722          0.745         1.3973         %C2
           57.13           0            964           12          16.323         12             0              %C3
          501.92           0              1.86e+06     0          373.72          0             0          ];  %C4

              %CO         %H2O         %CO2           %H2         %Inert        %H2perm       %sweep 

    % Ideal gas enthalpy of formation [kJ/kmol]  Trif =298.15 K

    DH0 = 1e-3*[-1.1053e+08   -2.41814e+08    -3.9351e+08        0            0            0       -2.41814e+08 ];
                    %CO            %H2O           %CO2          %H2        %Inert       %H2perm       %sweep 

    % Reference Temperature [K]

    Trif = 298.15;

    % Wall Temperature [K]

    Tw = 628.15;

    % Initial Sweep Temperature [K]

    TinSweep = 628.15;

    % Initial Temperature [K]

    Tin = 628;

    % Initial Sweep Pressure [kPa]

    PinSweep = 101;

    % Initial Pressure [kPa]

    Pin = 2376.53;

    % Initial Sweep Flow Rate [kmol/h] 

    FinSweep = 8;

    % Initial Flow Rate [kmol/h]

    Ftot0 = 14;

    % Inlet Compositions at Permeation Side [-]

    F0Perm = 0;

    % Inlet Compositions of reacting Mixture [-]

    x0 = [0.0521   0.4453   0.0272   0.4649   0.0105];
            %CO     %H2O     %CO2      %H2    %Inert

    % MW [kg/kmol]

    MW = [ 28.010        18.015       44.010      2.016       28.013       2.016      18.015 ];
             %CO          %H2O         %CO2        %H2        %Inert      %H2perm      %sweep 

    % Normal Boiling Temperature and corresponding Parameters [K]

    Tb = [81.15   373.15   194.65   79    77.35      79      373.15 ];
           %CO     %H2O     %CO2    %H2   %Inert   %H2perm   %sweep       
    S = 1.5.*Tb; 

    % Global Heat Exchange coefficient for permeate side [kJ/(m^2*h*K)]

    Um = 8.64;

    % Particle Diameter [m]

    Dp = 0.0053; 

    % Thermal Conductivity of packing Material [kJ/(m*h*K)]

    ls = 1.256;

    % Solid Surface Emissivity [-]

    p = 0.8;

    % Inlet Conditions 

    In = [Ftot0*x0, F0Perm, FinSweep, Pin, Tin, TinSweep - 100]; 

    % Sistem Solution

    zspan = [0 L]; % [m]
    options = odeset('Reltol',1e-11,'Abstol',1e-12);
    [Z,X] = ode15s(@Tronky,zspan,In,options);

    function f = Tronky(z,x)

    % Variable Definition

    global eta eps rext rint delta rhoc cpi PinSweep
    global DH0 Trif mui ki MW S R Ac Dp p ls Tw

    Fi = x(1:5);
    FiPS = x(6:7);
    P = x(8);
    T = x(9);
    Tperm = x(10);

    % Other Variables

    Ftot = sum(Fi); % [kmol/h]
    xi = Fi./Ftot; % [-]
    Pi = P*xi; % [kPa]
    FtotPS = sum(FiPS); % [kmol/h]
    xiPS = FiPS./FtotPS; % [-]
    PiPS = PinSweep*xiPS; % [kPa]
    Tmem = (T + Tperm)/2; % [K]

    % Reaction, Adsorption, Equilibrium Constants

    k = exp(-29364/(1.987*T) + 40.32/1.987); % [kmol/(kgcat*h)]
    Kco = exp(3064/(1.987*T) - 6.74/1.987)*1/(101325); % [1/Pa]
    Kh2o = exp(-6216/(1.987*T) + 12.77/1.987); % [1/Pa]
    Kco2 = exp(12542/(1.987*T) - 18.45/1.987); % [1/Pa]
    Keq = exp(4400/T - 4.063); % [-]

    % Reaction Rate [kmol/m^3*h]

    rco = k*Kco*Kh2o*(Pi(1)*Pi(2) - Pi(3)*Pi(4)/Keq)*(1 + Kco*Pi(1) + Kh2o*Pi(2) + Kco2*Pi(4))^-2*rhoc; 
    Rco = rco*eta*(1 - eps)*pi*(rext^2 - rint^2); 

    % Membrane Permeability [kmol/(m*h*Pa^0.5)]

    Bh = 2.95*1e-4*3600/1e3*exp(-5833.5/Tmem);

    % Hydrogen Flux [kmol/(h*m^2)]

    Jh2perm = Bh/delta*(Pi(4)^0.5 - PiPS(2)^0.5);
    Rh2 = Jh2perm*2*pi*rint;

    % Specific Heats, Enthalpies and molecular Weights of reacting Mixture

    for i = 1:5
        cp(i) = cpi(1,i) + cpi(2,i)*((cpi(3,i)/T)/sinh(cpi(3,i)/T))^2 +...
            cpi(4,i)*((cpi(5,i)/T)/cosh(cpi(5,i)/T))^2; % [kJ/(kmol*K)]
        DHi(i) = DH0(i) + (2*cpi(2,i)*cpi(3,i))*((1/(exp(2*cpi(3,i)/T) - 1) - 1/(exp(2*cpi(3,i)/Trif) - 1)) + ...
            (1/(exp(2*cpi(3,i)/T) + 1) - 1/(exp(2*cpi(3,i)/Trif) + 1))); % [kJ/kmol]
        cpmixi(i) = cp(i)*xi(i); % [kJ/(mol*K)]
        MWi(i) = MW(i)*xi(i);
    end

    % Mixture Molecular Weight [kg/kmol]

    MWm = sum(MWi);

    % Mixture Density [kg/m^3]

    rhom = P/(R*T)*MWm;

    % Mixture Specific Heat [kJ/(mol*K)]

    cpm = sum(cpmixi(i));

    % Reaction Heat [kJ/mol]

    DHr = -DHi(1) - DHi(2) + DHi(3) + DHi(4);

    % Viscosities and Conductivities 

    for i = 1:5

        mu = mui(1,i)*(T^(mui(2,i)))/(1 + mui(3,i)/T + mui(4,i)/(T^2)); % [Pa*s]
        k = 3600/1e3*(ki(1,i)*(T^(ki(2,i)))/(1 + ki(3,i)/T + ki(4,i)/(T^2))); % [kJ/(m*h*K)]
    end

    % Mixing Rules for Viscosities and Conductivities

    for i = 1:5
        for j= 1:5

            phi(i,j) = (1 + ((mu(i)/mu(j))^(1/2))*((MW(j)/MW(i))^(1/4)))^2/...
                (8*(1 + MW(i)/MW(j)))^(1/2); % [-]
            Ss(i,j) = 0.735*sqrt(S(i)*S(j)); % [K]
            Ai(i,j) = 0.25*(1 + (mu(i)/mu(j)*(MW(j)/MW(i))^(3/4)*((T + S(i))/(T + S(j))))^(1/2))^2*...
                ((T + Ss(i,j))/(T + S(i))); % [-]
            DENm(i) = phi(i,j)*xi(j); % [-]
            DENk(i) = Ai(i,j)*xi(j); % [-]
        end
        DENmt(i) = sum(DENm); % [-]
        DENkt(i) = sum(DENk); % 
        NOMm(i) = mu(i)*xi(i); % [Pa*s]
        NOMk(i) = k(i)*xi(i); % [kJ/(m*h*K)]

    end
    mum = sum(NOMm)/sum(DENmt);
    km = sum(NOMk)/sum(DENkt);

    % Mixture Velocity [m/s]

    v = Ftot*MWm*1/(rhom*Ac*3600);

    % Adimensional Numbers [-]
    Re = rhom*v*Dp/mum;
    Pr = mum*cpm/km;

    % Parameters for global Heat Exchange Coefficient

    alpharu = (0.8171*(T/100)^3)/(1 + eps/(2*(1 - eps))*((1 - p)/p));
    alphars = 0.8171*(p/(2 - p))*(T/100)^3;
    ler0 = eps*(km + 0.95*alpharu*Dp) + ((0.95*(1 - eps))/(2/(3*ls) + 1/(10*km + alphars*Dp))); % [kJ/(m*K*h)]
    ler = ler0 + 0.111*km*(Re*Pr^(1/3))/(1 + 46*(Dp/(2*rext))^2); % [kJ/(m*K*h)]
    alphaw = 8.694*(ler0/(2*rext)^(4/3)) + 0.512*km*2*rext*Re*Pr^(1/3)/Dp;
    Bi = alphaw*rext/ler; % [-]
    U = (1/alphaw + rext/(3*ler)*(Bi + 3)/(Bi + 4))^(-1); % [kJ/(m^2*K*h)]

    % Specific Heats and Enthalpies of Permeation Side

    for i = 6:7

        cp(i) = cpi(1,i) + cpi(2,i)*((cpi(3,i)/T)/sinh(cpi(3,i)/T))^2 +...
            cpi(4,i)*((cpi(5,i)/T)/cosh(cpi(5,i)/T))^2; % [kJ/(kmol*K)]
        DHi(i) = DH0(i) + (2*cpi(2,i)*cpi(3,i))*((1/(exp(2*cpi(3,i)/T) - 1) - 1/(exp(2*cpi(3,i)/Trif) - 1)) + ...
            (1/(exp(2*cpi(3,i)/T) + 1) - 1/(exp(2*cpi(3,i)/Trif) + 1))); % [kJ/kmol]
        DHTpermi(i) = DH0(i) + (2*cpi(2,i)*cpi(3,i))*((1/(exp(2*cpi(3,i)/Tperm) - 1) - 1/(exp(2*cpi(3,i)/Trif) - 1)) + ...
            (1/(exp(2*cpi(3,i)/Tperm) + 1) - 1/(exp(2*cpi(3,i)/Trif) + 1))); % [kJ/kmol]
    end

    % Ergun Expression

    fv = 150 + 1.75*Re/(1 - eps);
    Pressure = fv*v*mum/Dp^2*((1 - eps)^2/eps^3);

    % Balance Equations

    f(1) = -Rco;
    f(2) = -Rco;
    f(3) = Rco;
    f(4) = Rco - Rh2;
    f(5) = 0;
    f(6) = -Rh2;
    f(7) = 0;
    f(8) = Pressure;
    f(9) = 1/(Ftot*cpm)*(Rco*DHr + U*2*pi*rext*(Tw - T) - Um*2*pi*rint*(T - Tperm));
    f(10) = 1/(FiPS(1)*cp(6) + FiSweep*cp(7))*(Um*2*pi*rint*(T - Tperm) + Rh2*(DHTpermi(6) - DHi(6)));

    f = f(:);

    end

[Geoff Hayes]

Riccardo - please copy and paste the full error message so that we can see which line is generating this error message. The error message is clear - the code is trying to access (get) an element in a matrix using an index that is greater than the dimensions of the matrix.

[Walter Roberson]

Duplicates https://www.mathworks.com/matlabcentral/answers/364780-index-exceeds-matrix-dimensions