## Coupled ODEs by ode45 and could not get the figures

Hi everyone!

I’ve been learning MATLAB for about six months, so I apologize if I ask any basic questions. In part of my project, I have a system of six coupled ODEs (all functions of t). You can refer to the attached picture for more details.

To obtain the T profile in terms of t, I decided to define the function similarly to how many others do:

function dYdt = complexSystem(t, Y)

Aa = 1e14; % Pre-exponential factor for xa

Ea = 2.24e-19; % Activation energy for xa

Ac = 4e13; % Pre-exponential factor for xc

Ec = 2.03e-19; % Activation energy for xc

As = 1e15; % Pre-exponential factor for xs

Es = 2.24e-19; % Activation energy for xs

Aec = 1e12; % Pre-exponential factor for xec

Eec = 1.4e-19; % Activation energy for xec

Kb = 1.38e-23; % Boltzmann constant

z0 = 0.033; % Initial value for z

ma = 0.13;

ha = 1.714e6;

mc = 0.29;

hca = 3.14e5;

hs = 2.57e5;

C = -3600 * 1.1 * 3.0 * (1-0.12-0.51);

T = Y(1);

Xa = Y(2);

Xc = Y(3);

Xs = Y(4);

Z = Y(5);

Soc = Y(6);

dXadt = -Aa * Xa * exp(-Ea/(Kb*T));

dXcdt = Ac * Xc * exp(1 – Xc) * exp(-Ec / (Kb * T));

dXsdt = -As * Xs * exp((-Z / z0) * (-Es / (Kb * T)));

dZdt = Aa * Xs * exp((-Z / z0) * (-Es / (Kb * T)));

dSocdt = -Aec * (1 – Xc) * Xa * exp(-Eec / (Kb * T));

Qa = -ma * ha * dXadt;

Qc = -mc * hca * dXcdt;

Qs = -ma * hs * dXsdt;

Qec = C * dSocdt;

Q = Qa + Qc + Qs + Qec;

dTdt = Q;

dYdt = [dTdt; dXadt; dXcdt; dXsdt; dZdt; dSocdt];

end

Then, I defined the initial conditions and solved the ODEs using ode45.

% Initial conditions

Y0 = [298; 0.75; 0.04; 0.015; 0.0331; 0]; % [T0; Xa0; Xc0; Xs0; Z0; Soc0]

% Time span

tspan = [0 100000];

% Solve ODEs

[T, Y] = ode45(@(t, Y) complexSystem(t, Y), tspan, Y0);

Lastly, I attempted to plot the results. However, the code took forever to run, and I didn’t get any figures from it.

Could my initial conditions be incorrect? Or are my ODEs not converging?

Thank you so much for your time and help!

figure;

plot(T, Y)

legend(‘T’, ‘X_a’, ‘X_c’, ‘X_s’, ‘Z’, ‘S_{oc}’)

title(‘Solution of Corrected Complex System’)

xlabel(‘Time t’)

ylabel(‘State Variables’)Hi everyone!

I’ve been learning MATLAB for about six months, so I apologize if I ask any basic questions. In part of my project, I have a system of six coupled ODEs (all functions of t). You can refer to the attached picture for more details.

To obtain the T profile in terms of t, I decided to define the function similarly to how many others do:

function dYdt = complexSystem(t, Y)

Aa = 1e14; % Pre-exponential factor for xa

Ea = 2.24e-19; % Activation energy for xa

Ac = 4e13; % Pre-exponential factor for xc

Ec = 2.03e-19; % Activation energy for xc

As = 1e15; % Pre-exponential factor for xs

Es = 2.24e-19; % Activation energy for xs

Aec = 1e12; % Pre-exponential factor for xec

Eec = 1.4e-19; % Activation energy for xec

Kb = 1.38e-23; % Boltzmann constant

z0 = 0.033; % Initial value for z

ma = 0.13;

ha = 1.714e6;

mc = 0.29;

hca = 3.14e5;

hs = 2.57e5;

C = -3600 * 1.1 * 3.0 * (1-0.12-0.51);

T = Y(1);

Xa = Y(2);

Xc = Y(3);

Xs = Y(4);

Z = Y(5);

Soc = Y(6);

dXadt = -Aa * Xa * exp(-Ea/(Kb*T));

dXcdt = Ac * Xc * exp(1 – Xc) * exp(-Ec / (Kb * T));

dXsdt = -As * Xs * exp((-Z / z0) * (-Es / (Kb * T)));

