%% Copyright (c) 2025 Duong Ngoc Son <son@duongnotes.com>

%% THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
%% "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
%% LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
%% A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
%% HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
%% SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
%% LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
%% DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
%% THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
%% (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
%% OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.

%' @title Lagrange Interpolation Formula
%'
%' @description
%' Compute the value of the interpolating polynomial at a point using 
%' Lagrange Interpolation formula.

function w = LagrangeInterpolation(x, y, xi)

	% Tham số: x và y là hai dãy dữ liệu cần nội suy, xi là điểm cần
	% tính giá trị của đa thức nội suy.
	% Kết quả: w là giá trị của đa thức nội suy tại x = xi.
	
    % Tính độ dài của các dãy x và y

	n = length(x);
	m = length(y);
	
	if (m ~= n)
		t = NaN;	% Hai dãy x và y khác độ dài
	elseif (n == 1)
		t = 1;
	else
		t = ones(1, n)
		for i = 1:1:n
			for j = 1:1:n
				if j ~= i 
					p = (xi - x(j)) / (x(i) - x(j));
					t(i) = t(i) * p;
				end
			end
		end
	end

	w = 0;
	for i = 1:1:n
		w = w + y(i) * t(i);
	end
end

% Ví dụ áp dụng: x = (0, 2, 4, 5), y = (-2, 8, 115, 150)
% p(x) = -2 - (2911/60) * x + (1363/40) * x^2 - (439/120) * x^3

x = [0 2 4 5];
y = [-2 8 115 150];
xi = 1.2;

% Tính giá trị của đa thức nội suy đối với ví dụ trên

w = LagrangeInterpolation(x, y, xi);

% In ra giá trị của đa thức nội suy tại x = xi.

fprintf('%d\n', w);
