function A = Jpop(x,t,pop)
% Jacobian matrix
a11 = sum(exp(2*x(2,1)*t));
a12 = sum(2*x(1,1).*t.*exp(2*x(2,1)*t) - pop.*t.*exp(x(2,1)*t));
a21 = a12;
a22 = sum(2*x(1,1)^2*(t.^2).*exp(2*x(2,1)*t) - ...
    x(1,1).*(t.^2).*exp(x(2,1)*t).*pop);
A = [a11,a12;a21,a22];
end

