How to find/build implementation of protocol described in core.matrix.impl.defaults?

I’ve run across some interesting work in the default implementation of householder-qrj protocol (function?) found here:

clojure.core.matrix.impl.defaults/householder-qr (defaults.cljc:2501)

See also this interesting discussion of these implementations, from 2013!

Does anyone know how could I find/build or otherwise use an implementation to experiment with this?

use an implementation to experiment with this

Since this is a specific helper function used to implement the actual protocol PQRDecomposition/qr, and the helper itself is written in terms of core.matrix protocols apparently, then any core.matrix implementation should work. So even the default core.matrix implementation (I think vectorz-clj is included by default, as well as the protocols working with double arrays and persistent vectors) should work with this as-is.

Thanks @joinr ! I must just be doing it wrong then, using the :vectorz implementation somehow:

(require '[clojure.core.matrix :as m])
=> nil
clojure.core.matrix.implementations/*matrix-implementation*
=> :vectorz

Then attempting to define this test function:

(defn tstHouse [g]
    (let [n (first (m/shape g))
          us (double-array n)
          gamma (double-array n)]
         (def rv (m/householder-qr lap3 0 n n us gamma))
         (prn rv)
        )
    )

generates this compile error:

Syntax error compiling at (/private/var/folders/v0/f06jgp6n6rq140nhcb3_9xpw0000gp/T/form-init6998699632090313895.clj:4:13).
No such var: m/householder-qr

Here’s a test matrix I was going to use, FYI:

(def lap3 '[[4 -1 -1 -1 0 0 0 -1] [-1 4 -1 -1 0 0 0 -1] [-1 -1 3 -1 0 0 0 0] [-1 -1 -1 5 -1 0 0 -1] [0 0 0 -1 4 -1 -1 -1] [0 0 0 0 -1 3 -1 -1] [0 0 0 0 -1 -1 3 -1] [-1 -1 0 -1 -1 -1 -1 6]] )

I don’t know the specifics of this particular case but I had similar issues. Every implementations of ‘core.matrix’ seemed to have a different broken/unimplemented subset. For instance the default backend doesnt implement the LU decomposition, but it you switch to ‘vectorz’ then it works :man_shrugging:

Fortunately it’s really easily to swap backends, so I’d change backends depends on what I was doing

I’ve tried :vectorz, :clatrix, :clojure, :ndarray with (set-current-implementation), all fail with the same

No such var: m/householder-qr

Here is the poorly documented solution I figured out after inspecting the code. The point you are missing is that, while you “can” invoke the householder function, you really aren’t meant to. Instead, you leverage the qr protocol function, which in turn uses the householder helper function in its default implementation (the PQRDecomposition protocol is extended to the base Object type). So then it “just works” if you provide it a matrix implementation that has all the requisite protocols. The other implicit/undocumented trick is that you have to pass an options map on what to return. It looks like it takes a sequence of keys, expecting [:Q :R] e.g. if you want both back. There is also an option for the key :compact, which is boolean supplied to the helper compute-r, which I have no idea what it’s supposed to do.

(ns qrtest
  (:require [clojure.core.matrix :as m]
             ;;qr lives here
             [clojure.core.matrix.protocols :as mp]))

(def lap3 '[[4 -1 -1 -1 0 0 0 -1]
            [-1 4 -1 -1 0 0 0 -1]
            [-1 -1 3 -1 0 0 0 0]
            [-1 -1 -1 5 -1 0 0 -1]
            [0 0 0 -1 4 -1 -1 -1]
            [0 0 0 0 -1 3 -1 -1]
            [0 0 0 0 -1 -1 3 -1]
            [-1 -1 0 -1 -1 -1 -1 6]] )

(defn qr [g]
  (let [n (first (m/shape g))
        us (double-array n)
        gamma (double-array n)]
   (mp/qr g {:return [:Q :R]})))

#_
(qr lap3)
#_
{:Q
 [[-0.8944271909999157
   -2.629045857405082E-17
   1.9717843930538115E-17
   -0.054110331064667484
   0.08186963585930283
   0.10853774097467565
   0.23149296658959517
   -0.3535533905932737]
  [0.22360679774997896
   -0.8660254037844384
   -2.327375768669318E-18
   -0.0541103310646675
   0.08186963585930282
   0.10853774097467567
   0.23149296658959512
   -0.3535533905932737]
  [0.22360679774997896
   0.2886751345948128
   -0.792593923901217
   -0.2414153232115937
   0.028267691860215334
   0.09899596154833055
   0.21114193655974067
   -0.3535533905932737]
  [0.22360679774997896
   0.2886751345948128
   0.5661385170722979
   -0.616025307505446
   -0.0789361961379597
   0.07991240269564033
   0.17043987650003162
   -0.3535533905932736]
  [0.0
   0.0
   0.0
   0.2705516553233378
   -0.8829653561242725
   -0.06321428869953637
   -0.13482557394778624
   -0.35355339059327373]
  [0.0
   0.0
   0.0
   0.0
   0.23680858841387922
   -0.8799663488334247
   -0.21114193655974065
   -0.3535533905932738]
  [0.0
   0.0
   0.0
   0.0
   0.23680858841387922
   0.4004919326595827
   -0.811497322440449
   -0.35355339059327373]
  [0.22360679774997896
   0.2886751345948128
   0.22645540682891913
   0.6951096375230371
   0.29627741185565293
   0.14670485868005614
   0.31289708670901334
   -0.35355339059327384]],
 :R
 [[-4.47213595499958
   1.1180339887498945
   1.1180339887498945
   1.3416407864998736
   -0.4472135954999579
   -0.22360679774997896
   -0.22360679774997896
   1.788854381999831]
  [0
   -4.330127018922194
   1.443375672974064
   1.7320508075688765
   -0.5773502691896257
   -0.2886751345948128
   -0.2886751345948128
   2.3094010767585025]
  [0
   0
   -2.943920288775949
   3.3968311024337874
   -0.792593923901217
   -0.2264554068289191
   -0.2264554068289191
   0.7925939239012167]
  [0
   0
   0
   -3.6961518450326776
   1.00312229127576
   -0.9656612928463746
   -0.9656612928463746
   4.624352139449664]
  [0
   0
   0
   0
   -4.222819817042542
   1.0603051210963779
   1.0603051210963779
   2.102209574849785]
  [0 0 0 0 0 -3.123881549140376 1.9979515768316531 1.1259299723087235]
  [0 0 0 0 0 0 -2.4014215435228343 2.4014215435228334]
  [0 0 0 0 0 0 0 -4.440892098500626E-16]]}
2 Likes

Hi @joinr thanks very much for these notes, and I apologize for my delayed response; RealLife gets in the way. I’ll ponder…

Just realized in re-reading this, the stuff in the let binding is not necessary:

(defn qr [g] (mp/qr g {:return [:Q :R]})))

Since you already have a valid matrix, just use the qr protocol function, e.g.

(mp/qr g {:return [:Q :R]})

If you want to submit your own args for the lower level householder function, then the stuff you supplied would come into play (although apparently it’s not in the interface).

2 Likes