Mojo Matrix Operations

Matrices are two-dimensional grids of numbers. They appear everywhere in AI — weight tables in neural networks, image pixel grids, transformation operations in graphics, and covariance tables in statistics. Mojo lets you implement matrix operations with full control over memory layout and hardware utilization, reaching speeds that match specialized libraries.

Matrix Representation

  A 3×4 matrix has 3 rows and 4 columns:

       col0  col1  col2  col3
  row0 [  1    2    3    4  ]
  row1 [  5    6    7    8  ]
  row2 [  9   10   11   12  ]

  Row-major storage in memory (C order):
  [1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12]
   ←── row0 ──→←── row1 ──→←──── row2 ────→

  Element [row, col] is at index: row × num_cols + col
  Element [1, 2] = index 1×4 + 2 = 6 → value 7 ✓

Matrix Struct

from memory import UnsafePointer

struct Matrix:
    var data: UnsafePointer[Float32]
    var rows: Int
    var cols: Int

    fn __init__(inout self, rows: Int, cols: Int):
        self.rows = rows
        self.cols = cols
        self.data = UnsafePointer[Float32].alloc(rows * cols)
        for i in range(rows * cols):
            self.data.init_pointee_copy(0.0)

    fn __del__(owned self):
        for i in range(self.rows * self.cols):
            (self.data + i).destroy_pointee()
        self.data.free()

    fn get(self, row: Int, col: Int) -> Float32:
        return self.data[row * self.cols + col]

    fn set(inout self, row: Int, col: Int, value: Float32):
        self.data[row * self.cols + col] = value

    fn print_matrix(self):
        for r in range(self.rows):
            for c in range(self.cols):
                print(self.get(r, c), end="\t")
            print("")

Matrix Addition

Element-wise addition: add corresponding elements from two matrices of the same shape.

fn mat_add(a: Matrix, b: Matrix, result: inout Matrix):
    for r in range(a.rows):
        for c in range(a.cols):
            result.set(r, c, a.get(r, c) + b.get(r, c))

fn main():
    var a = Matrix(2, 3)
    var b = Matrix(2, 3)
    var result = Matrix(2, 3)

    a.set(0, 0, 1.0); a.set(0, 1, 2.0); a.set(0, 2, 3.0)
    a.set(1, 0, 4.0); a.set(1, 1, 5.0); a.set(1, 2, 6.0)

    b.set(0, 0, 10.0); b.set(0, 1, 20.0); b.set(0, 2, 30.0)
    b.set(1, 0, 40.0); b.set(1, 1, 50.0); b.set(1, 2, 60.0)

    mat_add(a, b, result)
    result.print_matrix()

Output:

11.0    22.0    33.0
44.0    55.0    66.0

Scalar Multiplication

fn mat_scale(m: Matrix, factor: Float32, result: inout Matrix):
    for r in range(m.rows):
        for c in range(m.cols):
            result.set(r, c, m.get(r, c) * factor)

Matrix Transpose

Transposing swaps rows and columns. A 3×2 matrix becomes a 2×3 matrix. Element [r, c] moves to position [c, r].

fn transpose(m: Matrix, result: inout Matrix):
    # result must be (m.cols × m.rows)
    for r in range(m.rows):
        for c in range(m.cols):
            result.set(c, r, m.get(r, c))
Transpose diagram:

  Original (2×3):          Transposed (3×2):
  [ 1  2  3 ]              [ 1  4 ]
  [ 4  5  6 ]              [ 2  5 ]
                            [ 3  6 ]

  [r, c] → [c, r]

Matrix Multiplication

Matrix multiplication combines two matrices into a third. To multiply A (M×K) by B (K×N), each output element C[r, c] is the dot product of row r of A with column c of B.

fn matmul(a: Matrix, b: Matrix, c: inout Matrix):
    # a: M×K, b: K×N, c: M×N
    for r in range(a.rows):
        for col in range(b.cols):
            var sum: Float32 = 0.0
            for k in range(a.cols):
                sum += a.get(r, k) * b.get(k, col)
            c.set(r, col, sum)
