the degradation of ozone in aqueous solutions (water). When I plot the all three components, the plot stays at the 0 value for all components.

function ProjectOzone

clc

ki = 7e-3;

kj = 2e-5;

k1 = 1e3;

k2 = 1e3;

k3 = 5e7;

OH = 1e-6;

%constants

t0=0;

tn=1000000;

% time iteration interval is set to a high value to reach steady state

Ca0 = 0.002292;

Cb0 = 0;

Cc0 = 0;

% initial concentrations

[tsoln, Csoln] = ode45(@Ozone,[t0 tn],[Ca0 Cb0 Cc0]);

plot(tsoln, Csoln)

xlabel('Time')

ylabel('Concentration')

legend('Cb')

function dC = Ozone(~, C)

dC = zeros(3,1);

dC(1) = -C(1)*(3*ki + kj*OH + 2*k2*sqrt((C(1)/(2*k3))*(2*ki+kj*OH)));

dC(2) = (2*ki)*C(1) - k1*C(1)*C(2) + k2*C(1)*C(3);

dC(3) = kj*C(1)*OH + k1*C(1)*C(2) - k2*C(1)*C(3) - (2*k3)*(C(3)^2);

end %dC(1) is evaluated by steady state approximation method

end

%HERE IS THE CODE, I POST IT IF YOU WANT TO CHECK

Answer 1

[Torsten]

There seems to be a singularity in the solution at t approximately 320.

ProjectOzone()

function ProjectOzone
clc
ki = 7e-3;
kj = 2e-5;
k1 = 1e3;
k2 = 1e3;
k3 = 5e7;
OH = 1e-6;
%constants
t0=0;
tn=321;%1000000;
% time iteration interval is set to a high value to reach steady state
Ca0 = 0.002292;
Cb0 = 0;
Cc0 = 0;
% initial concentrations
[tsoln, Csoln] = ode15s(@Ozone,[t0 tn],[Ca0 Cb0 Cc0]);
plot(tsoln, Csoln)
xlabel('Time')
ylabel('Concentration')
legend('Cb')
    function dC = Ozone(~, C)
        dC = zeros(3,1);
        dC(1) = -C(1)*(3*ki + kj*OH + 2*k2*sqrt((C(1)/(2*k3))*(2*ki+kj*OH)));
        dC(2) = (2*ki)*C(1) - k1*C(1)*C(2) + k2*C(1)*C(3);
        dC(3) = kj*C(1)*OH + k1*C(1)*C(2) - k2*C(1)*C(3) - (2*k3)*(C(3)^2);
    end %dC(1) is evaluated by steady state approximation method
end