Skip to content

Instantly share code, notes, and snippets.

@angellicadearaujo
Created October 24, 2016 13:40
Show Gist options
  • Select an option

  • Save angellicadearaujo/9770038bf80f69f5e05ff3146cb73fb2 to your computer and use it in GitHub Desktop.

Select an option

Save angellicadearaujo/9770038bf80f69f5e05ff3146cb73fb2 to your computer and use it in GitHub Desktop.
GMRS condicionado
import math
# @description Multiplies a matrix by a vector (one column matrix)
def multiply(A, x):
m=len(A)
n=len(x)
y=[0 for i in range(m)]
for i in range(m):
for j in range(n):
# print('Multiplicando A[' + str(i) + '][' + str(j) + '] por x['+ str(i) + ']')
y[i] += A[i][j]*x[j]
return y
# @description Subtract y from x (x - y) where both x and y are vectors
def subtract(x, y):
n=len(y)
z=[0 for i in range(n)]
for i in range(n):
z[i]= x[i] - y[i]
return z
# @description Performs the internal product between two vectors
def internProduct(u,v):
l=len(u)
p=0
for i in range(l):
p += u[i]*v[i]
return p
# @description Multiplies a scalar by a vector
def multiplyScalarByVector(s, v):
l=len(v)
x=[1 for i in range(l)]
for i in range(l):
x[i]=s*v[i]
return x
# @description Subtract y from x (x - y) where both x and y are vectors
def subtract(x, y):
n=len(y)
z=[0 for i in range(n)]
for i in range(n):
z[i]= x[i] - y[i]
return z
# @description Get the lenght of the vector
def norm(x):
s=0;
l=len(x)
for i in range(l):
s += x[i]**2
return math.sqrt(s)
def GMRES(A, b, guess, jmax):
x=list(guess)
w=[]
j=0
while j < jmax:
x= multiply(A,x)
if j==0:
x=subtract(b,x)
# implementar pre cond aqui (right)
for i in range(j):
t=internProduct(x,w)
tw=multiplyScalarByVector(t,w)
x=subtract(x, tw)
alpha=norm(x)
x=multiplyScalarByVector(1/float(alpha),x)
w=list(x)
# implementar pre cond aqui também (left)
j += 1
return x
print(GMRES([[3,-1,1], [-1,3,-1],[1,-1,3]],[-1,7,-7],[0,0,0], 10))
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment