I want to implement a conjugate gradient algorithm in python. However, when I run the code no results are printed/displayed. When I stop the program python shows me this:

def V(r):
return (r[0]-r[1])**4+2*r[0]**2+r[1]**2-r[0]+2*r[1]

def g(r):
return np.array([4*(r[0]-r[1])**3+4*r[0]-1, -4*(r[0]-r[1])**3+2*r[1]+2])  

def Conj_Grad(r):
iterations = 0
maxiterations = 1000
precision = 10**-10
a = 10**-4
b = 0.9
alpha = 0.09
r0 = r
g0 = g(r)
p0 = -g(r)
V0 = V(r)
Coord_List= []
S = np.dot(p0,g0)
while iterations < maxiterations:
    while True:
        r1 = r0 + alpha*p0 
        g1 = g(r1) 
        V1 = V(r1)
        if V1 <= V0 + a*alpha*S:
            g1 = g(r1)
            if abs(np.dot(g1,p0)) <= b*abs(np.dot(p0,g0)):
                alpha = alpha*(1.5)
            alpha = alpha/2
    if abs(V0-V1) < precision:
        iterations +=1
# Generating beta by Fletcher Reeves
    beta = np.dot(g1,g1)/np.dot(g0,g0)        
    p0 = -g1 + beta*p0
    g0 = g1
    V0 = V1
    iterations +=1
return (Coord_List, i, V1-V0)
<ipython-input-57-f26ac2ab8296> in Conj_Grad(r)
     29         Coord_List.append(r0)
     30         while True:
---> 31             r1 = r0 + alpha*p0
     32             g1 = g(r1)
     33             V1 = V(r1)

When I input an r, for example [1,1]; no matter how long or short I let it run it always seems to be stuck at line 31. What is wrong with this code?

1 Answers

Louis Matthijssen On

Your while loop will run forever if you don't exit it. The 3 lines in your while loop will be executed very fast so when you pause your code there's a big chance it will pause on that line, however that line may have been successfully executed thousands of times already.

You should include some condition to stop the while loop when you're done:

while True:
    # Computations
    if enough_iterations:

Edit: also check the indentation of your if statement. It should be "in" the while loop, but it's not at the moment. So probably the first 3 lines of your while loop are being executed.