%%% rate of convergence comparison

datasd1 = load('Xstoresd1');
Xstoresd1 = datasd1.Xstore;
datasd2 = load('Xstoresd2');
Xstoresd2 = datasd2.Xstore;
datasd3 = load('Xstoresd3');
Xstoresd3 = datasd3.Xstore;

N = length(Xstoresd1);

%true limit points
Xsd1 = zeros(4);
Xsd1(1) = 1;
Xsd2 = zeros(4);
Xsd2(1,1) = 0.6;
Xsd2(2,2) = 0.4;
Xsd3 = eye(4)/3;
Xsd3(4,4) = 0;

%rates of convergence
ratesd1 = zeros(N-1,1);
ratesd2 = zeros(N-1,1);
ratesd3 = zeros(N-1,1);
for i=1:N-1
    ratesd1(i) = norm(Xstoresd1(:,:,i+1) - Xsd1,'fro')/norm(Xstoresd1(:,:,i) - Xsd1,'fro');
    ratesd2(i) = norm(Xstoresd2(:,:,i+1) - Xsd2,'fro')/norm(Xstoresd2(:,:,i) - Xsd2,'fro');
    ratesd3(i) = norm(Xstoresd3(:,:,i+1) - Xsd3,'fro')/norm(Xstoresd3(:,:,i) - Xsd3,'fro');
end

figure(1)
plot(ratesd1(70:98),'k')
hold on
plot(ratesd2(70:98),'r')
plot(ratesd3(70:98),'b')
legend('sd1','sd2','sd3','Location','best')
hold off




