%% Runge Function Example
%
% Author(s): Paul R. Miles & Ralph C. Smith
%
% Date Created: April 10, 2019
%
% Exampled created using MATLAB2018a.  Certain features may not be
% compatible with previous versions of MATLAB.
%
% The Lagrange Polynomial toolbox is available for download on the File Exchange
% of MATLAB Central:
% <https://www.mathworks.com/matlabcentral/fileexchange/13151-lagrange-interpolator-polynomial?focused=5083678&tab=function>

% Setup workspace
close all; clear; clc;

% Add path to Lagrange Polynomial toolbox
% Note: I assume the uncompressed folder has been placed in the current
% working directory.
addpath('lagrangepoly');

%% High-Fidelity Model or "Data"
%
% $f = \frac{1}{1+25q^2}$
f = @(q) 1./(1+25*q.^2);

%% Grid Choice 1
% $q^j = -1 + (j-1)\frac{2}{M}$
M1 = 11;
j = 1:1:M1;
x1 = -1 + (j-1)*2/M1;
y1 = f(x1);
[P1, R1, S1] = lagrangepoly(x1, y1);

%% Grid choice 2: Clenshaw-Curtis
% $q^j = -\cos\Big(\frac{\pi(j-1)}{M-1}\Big)$
M2 = 10;
j = 1:1:M2;
x2 = -cos(pi*(j-1)/(M2-1));
y2 = f(x2);
[P2, R2, S2] = lagrangepoly(x2, y2);

%% Plot
% * L - Lagrange Polynomial
% * QP - Quadrature Points
% * HFM - High Fidelity Model
xx = -1:0.01:1;

% Plot Grid Choice 1
subplot(1, 2, 1);
hold on
h3 = plot(xx, f(xx), '-k', 'linewidth', 2);
h1 = plot(xx, polyval(P1, xx), '--r', 'linewidth', 2);
h2 = plot(x1, y1, 'or', 'linewidth', 2, 'markersize', 10);
hold off
xlabel('q');
ylabel('f');
legh = legend([h1, h2, h3], ...
    {'L', 'QP', 'HFM'});
set(legh, 'location', 'best');
box on
title(sprintf('Grid Choice 1: M = %i', M1));
set(gca, 'FontSize', 22);

% Plot Grid Choice 2
subplot(1, 2, 2);
hold on
h3 = plot(xx, f(xx), '-k', 'linewidth', 2);
h1 = plot(xx, polyval(P2, xx), '--b', 'linewidth', 2);
h2 = plot(x2, y2, 'sb', 'linewidth', 2, 'markersize', 10);
hold off
xlabel('q');
ylabel('f');
legh = legend([h1, h2, h3], ...
    {'L', 'QP', 'HFM'});
set(legh, 'location', 'best');
box on
title(sprintf('Grid Choice 2: M = %i\nClenshaw-Curtis', M2));
set(gca, 'FontSize', 22);
set_figure_dimensions(12, 6)

%% Support Functions
function set_figure_dimensions(width, height)
% adjust figure dimensions
set(gcf,'units','inches');
pos = get(gca, 'pos');
set(gcf, 'Position', [pos(1) pos(2) width height]);
end