clc;
clear;
close all;

mesh_plot = 1;
pause_plot = 0;
save_figures = 0;
print_titles = 0;

N = 1;
[x,y] = OptimalNodes2D(N); [r,s] = xytors(x,y);

x_axis_limits = [-1, 1]; y_axis_limits = [-1, 1];
% not very fine mesh
mesh_size_parameter_string = '0.25';
mesh_size_parameter_string_comma = strrep(mesh_size_parameter_string,'.',',');
% Gmsh - Delaunay
method_string = 'Gmsh - Delaunay';
algorithm_string = 'Gmsh_Delaunay';
[P,E] = MshToMatlab2D('Mesh_quality_analysis_h=0.25_Delaunay.msh');
MeshQualityAnalysis;
if mesh_plot == 1
    figure;
    trimesh(E',P(1,:),P(2,:),'Color','black');
    axis equal;
    xlim(x_axis_limits); ylim(y_axis_limits);
    if print_titles == 1
        title(sprintf('Mesh for %s with h = %s',method_string,mesh_size_parameter_string));
    end
end
if save_figures == 1
    filename_string = ['Trimesh_',algorithm_string,'_for_h=',mesh_size_parameter_string_comma];
    saveas(gcf,filename_string,'pdf');
    saveas(gcf,filename_string,'png');
end
if pause_plot == 1
    pause;
end
% Gmsh - MeshAdapt
method_string = 'Gmsh - MeshAdapt';
algorithm_string = 'Gmsh_MeshAdapt';
[P,E] = MshToMatlab2D('Mesh_quality_analysis_h=0.25_MeshAdapt.msh');
MeshQualityAnalysis;
if mesh_plot == 1
    figure;
    trimesh(E',P(1,:),P(2,:),'Color','black');
    axis equal;
    xlim(x_axis_limits); ylim(y_axis_limits);
    if print_titles == 1
        title(sprintf('Mesh for %s with h = %s',method_string,mesh_size_parameter_string));
    end
end
if save_figures == 1
    filename_string = ['Trimesh_',algorithm_string,'_for_h=',mesh_size_parameter_string_comma];
    saveas(gcf,filename_string,'pdf');
    saveas(gcf,filename_string,'png');
end
if pause_plot == 1
    pause;
end
% Gmsh - Frontal
method_string = 'Gmsh - Frontal';
algorithm_string = 'Gmsh_Frontal';
[P,E] = MshToMatlab2D('Mesh_quality_analysis_h=0.25_Frontal.msh');
MeshQualityAnalysis;
if mesh_plot == 1
    figure;
    trimesh(E',P(1,:),P(2,:),'Color','black');
    axis equal;
    xlim(x_axis_limits); ylim(y_axis_limits);
    if print_titles == 1
        title(sprintf('Mesh for %s with h = %s',method_string,mesh_size_parameter_string));
    end
end
if save_figures == 1
    filename_string = ['Trimesh_',algorithm_string,'_for_h=',mesh_size_parameter_string_comma];
    saveas(gcf,filename_string,'pdf');
    saveas(gcf,filename_string,'png');
end
if pause_plot == 1
    pause;
end
% Fluent Gambit - pouziva Hesthaven
method_string = 'FLUENT Gambit';
algorithm_string = 'FLUENT_Gambit';
[Nv, VX, VY, K, EToV] = MeshReaderGambit2D('Hesthaven_Maxwell_h=0.25.neu');
P = [VX; VY];
E = EToV';
MeshQualityAnalysis;
if mesh_plot == 1
    figure;
    trimesh(E',P(1,:),P(2,:),'Color','black');
    axis equal;
    xlim(x_axis_limits); ylim(y_axis_limits);
    if print_titles == 1
        title(sprintf('Mesh for %s with h = %s',method_string,mesh_size_parameter_string));
    end
end
if save_figures == 1
    filename_string = ['Trimesh_',algorithm_string,'_for_h=',mesh_size_parameter_string_comma];
    saveas(gcf,filename_string,'pdf');
    saveas(gcf,filename_string,'png');
end
if pause_plot == 1
    pause;
end


% finer mesh
mesh_size_parameter_string = '0.125';
mesh_size_parameter_string_comma = strrep(mesh_size_parameter_string,'.',',');
% Gmsh - Delaunay
method_string = 'Gmsh - Delaunay';
algorithm_string = 'Gmsh_Delaunay';
[P,E] = MshToMatlab2D('Mesh_quality_analysis_h=0.125_Delaunay.msh');
MeshQualityAnalysis;
if mesh_plot == 1
    figure;
    trimesh(E',P(1,:),P(2,:),'Color','black');
    axis equal;
    xlim(x_axis_limits); ylim(y_axis_limits);
    if print_titles == 1
        title(sprintf('Mesh for %s with h = %s',method_string,mesh_size_parameter_string));
    end
end
if save_figures == 1
    filename_string = ['Trimesh_',algorithm_string,'_for_h=',mesh_size_parameter_string_comma];
    saveas(gcf,filename_string,'pdf');
    saveas(gcf,filename_string,'png');
end
if pause_plot == 1
    pause;
end
% Gmsh - MeshAdapt
method_string = 'Gmsh - MeshAdapt';
algorithm_string = 'Gmsh_MeshAdapt';
[P,E] = MshToMatlab2D('Mesh_quality_analysis_h=0.125_MeshAdapt.msh');
MeshQualityAnalysis;
if mesh_plot == 1
    figure;
    trimesh(E',P(1,:),P(2,:),'Color','black');
    axis equal;
    xlim(x_axis_limits); ylim(y_axis_limits);
    if print_titles == 1
        title(sprintf('Mesh for %s with h = %s',method_string,mesh_size_parameter_string));
    end
end
if save_figures == 1
    filename_string = ['Trimesh_',algorithm_string,'_for_h=',mesh_size_parameter_string_comma];
    saveas(gcf,filename_string,'pdf');
    saveas(gcf,filename_string,'png');
end
if pause_plot == 1
    pause;
end
% Gmsh - Frontal
method_string = 'Gmsh - Frontal';
algorithm_string = 'Gmsh_Frontal';
[P,E] = MshToMatlab2D('Mesh_quality_analysis_h=0.125_Frontal.msh');
MeshQualityAnalysis;
if mesh_plot == 1
    figure;
    trimesh(E',P(1,:),P(2,:),'Color','black');
    axis equal;
    xlim(x_axis_limits); ylim(y_axis_limits);
    if print_titles == 1
        title(sprintf('Mesh for %s with h = %s',method_string,mesh_size_parameter_string));
    end
end
if save_figures == 1
    filename_string = ['Trimesh_',algorithm_string,'_for_h=',mesh_size_parameter_string_comma];
    saveas(gcf,filename_string,'pdf');
    saveas(gcf,filename_string,'png');
end
if pause_plot == 1
    pause;
end
% Fluent Gambit - pouziva Hesthaven
method_string = 'FLUENT Gambit';
algorithm_string = 'FLUENT_Gambit';
[Nv, VX, VY, K, EToV] = MeshReaderGambit2D('Hesthaven_Maxwell_h=0.125.neu');
P = [VX; VY];
E = EToV';
MeshQualityAnalysis;
if mesh_plot == 1
    figure;
    trimesh(E',P(1,:),P(2,:),'Color','black');
    axis equal;
    xlim(x_axis_limits); ylim(y_axis_limits);
    if print_titles == 1
        title(sprintf('Mesh for %s with h = %s',method_string,mesh_size_parameter_string));
    end
end
if save_figures == 1
    filename_string = ['Trimesh_',algorithm_string,'_for_h=',mesh_size_parameter_string_comma];
    saveas(gcf,filename_string,'pdf');
    saveas(gcf,filename_string,'png');
end
if pause_plot == 1
    pause;
end


% finest mesh
mesh_size_parameter_string = '0.0625';
mesh_size_parameter_string_comma = strrep(mesh_size_parameter_string,'.',',');
% Gmsh - Delaunay
method_string = 'Gmsh - Delaunay';
algorithm_string = 'Gmsh_Delaunay';
[P,E] = MshToMatlab2D('Mesh_quality_analysis_h=0.0625_Delaunay.msh');
MeshQualityAnalysis;
if mesh_plot == 1
    figure;
    trimesh(E',P(1,:),P(2,:),'Color','black');
    axis equal;
    xlim(x_axis_limits); ylim(y_axis_limits);
    if print_titles == 1
        title(sprintf('Mesh for %s with h = %s',method_string,mesh_size_parameter_string));
    end
end
if save_figures == 1
    filename_string = ['Trimesh_',algorithm_string,'_for_h=',mesh_size_parameter_string_comma];
    saveas(gcf,filename_string,'pdf');
    saveas(gcf,filename_string,'png');
end
if pause_plot == 1
    pause;
end
% Gmsh - MeshAdapt
method_string = 'Gmsh - MeshAdapt';
algorithm_string = 'Gmsh_MeshAdapt';
[P,E] = MshToMatlab2D('Mesh_quality_analysis_h=0.0625_MeshAdapt.msh');
MeshQualityAnalysis;
if mesh_plot == 1
    figure;
    trimesh(E',P(1,:),P(2,:),'Color','black');
    axis equal;
    xlim(x_axis_limits); ylim(y_axis_limits);
    if print_titles == 1
        title(sprintf('Mesh for %s with h = %s',method_string,mesh_size_parameter_string));
    end
end
if save_figures == 1
    filename_string = ['Trimesh_',algorithm_string,'_for_h=',mesh_size_parameter_string_comma];
    saveas(gcf,filename_string,'pdf');
    saveas(gcf,filename_string,'png');
end
if pause_plot == 1
    pause;
end
% Gmsh - Frontal
method_string = 'Gmsh - Frontal';
algorithm_string = 'Gmsh_Frontal';
[P,E] = MshToMatlab2D('Mesh_quality_analysis_h=0.0625_Frontal.msh');
MeshQualityAnalysis;
if mesh_plot == 1
    figure;
    trimesh(E',P(1,:),P(2,:),'Color','black');
    axis equal;
    xlim(x_axis_limits); ylim(y_axis_limits);
    if print_titles == 1
        title(sprintf('Mesh for %s with h = %s',method_string,mesh_size_parameter_string));
    end
end
if save_figures == 1
    filename_string = ['Trimesh_',algorithm_string,'_for_h=',mesh_size_parameter_string_comma];
    saveas(gcf,filename_string,'pdf');
    saveas(gcf,filename_string,'png');
end
if pause_plot == 1
    pause;
end
% Fluent Gambit - pouziva Hesthaven
method_string = 'FLUENT Gambit';
algorithm_string = 'FLUENT_Gambit';
[Nv, VX, VY, K, EToV] = MeshReaderGambit2D('Hesthaven_Maxwell_h=0.0625.neu');
P = [VX; VY];
E = EToV';
MeshQualityAnalysis;
if mesh_plot == 1
    figure;
    trimesh(E',P(1,:),P(2,:),'Color','black');
    axis equal;
    xlim(x_axis_limits); ylim(y_axis_limits);
    if print_titles == 1
        title(sprintf('Mesh for %s with h = %s',method_string,mesh_size_parameter_string));
    end
end
if save_figures == 1
    filename_string = ['Trimesh_',algorithm_string,'_for_h=',mesh_size_parameter_string_comma];
    saveas(gcf,filename_string,'pdf');
    saveas(gcf,filename_string,'png');
end
if pause_plot == 1
    pause;
end

%[Nv, VX, VY, K, EToV] = MeshReaderGambit2D('advection_2D_100.neu');
%P = [VX; VY];
%E = EToV';
%MeshQualityAnalysis;
%if mesh_plot == 1
%    figure;
%    trimesh(E',P(1,:),P(2,:),'Color','black');
%    axis equal;
%    xlim([-5, 15]); ylim([-5, 15]);
%end