On Di, 2014-07-15 at 10:22 +0100, Neil Hodgson wrote:
Hi,
We came across this bug while using np.cross on 3D arrays of 2D vectors.
Hi, which numpy version are you using? Until recently, the cross product simply did *not* work in a broadcasting manner (3d arrays of 2d vectors), it did something, but usually not the right thing. This is fixed in recent versions (not sure if 1.8 or only now with 1.9) - Sebastian
The first example shows the problem and we looked at the source for np.cross and believe we found the bug - an unnecessary swapaxes when returning the output (comment inserted in the code).
Thanks Neil
# Example
shape = (3,5,7,2)
# These are effectively 3D arrays (3*5*7) of 2D vectors data1 = np.random.randn(*shape) data2 = np.random.randn(*shape)
# The cross product of data1 and data2 should produce a (3*5*7) array of scalars cross_product_longhand = data1[:,:,:,0]*data2[:,:,:,1]-data1[:,:,:,1]*data2[:,:,:,0] print 'longhand output shape:',cross_product_longhand.shape # and it does
cross_product_numpy = np.cross(data1,data2) print 'numpy output shape:',cross_product_numpy.shape # It seems to have transposed the last 2 dimensions
if (cross_product_longhand == np.transpose(cross_product_numpy, (0,2,1))).all(): print 'Unexpected transposition in numpy.cross (numpy version %s)'% np.__version__
# np.cross L1464 if axis is not None: axisa, axisb, axisc=(axis,)*3 a = asarray(a).swapaxes(axisa, 0) b = asarray(b).swapaxes(axisb, 0) msg = "incompatible dimensions for cross product\n"\ "(dimension must be 2 or 3)" if (a.shape[0] not in [2, 3]) or (b.shape[0] not in [2, 3]): raise ValueError(msg) if a.shape[0] == 2: if (b.shape[0] == 2): cp = a[0]*b[1] - a[1]*b[0] if cp.ndim == 0: return cp else: ## WE SHOULD NOT SWAPAXES HERE! ## For 2D vectors the first axis has been
## collapsed during the cross product return cp.swapaxes(0, axisc)
_______________________________________________ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion