function fig10_21

% solve z'=Bz

n=6;
m=4;

S0=rand(m,m);
D0=diag([1,0.5,-0.5,-1]);
B=S0*D0*inv(S0);
[VB,DB]=eig(B);

for ik=1:5
    z=zeros(n,4);
    z(1,:)=[100 0 0 0];
    beta=VB\z(1,:)';
    k=10^(2-ik);
    t=0;
    for it=2:n
        t=t+k;
        for im=1:m
            z(it,:)=z(it,:)+beta(im)*exp(DB(im,im)*t)*VB(:,im)';
        end
    end

    XD=z(1:n-1,:);
    YD=z(2:n,:);

    % use the SVD
    [U,S,V]=svd(XD);
    SS=zeros(n-1,m);
    for i=1:m
        SS(i,i)=1/S(i,i);
    end
    P=YD'*U*SS*V';
    Ac=P;

    BB=(eye(m)-inv(P))/k;
    kk(ik)=k;
    norm_error(ik)=norm(B-BB,inf)/norm(B,inf)

    [VP,DP]=eig(P);
    [D0,I] = sort(diag(DP),'descend');
    D0;
    (1-1./D0)/k

end

figure(1)
clf
% get(gcf)
set(gcf,'Position', [25 1115 658 230])
subaxis(1,1,1,1,'MT',0.01,'MB',0.18,'MR',-0.02,'ML',0.05,'P',0.04)
co = [0 0 1;
    0 0.5 0;
    1 0 0;
    0 0.75 0.75;
    0.75 0 0.75;
    0.75 0.75 0;
    0.25 0.25 0.25];
set(groot,'defaultAxesColorOrder',co)
loglog(kk,norm_error,'--or','LineWidth',1.5,'MarkerSize',8)
hold on
box on
grid on
axis([1e-3 10 1e-4 300])
set(gca,'ytick',[1e-4 1e-2 1 100])
xlabel('k-axis')
ylabel('Relative Error')
%legend({' k = 1 ',' k = 2 ',' k = 3 '},'AutoUpdate','off','Location','SouthWest','FontSize',14,'FontWeight','bold')
set(gca,'FontSize',16,'FontWeight','bold')

%exportgraphics(gca,'/Users/mark/Desktop/modes.eps')


















