# Resolvendo o sistema linear Au=v pelo método dos gradientes.
x=c(1,2,3,4,5)
y=c(2,3,-1,3,2)


A=matrix(0,5,5)
for (i in 1:5){A[,i]=x^(i-1)}

print("matriz A"); A   # Entrada da matriz A
v=y; print("Vetor v"); v                # Entrada do vetor v

B=t(A)%*%A; print("Matriz A*A");B; b=t(A)%*%v; print("vetor A*v"); b      # Adaptando o sistema para A*Au=A*v, ou Bu=b, em que B=A*A é positiva definida e b=A*v.

n=100000
u=c(0,0,0,0,0)
t=0*(1:(n+1))  # Tempo inicial
fu=0*t

for ( j in 1:n){
w=B%*%u-b                      
h= t(w)%*%w /(t(B%*%w) %*%w)     # Escolha de passo ótimo para o
u=u-h[1]*w                       # Método de Euler
t[j+1]=t[j]+h[1]
fua=t(A%*%u-v)%*%(A%*%u-v)
fu[j+1]=fua[1]
}
print("Aproximação para a solução de Au=v"); u          # Aproximação para u(t), sendo t a soma dos passos h.
print(" tempo t estimado"); t[n+1]   
print(" Teste para erro de aproximação");fu[n+1]                        # Verificação do erro de aproximação.
plot(t,fu,'l')