Homec4science
Diffusion SPClibs (feature/refactoring) Browse

spclibs/matlab/pde1d_eig_zmumps.m
c723281d03b6feature/refactoring

Last Change
Raw File

pde1d_eig_zmumps.m
View Options

%
% @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
c4science · Help