Adams-Bashforth-Euler格式的稳定域。

figure(1), clf, hold on

    % Plot of stability domain
theta = 0:0.01:1; theta = theta*pi;
z = exp(2*i*theta) - exp(i*theta); z = z./(1.5*exp(i*theta) - 0.5); plot(z)

    % Plot of oval
b = 0.5^0.25; a = 0.5; z =  - a*(1-cos(theta)) + i*b*(sin(theta)).^0.5;
plot(z,'--')

    % Plot of half ellipse
a = 1.0; b = sqrt(2/3); theta = -pi/2 +0.5*theta;
z = -a*cos(theta) - i*b*sin(theta); plot(z,'.')

figure(1), clf, hold on

    % Plot of stability domain
theta = 0:0.01:1; theta = theta*pi;
z = exp(3*i*theta) - exp(2*i*theta);
z = 12*z./(23*exp(2*i*theta) - 16*exp(i*theta) + 5); plot(z)

    % Plot of half ellipse
a =6/11;b = (72/11)*sqrt(2/235); theta = pi/2 + 0.5*theta;
z = a*cos(theta) + i*b*sin(theta); plot(z,'--');

z = [0  0.8*i]; plot(z)%   Plot of imaginary axis

figure(1), clf, hold on

    % Plot of stability domain
theta = 0:0.005:0.8; theta = theta*pi; c = (2-cos(theta)).*(1+cos(theta));
a = 0.5*(1-cos(theta)).^2./c; b = - 2*sin(theta)./c; z = -(a+i*b); plot(z)

    % Plot of oval
a = 0.5; b = 0.75^0.25; theta = 0:0.002:1; theta = theta*pi;
z =  - a*(1-cos(theta)) + i*b*(sin(theta)).^0.5; plot(z,'--')

    % Plot of parabola
yy = 0:0.02:2; par = -0.75*(yy).^2; plot(par,yy,'.')

figure(1), clf, hold on

    % Plot of stability domain
theta = 0:0.01:1; theta = theta*pi;
z = - 2*(1-cos(theta)).^2./(3-cos(theta)) + i*2*sin(theta)./(3-cos(theta));
plot(z)

    % Stability domain of second order Adams-Bashforth scheme
z = exp(2*i*theta) - exp(i*theta); z = z./(1.5*exp(i*theta) - 0.5); plot(z,'-.')

    % Plot of oval
b = 1; a = 1; z =  - a*(1-cos(theta)) + i*b*(sin(theta)).^0.5; plot(z,'--')

    % Plot of half ellipse
a = 2.0; b = sqrt(1/3); theta = -pi/2 +0.5*theta;
z = - a*cos(theta) - i*b*sin(theta); plot(z,'.')

figure(1), clf

    % Plot of stability domain
theta = 0:0.01:1; theta = theta*pi;
z = 0.5*exp(-2*i*theta) - 2*exp(-i*theta) +3/2; plot(z)

figure(1), clf, hold on

    % Plot of stability domain
theta = 0:0.01:0.5; theta = theta*pi;
a = (6*cos(theta)-9)./(2*cos(2*theta) - 4*cos(theta)) - 3/2;
b = (sin(theta).*(2-cos(theta)))./(cos(2*theta) - 2*cos(theta));
z = - a - i*b; plot(z)

    % Plot of oval
a = 1; b =(2*a/4)^(1/4); theta = 0:0.002:1; theta = theta*pi;
z =  - a*(1-cos(theta)) + i*b*(sin(theta)).^0.5; plot(z,'--')

    % Plot of parabola
yy = 0:0.02:1.8; b = 1; par = -(yy/b).^2; plot(par,yy,'.')