Why did not numpy copy the J rank concept?

14 points by jrank a day ago

Recently there was a post [1] about numpy being not easy to work with when using arrays of shape greater than 2.

  One of the problems mentioned in [1] is that  you can not use python loops because they are slow.  In J you can solve for example 100 equations using the rank concept, this is a simple  example:

   a=: 2 2 $ 1 1 1 _1
   b=: 10 2 $ ? 20 # 10
   solutions =: b %. "(1 _) a
That code solve the systems a * v_i = b_i for ten random vectors. I think a similar concept could be developed in numpy. The syntax "(1 _)" indicates to take the rows from the left operator and all (_ is infinite) of a apply solve (that is %. in J). In this case the system is x+y=y0, x-y=y1.

So I would suggest somthing like numpy.linalg.solve(a,b,rank=(1,inf))

[1] https://news.ycombinator.com/item?id=43996431

tester89 a day ago

The question you're asking seems interesting, but I don't understand J code so I don't know what you're talking about. Expanding the explanation would be helpful!

  • jrank a day ago

    If you have a a Numpy function that takes for example two arguments, my proposal is to add a optional argument that allows that function to be applied to cells of dimensions i and j by function_name(a,b,range=(i,j)) so that the function is applied to subarrays of dimension i of a and subarrays of dimension j of b to create a new array. The broadcast operation and the axis arguments are not a general solution. I J you have such mechanism as the example shows.

    The 1 cell of a matrix are the rows, etc.

    • gus_massa 13 hours ago

      I'm still confused. Can you type the complete slow python code?

      We had a weird multiplication last year and the solution in python was to use @guvectorize https://numba.pydata.org/numba-doc/dev/user/vectorize.html and broadcasting. I don't know if that's possible with your problem.

      • jrank 5 hours ago

        guvectorize in numba seems to be a good approximation to the rank concept I mentioned, but is not a complete solution. Unfortunately I don't have the time now to study it and make a complete comparison, but guvectorize is in the good direction. Thanks for providing that information.

Pompidou a day ago

Maybe the internal broadcasting mecanism in numpy don't allow this nativelly ?

  • jrank a day ago

    The post I referenced before shows that broadcast is not a general solution.