Question

Write a python function to compute the solution to a linear system using the conjugate gradiant algorithm.

Here is the sudo code, but I can't figure out what T is:

1x(0)-initial guess 4for k-0,1,2, 5 α(k) 6 X(k+1)-X(k)+α(k)s(k)

Answer 1

Please find below the python code to implement the above pseudocode along with necessary details in comment. Also, note that 'T' in the above code refers to the transpose of a matrix.

import numpy as np

# define the matrices
A=np.array([[3,2],[2,6]]) # should be a positive definite matrix
b=np.array([[2],[-8]])

# initilaize the method
x0=np.array([[0],[1]])
r0=b-(A.dot(x0))
s0=r0
imax=30
iter=1
diffx=100

# iteratively update the solution
while (iter<=imax) and (diffx>0.01):
alpha=((np.transpose(r0).dot(r0))/((np.transpose(s0).dot(A)).dot(s0)))
x=x0+alpha*s0
r=r0-alpha*(A.dot(s0))
beta=(np.transpose(r).dot(r))/(np.transpose(r0).dot(r0))
s=r+beta*(s0)
r0=r
s0=s
iter=iter+1
diffx=np.linalg.norm(np.absolute(x-x0))
x0=x

  
print(x)

====================== SCREENSHOT OF CODE

1 import numpy as np 3# define the matrices 4A-np. array([ [3,2], [2,6]]) # should be a positive definite matrix 5 b-np.array ([[2], -8]]) 7 # initilaize the method 8 x0=np. array([ [e], [1] ]! 9 ro-b-(A.dot(xe)) 10 se-re 11 imax-30 12 iter=1 13 diffx-100 15# iteratively update the solution 16 while (iter<-imax) and (diffx0.01): 17 alpha-((np.transpose (r0).dot (r®))/((np.transpose (s0).dot(A)).dot (s0))) x-xe+alpha s0 19 r ro-alpha (A.dot (s0)) 20 beta (np.transpose(r).dot(r))/(np.transpose(re).dot (r®)) 21 r+beta(s0) 23 24 25 26 27 ro-r s@=s iter=iter+1 diffx-np.linalg.norm(np.absolute(x-xe)) 29 print (x)

=========================== SAMPLE OUTPUT

L 2.1 2 In [32]: