function Graph
addpath('/Users/yuhuazhu/Documents/MATLAB/VPFP_linear')
addpath('/Users/yuhuazhu/Documents/MATLAB/Reduced basis')
addpath('/Users/yuhuazhu/Documents/MATLAB/General')
addpath('/Users/yuhuazhu/Documents/MATLAB/Interpolation')

%% Setting up
Setup_mesh([0,2*pi,32],[-6,6,64],0.5,1,1);
global a; global b;  a= 1; b= 5; global dim; dim =1;
%get the exact solution
% %%%%%%%%%%%%only need to run once%%%%%%%%%%%%
approx_n = 50;
[z, f_ext, f_z] = Test_f(approx_n);
save('exact solu eq1', 'z', 'f_ext', 'f_z');
% %%%%%%%%%%%%only need to run once%%%%%%%%%%%%
load('exact solu eq1.mat');
[~,~,approx_n] = size(f_ext);

%% Get reduced basis
%%%%%%%%%%%%only need to run once%%%%%%%%%%%%
Nz = 101; Z = linspace(a,b,Nz)'; criteria = 10^(-8);
[f_0_select, Z_select, f_T_select] = RB_init(Z, criteria);
save('reduced basis eq1', 'f_0_select','Z_select','f_T_select')
%%%%%%%%%%%%only need to run once%%%%%%%%%%%%
load('reduced basis eq1.mat')


%% Greedy
[Nx,Nv,n]=size(f_T_select); N = Nx*Nv;
f_T_select = reshape(f_T_select, N,n);
f_0_select = reshape(f_0_select, N,n);
f_z = reshape(f_z, N,approx_n);
error_avg = zeros(n,1);
for i = 1:n
    f_approx = Proj_solu(f_z,f_0_select(:,1:i),f_T_select(:,1:i));
    f_approx = reshape(f_approx,Nx,Nv,approx_n);
    error = reshape(L2_xv(f_ext - f_approx),approx_n,1);
    error_avg(i) = sum(error)/approx_n;
end

%plot dim-error
semilogy(1:n, error_avg, '-r*', 'MarkerSize',8);
title('Error of Greedy')
xlabel('Dimension') % x-axis label
ylabel('Averaged Error') % y-axis label
hold on;


%% Polynomial interpollation on quadratic points
if n > 10
    n = 10;
end
error_avg = zeros(n,1);
for i = 1:n
    Z_select = Guadrature_pt(i);
    f_0 = ID_EQ(Z_select);
    [~,f_T_select] = VPFP_linear(f_0, Z_select);
    [f_approx,~] = Poly_intp(z,Z_select,reshape(f_T_select,N,i));
    f_approx = reshape(f_approx,Nx,Nv,approx_n);
    error = reshape(L2_xv(f_ext - f_approx),approx_n,1);
    error_avg(i) = sum(error)/approx_n;
end

%plot dim-error
semilogy(1:n, error_avg, '-k*', 'MarkerSize',8);
title('Error of Collocation')
xlabel('Dimension') % x-axis label
ylabel('Averaged Error') % y-axis label
legend('Greedy','Quad pts interpolation')

1;