Clojure Numerics: Orthogonalization and Least Squares | Hacker News Reader

Clojure Numerics: Orthogonalization and Least Squares

dragan.rocks

118 pointsdragandj20 comments

Source code of the library is at: https://github.com/uncomplicate/neanderthal

What’s the benefit of writing this in Clojure, instead of in Java and designing it such that calling it from Clojure feels natural? Does it use any Clojure specific features or have strict Clojure dependencies?

In Java, it would require at least an order of magnitude more code. Additionally, you'd have to know exactly how to do it from the beginning, since it is much more difficult to experiment. High performance computing is not something that is fabulously well documented, so you have to poke around a lot to discover how things work. To see what I'm talking about, you should compare the codebase of neanderhtal (20.000 lines of code) with roughly similar libraries in Java (ND4J is something like 400.000 lines of code), or even Python (numpy is roughly 250.000).

Bit offtopic. Would you recommend clojure as a web development programming language? I am from a python background and want to learn a new language. Previously thought of Rust and Haskell, rust seems to be too low level and carries mental burden of manual memory management and haskell seems too complicated for me. Thanks in advance

I would. Web dev is (for better or worse :) one of the main applications of Clojure. See further at https://clojurescript.org/ and https://pragprog.com/book/dswdcloj2/web-development-with-clo...

Clojurescript today is probably the best, most pragmatic choice to do web dev in. It's not difficult to learn, FFI isn't cumbersome, compiled code is extremely optimized, tooling is pretty nice, syntax is a joy - writing Clojure is like poetry. I would never ever touch vanilla JS or TS (unless I really have to).

Does neanderthal claim feature parity with numpy and ND4J?

Do numpy and ND4J have feature parity among themselves?

Fine...

I should have said "or" instead of "and". I think it's a valid question in the context of your code size comparison.

The point is that Neanderthal is not a port of either of those libraries. Sure, it does not have all the features of numpy or ND4J, but, on the other hand, it has many features that neither numpy nor ND4J have.

I'm not very familiar with Clojure, but are the numerics done with unboxed numbers? In the end, what every language that wants to do scientific programming right is expose some sort of c-array like data structure that you can pass onto c/fortran code.

Do not worry. Neanderthal works with barest primitives, and uses Intel MKL, cuBLAS, CLBlast, and custom kernels at the low level. This is practicaly as fast as you can get in a general library.

> Neanderthal works with barest primitives

I see what you did there with naming.

How do you use MKL and friends on objects in the Java heap? It seems like you'd have to copy to the native heap, do your linear algebra operation, and then copy the result back from the native heap.

It is not used on objects on the java heap. No copying occurs.

Aha, makes sense. (Now I see that at http://neanderthal.uncomplicate.org/articles/tutorial_native...) Doesn't this make it very painful to do unusual operations on matrices? For example, solving linear systems Ax=b where A is of the form

  [ 1                                             ]
  [ a_2 b_1  1                                    ]
  [ a_3 b_1  a_3 b_2  1                           ]
  [ ...      ...      ...      ...                ]
  [ a_n b_1  a_n b_2  a_n b_3  ...  a_n b_{n-1} 1 ].

I do not understand what the form of that matrix has to do with java or native heap. This seems to me as a completely orthogonal issue. As for the triangular form, this is supported in Neanderthal, as well as the rest of special structural sparse shapes.

> I do not understand what the form of that matrix has to do with java or native heap. This seems to me as a completely orthogonal issue. As for the triangular form, this is supported in Neanderthal, as well as the rest of special structural sparse shapes.

You support triangular matrices. However, you can solve a linear system of the form I gave in linear time, while it takes quadratic space and time to form the corresponding triangular matrix and do a triangular solve.

Not all special matrices are supported by BLAS/LAPACK. Other common examples might be block Toeplitz/Hankel matrices for which fast multiplication and fast solvers are available. In order to support special (not in BLAS/LAPACK) matrix operations naturally, you'd want natural, no-extra-copying access to the vectors within Java or Clojure so that you can write the good algorithm manually, as you'd do in C or Fortran.

Sure! If you need to efficiently implement that yourself, you can use OpenCL to implement the kernels on the CPU, or OpenCL & CUDA on the GPU (If that makes sense). Check out ClojureCUDA and ClojureCL.

Thanks for writing this series. I'm really interested in this stuff and I'm sorta writing my own version where I'm implementing the algorithms in lisp: http://geokon-gh.github.io/linearsystems/ (I'm sure it's quite amateurish compared to what you're doing)

Nothing quite forces you to understand the issues deeply like implementing it in code!

I haven't gotten to LS/QR yet (just getting the LU working correctly was quite tricky), but when I do I'll make sure to reference your series :)

Reply on news.ycombinator.com