ABCCHEM
ABCCHEM

Reputation: 35

My code doesnt seem to run past this line

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:
    Coord_List.append(r0)
    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)):
                break
            else:
                alpha = alpha*(1.5)
                continue
        else:
            alpha = alpha/2
            continue
    if abs(V0-V1) < precision:
        Coord_List.append(r1)
        iterations +=1
        break       
# 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?

Upvotes: 0

Views: 93

Answers (1)

Louis Matthijssen
Louis Matthijssen

Reputation: 566

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:
        break

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.

Upvotes: 3

Related Questions