%
% @file pde1d_eig_zmumps.m
%
% @brief
%
% @copyright
% Copyright (©) 2021 EPFL (Ecole Polytechnique Fédérale de Lausanne)
% SPC (Swiss Plasma Center)
%
% SPClibs is free software: you can redistribute it and/or modify it under
% the terms of the GNU Lesser General Public License as published by the Free
% Software Foundation, either version 3 of the License, or (at your option)
% any later version.
%
% SPClibs is distributed in the hope that it will be useful, but WITHOUT ANY
% WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
% FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.
%
% You should have received a copy of the GNU Lesser General Public License
% along with this program. If not, see <https://www.gnu.org/licenses/>.
%
% @author
% (in alphabetical order)
% @author Trach-Minh Tran <trach-minh.tran@epfl.ch>
%
file='pde1d_eig_zmumps.h5';

%%% Read mumps matrix and convert to Matlab sparse format
[mata,diag_ele]=zmumps_mat(file,'/MAT');
n=size(mata,1);

%spy(mata,12);
arpack_vals=h5Complex_ll(file,'/eig_vals');
nev=size(arpack_vals,1);
arpack_vecs = h5Complex_ll(file, '/eig_vecs');

%%% Compute eigen values and vectors, using EIGS
[V,D,FLAG]=eigs(mata,nev,'SM');
[eigs_vals,perm]=sort(diag(D));
eigs_vecs=V(:,perm);

fprintf('Eigenvalues from Arpack and Matlab eigs\n');
for i=1:nev
   fprintf('%i (%.5e %.5e), %.5e\n',i,real(arpack_vals(i)),imag(arpack_vals(i)),eigs_vals(i));
end
fprintf('Norm of difference %.3e\n', norm(arpack_vals-eigs_vals,Inf));

%%% Renormalize EIGS %%%
fprintf('\n\nDiff of Eigenvectors from Arpack and Matlab eigs\n');
for i=1:nev
    nrm=arpack_vecs(1,i);
    eigs_vecs(:,i) = (eigs_vecs(:,i)./eigs_vecs(1,i)).*nrm;
    diff_vecs = norm(arpack_vecs(:,i)-eigs_vecs(:,i),Inf);
    fprintf('%i %10.3e\n', i, diff_vecs);
end