2  Quickstart

A minimal introduction to La Linealinear algebra with tensor abstractions. Matrices are dtype-next tensors, backed by EJML for computation and shared via zero-copy interop.

(ns lalinea-book.quickstart
  (:require
   ;; La Linea (https://github.com/scicloj/lalinea):
   [scicloj.lalinea.linalg :as la]
   [scicloj.lalinea.elementwise :as el]
   [scicloj.lalinea.tensor :as t]
   ;; FFT bridge — Fastmath transforms <-> ComplexTensor:
   [scicloj.lalinea.transform :as ft]
   ;; Visualization annotations (https://scicloj.github.io/kindly-noted/):
   [scicloj.kindly.v4.kind :as kind]
   [clojure.math :as math]))

Creating matrices

Matrices are dtype-next tensors of shape \([r, c]\), backed by a flat double[] in row-major order.

(t/matrix [[1 2 3]
           [4 5 6]])
#la/R [:float64 [2 3]
       [[1.000 2.000 3.000]
        [4.000 5.000 6.000]]]

Identity and zero matrices:

(t/eye 3)
#la/R [:float64 [3 3]
       [[1.000 0.000 0.000]
        [0.000 1.000 0.000]
        [0.000 0.000 1.000]]]
(t/zeros 2 3)
#la/R [:float64 [2 3]
       [[0.000 0.000 0.000]
        [0.000 0.000 0.000]]]

Diagonal matrices:

(t/diag [1 2 3])
#la/R [:float64 [3 3]
       [[1.000 0.000 0.000]
        [0.000 2.000 0.000]
        [0.000 0.000 3.000]]]

Matrix arithmetic

All operations accept tensors and return tensors.

(la/mmul (t/matrix [[1 2] [3 4]])
         (t/matrix [[5 6] [7 8]]))
#la/R [:float64 [2 2]
       [[19.00 22.00]
        [43.00 50.00]]]

\(AB = \begin{pmatrix} 19 & 22 \\ 43 & 50 \end{pmatrix}\)

(la/transpose (t/matrix [[1 2] [3 4]]))
#la/R [:float64 [2 2]
       [[1.000 3.000]
        [2.000 4.000]]]

Scalar properties

(la/det (t/matrix [[1 2] [3 4]]))
-2.0
(la/trace (t/matrix [[1 2] [3 4]]))
5.0
(la/norm (t/matrix [[1 2] [3 4]]))
5.477225575051661

Solving linear systems

Solve \(Ax = b\) where \(A = \begin{pmatrix} 2 & 1 \\ 1 & 3 \end{pmatrix}\), \(b = \begin{pmatrix} 5 \\ 10 \end{pmatrix}\).

(la/solve (t/matrix [[2 1] [1 3]])
          (t/matrix [[5] [10]]))
#la/R [:float64 [2 1]
       [[1.000]
        [3.000]]]

\(x = \begin{pmatrix} 1 \\ 3 \end{pmatrix}\)

Decompositions

Eigendecomposition of a symmetric matrix:

(la/real-eigenvalues (t/matrix [[2 1] [1 2]]))
#la/R [:float64 [2] [1.000 3.000]]

Eigenvalues are \(3\) and \(1\).

SVD:

(:S (la/svd (t/matrix [[1 2] [3 4]])))
#la/R [:float64 [2] [5.465 0.3660]]

Complex tensors

ComplexTensors wrap a dtype-next tensor whose last dimension is 2 (interleaved re/im pairs).

The same la/ functions work for both real and complex inputs — el/+, la/mmul, la/transpose are polymorphic over the number field.

(t/complex-tensor [1.0 2.0 3.0] [4.0 5.0 6.0])
#la/C [:float64 [3] [1.000 + 4.000 i  2.000 + 5.000 i  3.000 + 6.000 i]]

Complex matrix multiply through the same mmul:

(let [A (t/complex-tensor [[1.0 0.0] [0.0 1.0]]
                          [[0.0 0.0] [0.0 0.0]])]
  (la/mmul A A))
#la/C [:float64 [2 2]
       [[1.000 + 0.000 i  0.000 + 0.000 i]
        [0.000 + 0.000 i  1.000 + 0.000 i]]]

FFT bridge

The forward FFT takes a real signal and returns a ComplexTensor spectrum — zero-copy from Fastmath’s interleaved output.

(ft/forward [1.0 0.0 1.0 0.0])
#la/C [:float64 [4] [2.000 + 0.000 i  0.000 + 0.000 i  2.000 + 0.000 i  0.000 + 0.000 i]]

\(\hat{f} = [2, 0, 2, 0]\) — a signal with energy at DC and Nyquist.

Round-trip:

(let [signal [1.0 2.0 3.0 4.0]
      recovered (ft/inverse-real (ft/forward signal))]
  recovered)
#la/R [:float64 [4] [1.000 2.000 3.000 4.000]]

Composing with dtype-next

Since matrices are tensors, all dtype-next operations work.

(el/sum (t/matrix [[1 2] [3 4]]))
10.0

Element-wise operations preserve tensor structure:

((el/scale (t/matrix [[1 2] [3 4]]) 2.0) 1 1)
8.0

Convenience wrappers

Common SVD-based analyses are available as one-liners:

(la/rank (t/matrix [[1 2] [2 4]]))
1
(la/condition-number (t/matrix [[2 1] [1 3]]))
2.618033988749894

Pseudoinverse:

(la/close? (la/mmul (t/matrix [[2 1] [1 3]])
                    (la/pinv (t/matrix [[2 1] [1 3]])))
           (t/eye 2))
true

Matrix power:

(la/mpow (t/matrix [[1 1] [0 1]]) 5)
#la/R [:float64 [2 2]
       [[1.000 5.000]
        [0.000 1.000]]]

Tagged literals

La Linea tensors print as #la/R and #la/C tagged literals, enabling round-trip through pr-str / read-string:

(pr-str (t/matrix [[1 2] [3 4]]))
"#la/R [:float64 [2 2]\n       [[1.000 2.000]\n        [3.000 4.000]]]\n"
(pr-str (t/column [5 6 7]))
"#la/R [:float64 [3 1]\n       [[5.000]\n        [6.000]\n        [7.000]]]\n"

Element-wise functions

la/ covers linear algebra; el/ covers element-wise math — both are tape-aware and work with both real and complex inputs.

scicloj.lalinea.elementwise provides tape-aware element-wise math functions (sqrt, sin, exp, …):

(el/exp (t/column [0.0 1.0 2.0]))
#la/R [:float64 [3 1]
       [[1.000]
        [2.718]
        [7.389]]]
(el/clip (t/column [-2 0.5 3]) -1 1)
#la/R [:float64 [3 1]
       [[-1.000]
        [0.5000]
        [ 1.000]]]

Computation tape

Record operations as a DAG:

(require '[scicloj.lalinea.tape :as tape])
(let [{:keys [entries]} (tape/with-tape
                          (la/mmul (t/matrix [[1 2] [3 4]])
                                   (t/column [1 0])))]
  (mapv :op entries))
[:t/matrix :t/column :la/mmul]

Automatic differentiation

Reverse-mode autodiff computes gradients via VJP rules:

(require '[scicloj.lalinea.grad :as grad])
(let [A (t/matrix [[1 2] [3 4]])
      tape-result (tape/with-tape
                    (el/sum (el/sq (el/- (la/mmul A A)
                                           (t/matrix [[1 0] [0 1]])))))]
  ((grad/grad tape-result (:result tape-result) A) 0 0))
154.0
source: notebooks/lalinea_book/quickstart.clj