# Resolvendo o sistema linear Au=v pelo método dos gradientes.

c=c(8, 1, 1, 0, 2, 20, 9, 2, 1, 1, 2, 1, 11,2, 0, 3, 2, 1, -8, 1, 2, 1, 2, 3, 10)
A=matrix(c,5,5,,byrow=TRUE); A   # Entrada da matriz A
b=c(1,2,3,4,5); b                # Entrada do vetor v

f<-function(u){A%*%u-b}          # Au=v quando f(u)=0.

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


gradg1<-function(w){u}
gradg2<-function(u){B%*%u-bb}

n=4
u=c(1,1,1,1,1)
                                 # Início do método
w=B%*%u-bb                        # w0
v=w                              # d1
alpha= t(w)%*%w /(t(B%*%w) %*%w); alpha=alpha[1]     # Escolha de h0
u=u-alpha*w                       # u1

for ( j in 1:n){
w=(B%*%u-bb)
beta=t(B%*%w) %*%v/(t(B%*%v) %*%v); beta=beta[1]
v=w-beta*v                   # Método de Euler
alpha= t(v) %*%w/(t(B%*%v) %*%v);alpha=alpha[1]
u=u-alpha*v                      # Método quase Euler (integral)
}
print(" Aproximação para u");u;v                               # Aproximação para u(20).
print(" Aproximação para f(u)") ;f(u)                         # Verificação do erro de aproximação.
print(" Número de iterações"); n