% Specify theta
theta = linspace(0,2*pi,301);

% Set growth factor
g = exp(i*theta);

% Calculate lambda dt location for each method and plot.

% 1st order
z1 = (g-1)./g;
plot(z1); hold on;

% Second order
z2 = (g.^2 - 4/3*g + 1/3)./(2/3*g.^2);
plot(real(z2),imag(z2),'r');

% Third order
z3 = (g.^3 - 18/11*g.^2 + 9/11*g - 2/11)./(6/11*g.^3);
plot(z3);

% Fourth order
z4 = (g.^4 - 48/25*g.^3 + 36/25*g.^2 - 16/25*g + 3/25)./(12/25*g.^4);
plot(z4);

axis('equal');
%axis([-3,3,-3,3]);
xlabel('Real \lambda \Delta t');
ylabel('Imag \lambda \Delta t');
grid on;