NumPy gives off NaN far too often compared with Octave. This is understandable from the CompSci origin. Python is heavily integer-oriented (1/2==0) while NumPy still favors real numbers much more than complex ones (sqrt(-2)==nan). These are very awkward for matrix computations because many moderately complicated matrix computations quickly leads to full of nans.
The short term solution is to add 0j to the operands to force complex computation.
Long term solution is have more robust implementations. The experience of Octave and its large code base (which dates back to some excellent Fortran packages) should be well utilized in this regard. It is not clear to me how, but it could be extremely useful, to directly implement many operations in terms of the Octave source code base.