1

I'm trying to figure out why the above error came up in my code. The program is suppose to use the Newton Raphson Method to find the position of links in a 4 bar linkage given two input parameters.

The error is happening in this line

g1 = L1 * np.cos(theta) + L2 * np.cos(alpha) - L3

Thanks in advance for your help.

import numpy as np

L1=1
L2=1.5 * L1
theta = 40 * np.pi / 180

#initial guesses

L3 = 1.5
alpha = 30 * np.pi / 180
epsilon = 1
n = 0

while epsilon > 0.0001:
    g1 = L1 * np.cos(theta) + L2 * np.cos(alpha) - L3
    dg1dalpha = -L2 * np.sin(alpha)
    dg1dL3 = -1;

    g2 = L1 * np.sin(theta) - L2 * np.sin(alpha)
    dg2dalpha = -L2 * np.cos(alpha);
    dg2dL3 = 0

    J = np.array([[dg1dalpha, dg1dL3], [dg2dalpha, dg2dL3]])

    s = np.array([[alpha], [L3]]) - J/np.array([[g1], [g2]])

    epsilon_alpha = abs(s[0] - alpha)

    epsilon_L3 = abs(s[1] - L3)

    epsilon = max(epsilon_alpha.all, epsilon_L3.all)

    alpha = s[0]

    L3 = s[1]

    n = n + 1

print(n, alpha, L3)
13
  • Is there any place where you assign an np.array(..) to a variable named np? Commented Sep 7, 2015 at 18:24
  • 3
    The error occurs because after a few iterations, alpha ends up being an array with dtype object. Commented Sep 7, 2015 at 18:32
  • @WarrenWeckesser: ah, good catch! (The code won't even get through one loop for me because I'm using Python 3 and so the max fails-- in 2 it'll just give an answer that I doubt the OP intends.) Commented Sep 7, 2015 at 18:34
  • @JustinWhite: Print alpha at the start of the loop. You'll see that it is probably not what you expect it to be. Commented Sep 7, 2015 at 18:34
  • I ran it and it failed with 'TypeError: unorderable types: builtin_function_or_method() > builtin_function_or_method()' on line 31 which is "epsilon = max(epsilon_alpha.all, epsilon_L3.all)" Commented Sep 7, 2015 at 18:52

1 Answer 1

3

In Python2.7 adding a print('alpha',alpha) at the start of the loop produces:

('alpha', 0.5235987755982988)
('alpha', array([ 1.85083849,  2.29325173]))
('alpha', array([[array([ 1.98227296,  1.95343536]), 1.7138098231972174],
       [array([ 1.81303794,  1.7604074 ]), 2.2932517265367176]], dtype=object))
Traceback (most recent call last):
  File "stack32444132.py", line 17, in <module>
    g1 = L1 * np.cos(theta) + L2 * np.cos(alpha) - L3
AttributeError: 'numpy.ndarray' object has no attribute 'cos'

So the error is caused by calling np.cos(alpha) where alpha is an object array. alpha is (2,2); the 1st column contains length 2 arrays; the 2nd contains floats.

So at some point in the loop you are appending or conconatenating arrays or lists of different length.

s = np.array([[alpha], [L3]]) - J/np.array([[g1], [g2]])
alpha = s[0]

adding some more prints (before s)

('J', (2, 2), dtype('float64'))
('alpha', 0.5235987755982988)
('L3', 1.5)
....
('J', (2, 2), dtype('O'))   
('alpha', array([ 1.85083849,  2.29325173]))
('L3', array([-10.61649234,   1.5       ]))

In the 2nd loop J changes from a 2x2 matrix of floats to a 2x2 of objects.

Python3 raises an error earlier, at the first encounter of the epsilon = max(epsilon_alpha.all, epsilon_L3.all) expression. epsilon_alpha.all is a method; epsilon_alpha.all() is a boolean. But even that produces an error when epsilon_alpha becomes an array.

OK, this loop runs (alpha remains a scalar); it doesn't stop, presumably because epsilon does not get small enough; but I'll leave that to you.

while epsilon > 0.0001:
    # print('alpha', alpha)
    g1 = L1 * np.cos(theta) + L2 * np.cos(alpha) - L3
    dg1dalpha = -L2 * np.sin(alpha)
    dg1dL3 = -1;

    g2 = L1 * np.sin(theta) - L2 * np.sin(alpha)
    dg2dalpha = -L2 * np.cos(alpha);
    dg2dL3 = 0

    J = np.array([[dg1dalpha, dg1dL3], [dg2dalpha, dg2dL3]])
    print('J', J.shape,J.dtype)  # (2,2) floats
    s = np.array([[alpha], [L3]]) - J/np.array([[g1], [g2]])
    s = s[:,0]  # fudge to turn (2,2) array into a (2,) array

    epsilon_alpha = abs(s[0] - alpha)
    epsilon_L3 = abs(s[1] - L3) 
    epsilon = max(epsilon_alpha, epsilon_L3)
    # max on 2 scalars is ok

    alpha = s[0] # scalar
    L3 = s[1]    # scalar
    n = n + 1

The root of the problem is in the

 s = np.array([[alpha], [L3]]) - J/np.array([[g1], [g2]])

If alpha and L3 are scalars, then np.array([[alpha], [L3]]) is (2,1). So is np.array([[g1], [g2]]). But since J is (2,2), s is also (2,2). But you've been using s[0] and s[1], evidently assuming s is `(2,).

 s = s[:,0] 

makes s a (2,) so the rest of the code works. Since epsilon does not converge it probably is the wrong fix.

May I stress - when developing numpy code, pay close attention to the array shapes. If the shapes are wrong you'll get errors like this. In my experience getting the shapes right is 80% of the debugging work.

Sign up to request clarification or add additional context in comments.

1 Comment

Wow, @hpaulj thanks for the detailed post. I need to learn how to do efficient debugging like you just did. I'll keep working on getting the shapes right as you suggested.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.