%% .. Polymerization (RHS)
% ...  
% ... 13.03.2014:
% ...
% ... Milan Bernauer, ICT, Prague 2014
%
function dydt = dcdt(t,c)
global kd kp kt n nm

% c(1)     ... cI
% c(2)     ... cM
% c(3-n+2) ... cP(n)

% ... Substitution
A = (2.0 .* kd .* c(1) ./ (kp .* c(2) + sqrt(2.0.* kt .* kd .* c(1)))).^2;
B = @(z)((kp.*c(2) ./ (kp .* c(2) + sqrt(2.0.* kt .* kd .* c(1)))).^(z-2));

% --- RHS 
dydt(1) = -kd .* c(1);
dydt(2) = -kp .* c(2) .* sqrt(2.0 .* kd .* c(1) ./ kt);
 
% --- Using loop
%for i = 1:n
%    dydt(i+2) = (nm(i) - 1) .* kt .*A .* B(nm(i)) ./ 2.0;
%end

% --- Without loop (faster execution in Matlab)
dydt(3:n+2) =  (nm(1:n) - 1) .* kt .*A .* B(nm(1:n)) ./ 2.0;
dydt = dydt';

end