I get this error even though both the initial conditions vector and the function output vector have 10 elements. Could you help me how to solve it?

% Condizioni iniziali di composizione (frazioni molari)

Y0_CH4_values = [100 / 1E6, 1000 / 1E6, 10000 / 1E6]; % ppm convertiti in frazioni molari

Y0_O2 = 0.21;

Y0_CO2 = 0;

Y0_H2O = 0;

% Range di temperature

T_range = 273.15:10:800; % Vettore di temperature [K]

% Preallocazione del vettore di conversioni

conversion_matrix = zeros(length(T_range), length(Y0_CH4_values));

% Definizione del sistema di equazioni differenziali

NEQ = 10;

M = zeros(NEQ, NEQ);

for i = 1:4

M(i, i) = 1;

end

M(NEQ, NEQ) = 1; % Pressione

options = odeset('Mass', M, 'AbsTol', 1e-06, 'RelTol', 1e-06);

% Risoluzione delle ODE per ogni temperatura e concentrazione iniziale di CH4

for j = 1:length(Y0_CH4_values)

Y0_CH4 = Y0_CH4_values(j);

conversion_vector = zeros(size(T_range));

for i = 1:length(T_range)

T = T_range(i);

% Sezione 1: Zona senza catalizzatore (inizio)

zspan_no_catalyst_1 = [0, z_start]; % da 0 a z_start

y0 = [C0_CH4; C0_O2; C0_CO2; C0_H2O; P; u];

[z1, y1] = ode45(@(z, y) reactor_model_no_reaction(z, y, T), zspan_no_catalyst_1, y0);

C0_CH4_sol = y1(end, 1);

C0_O2_sol = y1(end, 2);

C0_CO2_sol = y1(end, 3);

C0_H2O_sol = y1(end, 4);

% Sezione 2: Zona con catalizzatore

zspan_with_catalyst = [z_start, z_end]; % da z_start a z_end

y0_catalyst = [y1(end, 1:4)'; C0_CH4_sol; C0_O2_sol; C0_CO2_sol; C0_H2O_sol; y1(end, 5:6)'];

[z2, y2] = ode15s(@(z, y) reactor_model_with_catalyst(z, y, T), zspan_with_catalyst, y0_catalyst);

C_end_CH4 = y2(end, 1);

C_end_O2 = y2(end, 2);

C_end_CO2 = y2(end, 3);

C_end_H2O = y2(end, 4);

% Sezione 3: Zona senza catalizzatore (fine)

zspan_no_catalyst_2 = [z_end, L_total]; % da z_end a L_total

y0_no_catalyst_2 = [C_end_CH4; C_end_O2; C_end_CO2; C_end_H2O; y2(end, 5); y2(end, 6)];

[z3, y3] = ode45(@(z, y) reactor_model_no_reaction(z, y, T), zspan_no_catalyst_2, y0_no_catalyst_2);

C_end_CH4_final = y3(end, 1);

% Calcolo della conversione

conversion = (C0_CH4 - C_end_CH4_final) / C0_CH4;

conversion_vector(i) = conversion;

end

% Salva i risultati nella matrice

conversion_matrix(:, j) = conversion_vector;

end

% Funzione del modello senza reazione

function dydz = reactor_model_no_reaction(z, y, T)

global kc eps eta

C_CH4 = y(1);

C_O2 = y(2);

C_CO2 = y(3);

C_H2O = y(4);

u = y(6);

% Velocità di reazione per la zona senza catalizzatore

dydz = zeros(6, 1);

dydz(1) = -0; % Non c'è reazione senza catalizzatore

dydz(2) = -0;

dydz(3) = 0;

dydz(4) = 0;

dydz(5) = 0; % Pressione costante

dydz(6) = 0; % Velocità costante

end

% Funzione del modello con catalizzatore

function dydz = reactor_model_with_catalyst(z, y, T)

global kc eps eta

C_CH4 = y(1);

C_O2 = y(2);

C_CO2 = y(3);

C_H2O = y(4);

P = y(5);

u = y(6);

% Velocità di reazione

k = 4204.88 * exp(-39545.5 / (8.314 * T)); % [1/s]

r_CH4 = -k * C_CH4^0.792 * eta * (1 - eps); % velocità di reazione per CH4

r_O2 = 2 * r_CH4; % velocità di reazione per O2

r_CO2 = -r_CH4; % velocità di reazione per CO2

r_H2O = -2 * r_CH4; % velocità di reazione per H2O

% Bilanci di massa

dydz = zeros(6, 1);

dydz(1) = -kc * (C_CH4 - C_CH4) / (u * eps) - r_CH4; % Bilancio di CH4

dydz(2) = -kc * (C_O2 - C_O2) / (u * eps) - r_O2; % Bilancio di O2

dydz(3) = kc * (C_CO2 - C_CO2) / (u * eps) + r_CO2; % Bilancio di CO2

dydz(4) = kc * (C_H2O - C_H2O) / (u * eps) + r_H2O; % Bilancio di H2O

dydz(5) = 0; % Pressione costante

dydz(6) = 0; % Velocità costante

end

Error using odearguments

@(Z,Y)REACTOR_MODEL_WITH_CATALYST(Z,Y,T) returns a vector of length 6, but the length of initial

conditions vector is 10. The vector returned by @(Z,Y)REACTOR_MODEL_WITH_CATALYST(Z,Y,T) and the

initial conditions vector must have the same number of elements.

Error in ode15s (line 153)

odearguments(odeIsFuncHandle, odeTreatAsMFile, solver_name, ode, tspan, y0, options, varargin);

Error in section_heterogenous (line 120)

[z2, y2] = ode15s(@(z, y) reactor_model_with_catalyst(z, y, T), zspan_with_catalyst, y0_catalyst);

>>

Answer 1

[Aravind]

Hello @Luigi,

I noticed an issue in your code related to the initial conditions for the differential equations. You are using six state variables for the ODEs, and these initial values need to be correctly passed to the ODE solvers. While you have done this correctly for sections 1 and 3, which do not involve the catalyst, there is an error in section 2. You are passing ten variables as the initial condition instead of six, which is causing the problem.

Upon closer inspection, it appears that you are duplicating the initial concentrations of the gases in section 2. To fix this, you should modify the line:

y0_catalyst = [y1(end, 1:4)'; C0_CH4_sol; C0_O2_sol; C0_CO2_sol; C0_H2O_sol; y1(end, 5:6)'];

to:

y0_catalyst = y1(end,1:6)';

Since the initial condition vector already contains “y1(end, 1:4)”, there is no need to include the variables “C0_CH4_sol”, “C0_O2_sol”, “C0_CO2_sol”, and “C0_H2O_sol” again, as they hold the same values as “y1(end, 1:4)”.

Hope this helps!