Untitled

clear;
clc;
clf;

m=1
n=8
x_n = [-1.1, -0.83, -0.57, -0.28, -0.02, 0.16, 0.95, 1.23]

function y=f(x);
    y = exp(sin(x))
endfunction

function y=g(m,n,x);
    for j = 0:(2*m)
        s = 0
        for i = 0:n
            s = s + x(i+1)^j
        end
        y(j+1) = s
    end
endfunction

function y=rho(m,n,x);
    for j = 0:m
        s=0
        for i = 0:n
            s = s + f(x(i+1)) * x(i+1)^j
        end
        r(j+1) = s
    end
endfunction

function y=F(x,a,m);
    y = a(1)
    for i=1:m
        y(i+1) = y(i+1) + a(i+1)*x^(i+1)
    end
endfunction

function y =generate(lp,xt,a,m);
    for i=1:lp
        y(i) = F(xt(i),a,m);
    end
endfunction
    

y_n = f(x_n);
gw = g(m,n,x_n);
rw = rho(m,n,x_n);

select m
case 0 then
    a(1) = rw(1)/gw(1)
case 1 then
    l=-gw(3)*rw(1) + gw(2)*rw(2)
    m=gw(2)^2 - gw(1)*gw(3)
    a(1) = l/m
    
    l = gw(2)*rw(1) - gw(1)*rw(2)
    m = gw(2)^2 - gw(1)*gw(3)
    a(2) = l/m
case 2 then
    m = gw(3)^3 + gw(1)*gw(4)^2 + gw(2)^2*gw(5) - gw(3)*(2*gw(2)*gw(4) + gw(1)*gw(5))
    a(1) = (gw(4)^2*rw(1) - gw(3)*gw(5)*rw(2) + gw(3)^2*rw(3) - gw(4)*(gw(3)*rw(2) + gw(2)*rw(3)))/m
    a(2) = (gw(2)*gw(4)*rw(1) + gw(3)^2*rw(2) - gw(3)*(gw(4)*rw(1) + gw(2)*rw(3)) + gw(1)*(-gw(5)*rw(2) + gw(4)*rw(3))) / m
    a(3) = (gw(3)^2*rw(1) - gw(2)*gw(4)*rw(1) + gw(1)*gw(4)*rw(2) * gw(2)^2*rw(3) - gw(3)*(gw(2)*rw(2) + gw(1)*rw(3))) / m
end
lp=100
xt = linspace(x_n(1) - 0.1, x(n+1) + 0.1,lp)

Fn = generate(lp,xt,a,m)
plot(xt,Fn)