dZdt = Aa * Xs * exp((-Z / z0) * (-Es / (Kb * T)));

dSocdt = -Aec * (1 – Xc) * Xa * exp(-Eec / (Kb * T));

Qa = -ma * ha * dXadt;

Qc = -mc * hca * dXcdt;

Qs = -ma * hs * dXsdt;

Qec = C * dSocdt;

Q = Qa + Qc + Qs + Qec;

dTdt = Q;

dYdt = [dTdt; dXadt; dXcdt; dXsdt; dZdt; dSocdt];

end

Then, I defined the initial conditions and solved the ODEs using ode45.

% Initial conditions

Y0 = [298; 0.75; 0.04; 0.015; 0.0331; 0]; % [T0; Xa0; Xc0; Xs0; Z0; Soc0]

% Time span

tspan = [0 100000];

% Solve ODEs

[T, Y] = ode45(@(t, Y) complexSystem(t, Y), tspan, Y0);

Lastly, I attempted to plot the results. However, the code took forever to run, and I didn’t get any figures from it.

Could my initial conditions be incorrect? Or are my ODEs not converging?

Thank you so much for your time and help!

figure;

plot(T, Y)

legend(‘T’, ‘X_a’, ‘X_c’, ‘X_s’, ‘Z’, ‘S_{oc}’)

title(‘Solution of Corrected Complex System’)

xlabel(‘Time t’)

ylabel(‘State Variables’) Hi everyone!

I’ve been learning MATLAB for about six months, so I apologize if I ask any basic questions. In part of my project, I have a system of six coupled ODEs (all functions of t). You can refer to the attached picture for more details.

To obtain the T profile in terms of t, I decided to define the function similarly to how many others do:

function dYdt = complexSystem(t, Y)

Aa = 1e14; % Pre-exponential factor for xa

Ea = 2.24e-19; % Activation energy for xa

Ac = 4e13; % Pre-exponential factor for xc

Ec = 2.03e-19; % Activation energy for xc

As = 1e15; % Pre-exponential factor for xs

Es = 2.24e-19; % Activation energy for xs

Aec = 1e12; % Pre-exponential factor for xec

Eec = 1.4e-19; % Activation energy for xec

Kb = 1.38e-23; % Boltzmann constant

z0 = 0.033; % Initial value for z

ma = 0.13;

ha = 1.714e6;

mc = 0.29;

hca = 3.14e5;

hs = 2.57e5;

C = -3600 * 1.1 * 3.0 * (1-0.12-0.51);

T = Y(1);

Xa = Y(2);

Xc = Y(3);

Xs = Y(4);

Z = Y(5);

Soc = Y(6);

dXadt = -Aa * Xa * exp(-Ea/(Kb*T));

dXcdt = Ac * Xc * exp(1 – Xc) * exp(-Ec / (Kb * T));

dXsdt = -As * Xs * exp((-Z / z0) * (-Es / (Kb * T)));

dZdt = Aa * Xs * exp((-Z / z0) * (-Es / (Kb * T)));

dSocdt = -Aec * (1 – Xc) * Xa * exp(-Eec / (Kb * T));

Qa = -ma * ha * dXadt;

Qc = -mc * hca * dXcdt;

Qs = -ma * hs * dXsdt;

Qec = C * dSocdt;

Q = Qa + Qc + Qs + Qec;

dTdt = Q;

dYdt = [dTdt; dXadt; dXcdt; dXsdt; dZdt; dSocdt];

end

Then, I defined the initial conditions and solved the ODEs using ode45.

% Initial conditions

Y0 = [298; 0.75; 0.04; 0.015; 0.0331; 0]; % [T0; Xa0; Xc0; Xs0; Z0; Soc0]

% Time span

tspan = [0 100000];

% Solve ODEs

[T, Y] = ode45(@(t, Y) complexSystem(t, Y), tspan, Y0);

Lastly, I attempted to plot the results. However, the code took forever to run, and I didn’t get any figures from it.

Could my initial conditions be incorrect? Or are my ODEs not converging?

Thank you so much for your time and help!

figure;

plot(T, Y)

legend(‘T’, ‘X_a’, ‘X_c’, ‘X_s’, ‘Z’, ‘S_{oc}’)

title(‘Solution of Corrected Complex System’)

xlabel(‘Time t’)

ylabel(‘State Variables’) ode45 MATLAB Answers — New Questions