Homec4science
Diffusion SPClibs (feature/refactoring) Browse

spclibs/matlab/gs_fe.m
c723281d03b6feature/refactoring

Last Change
Raw File

gs_fe.m
View Options

%
% @file gs_fe.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>
%
clear all
N=1024;
L=100;
n0=1;
Narr=[128 256 512 1024];
for i=1:length(Narr);
N=Narr(i);
h=L/N; h2=h^2;
title_str=sprintf(['N=%d, L=%.1f, n0=%.1f'], N, L, n0);
v0=-2.+2.0*n0*h2/6.0;
u0=1.+n0*h2/6.0;
dom(i)= 2*abs(u0)/abs(v0);
v=v0.*ones(N,1);
u=u0*ones(N-1,1);
mata=diag(u,1) + diag(u,-1) + diag(v);
matl= tril(mata,0); % D+L
lambda = eig(-triu(mata,1),matl);
rho(i)=max(abs(lambda));
% $$$ figure
% $$$ plot(lambda,'o')
% $$$ xlabel('Real eigenvalue'), ylabel('Imag eigenvalue')
% $$$ grid on
% $$$ axis equal
% $$$ title_str=sprintf(['N=%d, L=%.1f, n0=%.1f, Spec. Radius=%.4f, dom=%.4f'], N, L, ...
% $$$ n0, rho(i),dom(i));
% $$$ title(title_str);
fprintf(1, 'Spectral Radius of GS relaxation matrix = %.4f, dom=%.4f\n', ...
rho(i),dom(i))
end
figure
plot(Narr, rho,'o-')
xlabel('N'); ylabel('GS spectral radius')
title_str=sprintf(['FE scheme, L=%.1f, n0=%.1f'], L, n0);
title(title_str);
grid on
c4science · Help