n=41
a=-5;b=5;x=0*1:n

for ( k in 1:n){  # Nós de Chebyshev reescalados.
x[k]=(a+b)/2+(b-a)*cos((2*k-1)/(2*n)*pi )/2}

y=1/(1+x^2)

n=length(x)
plot(x,y,col="red")

# Agora cálculos mostram que a função g(u)= |Au-y|^2 / 2
# em que A é a matriz e é um vetor, como montamos no que segue.
n=length(x)
A=matrix(0,n,n); 
for ( j in 1:n){A[,j]=x^(j-1)}

f<-function(u){A%*%u-y} # f(u)=Au-y

g<-function(u){  #  |Au-v|^2 / 2
w=A%*%u-y; t(w)%*%w/2}

# Agora é conhecido que o jacobiano de f(u) é A e o gradiente de g(u) é dado por A^*(Au-z).
#
# 
Jacf<-function(u){ A }   # função jacobiano de f(u).


gradg<-function(u){t(A)%*%f(u)}


#---- Início da rotin para o método de Euler.


ZeroEuler<-function(u0,m){
Y=u0
for ( i in 1:m){
B=Jacf(Y)
w=t(B)%*%f(Y)
bb=abs( sum( (( t(B)%*%B )%*%w )*w ) )
h=sum(w*w)/bb
if (bb!=0){Y=Y-h*w }
}
Y}

# Teste de aproximação de solução do problema de minimização. 
u0=0*x;print("Condicao inicial nula") 
m=1000
print("Teste de g(u0)");g(u0)                       # Teste de g(u0).

u=ZeroEuler(u0,m) ; print("Aproximacao"); u     # Aproximadamente u(t)
print("Teste de minimo aproximado");g(u)                        # Teste de minímo aproximado

# Definição da aproximação de f(u).

poli<-function(s){
    p=u[1]
    for ( i in 2:n){p=p+u[i]*s^(i-1)}
    return(p)
} 

curve(poli,-5,5,col="blue",ylim=c(0,1))
points(x,y,col="red")        # Teste de ajuste.

plot(x,y-poli(x),'l')

# Se quiser pode aumentar t e n ou escolher outra condição inicial 
# para se convencer do resultado obtido.