Matrix multiply (2×3) × (3×2) → (2×2):

  A:                B:             C = A×B:
  [ 1  2  3 ]   [ 7  8  ]      [ 1×7+2×9+3×11  1×8+2×10+3×12 ]
  [ 4  5  6 ]   [ 9  10 ]   =  [ 4×7+5×9+6×11  4×8+5×10+6×12 ]
                [ 11 12 ]

  C[0,0] = 1×7 + 2×9 + 3×11 = 7 + 18 + 33 = 58
  C[0,1] = 1×8 + 2×10 + 3×12 = 8 + 20 + 36 = 64
  C[1,0] = 4×7 + 5×9 + 6×11 = 28 + 45 + 66 = 139
  C[1,1] = 4×8 + 5×10 + 6×12 = 32 + 50 + 72 = 154

  Result: [ 58   64 ]
          [139  154 ]

Vectorized Matrix Multiply

The inner dot-product loop is the hot path in matrix multiplication. Vectorize it with SIMD for a direct performance gain.

fn matmul_simd(a: Matrix, b: Matrix, c: inout Matrix):
    alias width = 8
    for r in range(a.rows):
        for col in range(b.cols):
            var acc = SIMD[DType.float32, width].splat(0.0)
            var k = 0
            while k + width <= a.cols:
                var av = SIMD[DType.float32, width].load(a.data + r * a.cols + k)
                var bv = SIMD[DType.float32, width].load(b.data + k * b.cols + col)
                acc += av * bv
                k += width
            var sum = acc.reduce_add()
            while k < a.cols:
                sum += a.get(r, k) * b.get(k, col)
                k += 1
            c.set(r, col, sum)

Identity Matrix

The identity matrix is the matrix equivalent of the number 1. Multiplying any matrix by the identity matrix returns the original matrix unchanged. The identity has 1s on the main diagonal and 0s everywhere else.

fn identity(n: Int, result: inout Matrix):
    for r in range(n):
        for c in range(n):
            result.set(r, c, 1.0 if r == c else 0.0)
3×3 Identity matrix:
  [ 1  0  0 ]
  [ 0  1  0 ]
  [ 0  0  1 ]

  A × I = A  (for any matrix A with matching dimensions)

Element-wise Operations

Unlike matrix multiplication, element-wise operations apply a function to each corresponding pair of elements independently. These appear everywhere in AI — applying activation functions, scaling layers, computing losses.

fn element_multiply(a: Matrix, b: Matrix, result: inout Matrix):
    for r in range(a.rows):
        for c in range(a.cols):
            result.set(r, c, a.get(r, c) * b.get(r, c))
Element-wise multiply (Hadamard product):
  [ 2  3 ]   ⊙   [ 5  6 ]   =   [ 10  18 ]
  [ 4  5 ]       [ 7  8 ]       [ 28  40 ]

  Notation: ⊙ (not the same as matrix multiply ×)

Performance Summary

Operation             | Complexity | Key optimization
----------------------|------------|-------------------------------
Matrix add            | O(M×N)     | SIMD row scan
Matrix transpose      | O(M×N)     | Cache-friendly access pattern
Matrix multiply       | O(M×K×N)   | SIMD inner loop + tiling
Element-wise ops      | O(M×N)     | SIMD + vectorize
Row/column reduction  | O(M×N)     | SIMD reduce_add per row

Key Takeaways

Matrices store two-dimensional data in a flat array using row-major order. Element [r, c] lives at index r × num_cols + c. Matrix addition and element-wise multiplication operate independently on each position. Transpose swaps row and column indices. Matrix multiplication computes a dot product between each row of the left matrix and each column of the right — it is the central operation in neural networks. Vectorize the inner loop with SIMD and tile the outer loops for cache efficiency to reach peak performance on large matrices.

Leave a Comment

Your email address will not be published. Required fields are marked *