function Y = hw2a(k, initL, initR, initLR, initLD, initD)
f
%The equation for concentration of L to be solved it:
%
%dL(t)/dt = -k1*L(t)*R(t) + k2*LR(t) - k5*L(t)*D(t) + k6*LD(t);
%
%The equation for concentration of R to be solved it:
%
%dR(t)/dt = -k1*L(t)*R(t) + k2*LR(t) - k7*L(t)*R(t) + k8*L(t)*D(t);
%
%The equation for concentration of LR to be solved it:
%
%dLR(t)/dt = k1*L(t)*R(t) + k2*LR(t) - k7*L(t)*R(t) + k8*L(t)*D(t);
%
%The equation for concentration of LD to be solved it:
%
%dLD(t)/dt = k3*LR(t) - k4*LD(t) + k5*L(t)*D(t) - k6*LD(t)
%
%The equation for concentration of D to be solved it:
%
%dD(t)/dt = k7*L(t)*R(t) - K8*L(t)*D(t) + k6*LD(t) - k5*L(t)*D(t)
%
%----------------Set up time interval--------------------
%
h = 0.1;
h
%----------------Set up duration-------------------------
%
t = 0:h:100;
t
%----------------Initialize Vector Y---------------------
%
Y = zeros(5,1);
Y
%---------------Initial conditions-----------------------
%
Y(1,1) = initL;
Y(2,1) = initR;
Y(3,1) = initLR;
Y(4,1) = initLD;
Y(5,1) = initD;
Y
%---------------RK2 Method's Loop------------------------
%
for i=1:length(t)-1,
%Set up "for" loop.
k1L = h * (-k(1,1)*Y(1,i)*Y(2,i) + k(2,1)*Y(3,i) - k(5,1)*Y(1,i)*Y(5,i) +
k(6,1)*Y(4,1));
k2L = h * (-k(1,1)*(Y(1,i)+k1L/2)*(Y(2,i)+k1L/2) + k(2,1)*(Y(3,i)+k1L/2) -
k(5,1)*(Y(1,i)+k1L/2)*(Y(5,i)+k1L/2) + k(6,1)*(Y(4,1)+k1L/2));
Y(1,i+1) = Y(1,i) + k2L;
k1R = h * (-k(1,1)*Y(1,i)*Y(2,i) + k(2,1)*Y(3,i) - k(7,1)*Y(1,i)*Y(2,i) +