# Método de Gauss-Seidel

Sass<-function(A){ # Teste de convergêncial de Sassenfeld
n=length(A[1,])
b=matrix(1,n,1)
for ( i in 1:n){
b[i]=sum(abs(b[-c(i)])*abs(A[i,-c(i)]))/abs(A[i,i])
}
return(b)
}

Sassb<-function(A){ # Teste de convergêncial de Sassenfeld
n=length(A[1,])
cont=1:n
b=matrix(1,n,1)
gamma=0*A
gamma[1]=max(A[1,-c(1)]/A[1,1])
b[1]=gamma[1]
for ( i in 2:(n-1)){
s=sum(b[1:(i-1)]*abs(A[i,1:(i-1)]))
for (j in (i+1):n){
gamma[i,j]=(s+abs(A[i,j]))/abs(A[i,i])
}
b[i]=max(gamma[i,])
}
b[n]=sum(b[1:(n-1)]*abs(A[n,1:(n-1)]))/abs(A[n,n])
return(list(gamaij=gamma,gamma=b))
}

# Teste --------------------------------------------
C=matrix(c(25,2,1,2,10,2,1,1,1),3,3)
print("teste convergencia")
C 
Sass(C)
3/25;31/250; 3/25+62/250
print("Teste 1")
Sassb(C)

print("fim teste")
#------------------------------------------------------

GSeidel<-function(A,b,erro){  # Gauss-Seidel
m=length(b)
Erro=erro+1

c=0*v
B=0*A

for ( i in 1:m){
c[i]=b[i]/A[i,i]
B[i,]=A[i,]/A[i,i]
B[i,i]=0
}

x=0*b
cont=1
while (Erro>=erro){
y=x
mm=1:m
for (i in mm){
p=0
for ( j in mm[-c(i)]){
p=p+B[i,j]*x[j]
}
x[i]=c[i]-p
}
cont=cont+1
Erro=sum(abs(x-y))

if (cont>=10^6){Erro=0}}
x
}

# ------------------------------------------Um exemplo 

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

b=Sass(A);b;max(b)                # Teste de convergência teórica
b=Sassb(A);b;max(b$gamma) 
x=GSeidel(A,v,10^(-16))         # Gauss-Seidel
x                               # Solução aproximada     
A%*%x-v                         # Teste de erro na aproximação

# ------------------------------------------Um exemplo 

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

b=Sass(A);b;max(b)                # Teste de convergência teórica
b=Sassb(A);b;max(b$gamma) 
x=GSeidel(A,v,10^(-16))         # Gauss-Seidel
x                               # Solução aproximada     
A%*%x-v                         # Teste de erro na aproximação