Homec4science
Diffusion SPClibs (master) Browse

spclibs/matlab/pde2d_nh.m
253e659ccb68master

Last Change
Raw File

pde2d_nh.m
View Options

%
% @file pde2d_nh.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='pde2d_nh.h5';
%
% Get data from data sets
%
r=hdf5read(file,'/xgrid');
t=hdf5read(file,'/ygrid');
sol=hdf5read(file,'/sol')';
solexact=hdf5read(file,'/solana')';
err=hdf5read(file,'/errors')';
errx=hdf5read(file,'/errors_x')';
erry=hdf5read(file,'/errors_y')';
solr=hdf5read(file,'/derivx')';
solt=hdf5read(file,'/derivy')';
%
% Attributes
%
NR=hdf5read(file,'/NX'); NTH=hdf5read(file,'/NY');
NIDBAS1=hdf5read(file,'/NIDBAS1');
NIDBAS2=hdf5read(file,'/NIDBAS2');
MBESS=hdf5read(file,'/MBESS');
LABEL=sprintf('nr = %d, ntheta = %d, nidbas = (%d,%d), mbess = %d', NR, NTH, ...
NIDBAS1, NIDBAS2, MBESS);
attr=hdf5read(file,'/xgrid/title'); x_ann=attr.Data;
attr=hdf5read(file,'/ygrid/title'); y_ann=attr.Data;
attr=hdf5read(file,'/sol/title'); sol_ann=attr.Data;
attr=hdf5read(file,'/solana/title'); solexact_ann=attr.Data;
attr=hdf5read(file,'/errors/title');err_ann=attr.Data;
[R,T]=meshgrid(r,t);
x = R.*cos(T); y= R.*sin(T);
solx = cos(T).*solr - sin(T)./R.*solt;
soly = sin(T).*solr + cos(T)./R.*solt;
figure
subplot(221)
pcolor(double(r),double(t),double(sol));
shading interp
hold on, quiver(r,t,solr,solt)
xlabel(x_ann); ylabel(y_ann)
title(LABEL)
colorbar
subplot(222)
pcolor(double(x),double(y),double(sol))
shading interp
hold on, quiver(x,y,solx,soly)
hold off, axis image
xlabel('X'); ylabel('Y')
title('X-Y plane')
colorbar
subplot(223)
surfc(double(x),double(y),double(sol))
xlabel('X'); ylabel('Y');
title(sol_ann)
subplot(224)
surfc(double(x),double(y),double(err))
xlabel('X'); ylabel('Y');
title(err_ann)
figure
subplot(311)
plot(t,err(:,NR+1),'o-')
xlabel('\theta'); ylabel('Error on solution')
grid on
title('Error on Boundary r=1');
subplot(312)
plot(t,errx(:,NR+1),'o-')
xlabel('\theta'); ylabel('Error on x-derivative')
grid on
subplot(313)
plot(t,erry(:,NR+1),'o-')
xlabel('\theta'); ylabel('Error on y-derivative')
grid on
c4science · Help