%% Logistic Method - Predicator-Corrector
%
% We solve the logistic problem
%
% $$x'=(1-x)*x, x(0)=2$$
%
% using PECECE Predicator-Corrector (with ab1 and am2), Euler, and ode23

a = 1;
b = 1;
c = 0;
logistic = @(t,x)(a-b*x)*x-c;

figure;
axis([0 3 1 2]);
hold on;
xlabel('t');
ylabel('x');
[t,x]=pred_corr(logistic,0,3, 2, .1);
plot(t,x(:,1),'ro','DisplayName','Predictor-Corrector');

[t,x]=ode23(logistic,[0,3], 2);
plot(t,x(:,1),'k-','DisplayName','ode23');

legend('Location', 'NorthEast');


%% Linear Oscillator - Predicator-Corrector
%
% We solve linear oscillator with $a = 0$, $b = 9$, $c = 10$, $\omega = 1$
% and $x_0=(2,1)^T$ using PECECE Predicator-Corrector (with ab2 and am2),
% Euler, and ode23.

a = 0;
b = 9;
c = 10;
omega = 1;

oscillator = @(t,x) [0 1; -b -a]*x + [0; c*cos(omega*t)];

figure
axis([0 10 -4 4]);
hold on;
xlabel('t');
ylabel('x_1');
[t,x]=pred_corr(oscillator,0,10, [2;1], .1);
plot(t,x(:,1),'ro','DisplayName','Predictor-Corrector');

% compare with  ode23
[t,x]=ode23(oscillator,[0,10], [2;1]);
plot(t,x(:,1),'k-','DisplayName','ode23');

legend('Location', 'NorthEast');


%% Satellite ODE - Predicator-Corrector
%
% We solve satellite ODE using PECECE Predicator-Corrector (with ab2 and am2),
% Euler, and ode23.

figure
hold on;
xlabel('t');
ylabel('x_1');
zlabel('x_2');
xlim([0 7]);
ylim([-5 5]);
zlim([-5 5]);
T = 6.19216933131963970674;
x0 = [1.2; 0; 0; -1.04935750983031990726];
[t,x]=pred_corr(@sat_ode,0,T, x0, 1e-3);
plot3(t, x(:,1), x(:,2),'ro','DisplayName','Predictor-Corrector');

% compare with  ode23
[t,x]=ode23(@sat_ode,[0,T], x0);
plot3(t, x(:,1), x(:,2),'k-','DisplayName','ode23');

view(3);
legend('Location', 'NorthEast');

%% Convergence Analysis
% We perform convergence analysis for PECECE Predicator-Corrector
% (with ab1 and am2)

fprintf('Predicator Corrector (PECECE, ab1, am2): = ');
conv_analysis(@pred_corr);
