function[P] = lagrange(X, Y)
    n = length(X);
    x = poly(0,"x");
    P = 0;

    for i = 1:n, L=1;
        for j = [1:i-1,i+1:n]
            L = L*(x-X(j)) / (X(i)-X(j));
        end

        P = P + L*Y(i);
    end
endfunction

clear;
clc;

//x_in = input("Entre com o vetor x: ");
//y_in = input("Entre com o vetor y: ");
x_in = [0, 2, 4, 5];
y_in = [0, 4, 16, 25];
p = lagrange(x_in, y_in);

printf("Polinômio: ");
disp(p);


x = linspace(x_in(1), x_in($));
y = linspace(y_in(1), y_in($));

plot(x, y);