/*function res = PolDeg(x)
    res = 0
    for i = 1:length(x)
        if x(length(x)-i+1)~=0 then
            res = length(x) - i
            break
        end
    end

endfunction
function res = PolS(a,b)

    sizeA = PolDeg(a)
    sizeB = PolDeg(b)
    sizeRes = sizeA
    if sizeB > sizeA then
    c = b
    b = a
    a = c
    sizeRes = sizeB
    sizeA = PolDeg(a)
        sizeB = PolDeg(b)
    end
    sizeRes = sizeRes + 1
    res = [1:sizeRes] * 0
    for i =1:sizeRes
    res(i) = res(i) + a(i)
    if i <= sizeB+1
        res(i) = res(i) + b(i)
    end
    end


endfunction

function res = PolM(a,b)
    sizeA = PolDeg(a)
    sizeB = PolDeg(b)
    sizeRes = sizeA + sizeB
    res = [1:sizeRes+1] * 0
    for i =1:sizeRes+1
        for j = 1 :i
            if (j <= sizeA+1) & (i-j+1<=sizeB+1) then
                res(i) =res(i)+ b(i-j+1) * a(j)
            end
        end
    end
endfunction

function res = EvaluatePol1(a,x)
    sizeA = PolDeg(a)
    res = 0
    r = 1
    for i = 1: sizeA+1
        res = res + a(i) * r
        r = r * x
    end
endfunction
function res = EvaluatePol(a,x)
    for i = 1: length(x)
        x(i) = EvaluatePol1(a,x(i))
    end
    res = x
endfunction

function res = Lagrange(x,y)
        res = [0]
        for i = 1:length(x)
            yi = [y(i)]
            for j = 1:length(x)
                if i ~= j then
                    yi = PolM(yi ,( [-x(j)/ (x(i) - x(j)),1/ (x(i) - x(j))]  ))
                end
            end
            res = PolS(res ,yi)
        end
    endfunction
function res = DivDiff(x,y,n)
        m = rand(length(x),n) * 0
        for i =1:length(x)
            m(1,i) = x(i)
            m(2,i) = y(i)
        end
        for i=3:n+2
            for j=1:length(x)+2-i
                m(i,j) = (m(i-1,j) - m(i-1,j+1))/(x(j) - x(j+i-2))
            end
        end
    res = m
    endfunction

    function res = Newton(x,y)
        m= DivDiff(x,y,length(x))
        res =[y(1)]
        yi = [1]
        for i = 1:length(x)-1
            c = m(2+i,1)
            yi = PolM(yi ,[-x(i),1])
            res = PolS(res, PolM(yi,[c]) )
        end
    endfunction
/////////
*/
n=25
from=0
to=6
x=[from:(to-from)/n:to]
xc=[0]
for i=1:n+1
    xc(i)=(to+from)/2+((to-from)/2)*cos(((2*i-1)*%pi)/(2*(n+1)))
end
my_y=[0]
//for i=1:n+1
    my_y=exp(sin(x))
    //disp(my_y(i))
//end
my_yc=[0]
for i=1:n+1
    my_yc(i)=exp(sin(xc(i)))
    //disp(my_y(i))
end
function DD = DivDiff(x,y,n)
        m = rand(length(x),n) * 0
        for i =1:length(x)
            m(1,i) = x(i)
            m(2,i) = y(i)
        end
        for i=3:n+2
            for j=1:length(x)+2-i
                m(i,j) = (m(i-1,j) - m(i-1,j+1))/(x(j) - x(j+i-2))
            end
        end
    DD = m
endfunction
function Lagrange = Lag(x, y)
    a=1
    Lagrange=poly([0],"x","c")
    for i=1:length(x)
        c=1
        Li=poly([a],"x","c")
        for j=1:length(x)
            if i~=j then
                li=poly([-x(j),1],"x", "c")/(x(i)-x(j))
                li=li
                Li=li*Li
                //disp(Li)
            end
        end
        disp(y(i))
        Lagrange=y(i)*Li+Lagrange
    end
    Lagrange=Lagrange/a
    //disp(Lagrange)
endfunction
function Newtoon = Newt(x, y)
    N=poly([1], "x", "c")
    m= DivDiff(x,my_y,length(x))
    disp(m)
    Newtoon=poly([y(1)],"x","c")
    for i=1:length(x)-1
        Ni=poly(x(i), "x", "r")
        N=N*Ni
        Newtoon=Newtoon+N*m(2+i,1)
    end
endfunction

Newtoon=Newt(x,my_y)
PolLag = Lag(x, my_y)
PolLagCh = Lag(xc, my_yc)
_nodes=[from:1/1000:to]
value_calc = horner(PolLag,_nodes)
value_calcCh = horner(PolLagCh,_nodes)
//plot(_nodes,value_calc, "red")
//plot(_nodes,value_calcCh, "green")
plot(_nodes,exp(sin(_nodes)), "blue")
value_calcN = horner(Newtoon,_nodes)
plot(_nodes,value_calcN, "red")
///////