Array Morphisms

A unified backend for numerical computing combining fusion-based lazy evaluation, Mathematics of Arrays (MoA) index transformations, and automatic memory reuse.

Just as the Scheme programming language demonstrates that a small set of orthogonal primitives suffices for all of computation, array morphisms demonstrate that index transformations suffice for most of what we call "matrix operations."

The proliferation of reshape, transpose, slice, and broadcast in numerical multi-dimensional array libraries such as NumPy is an indication that the underlying abstraction, the coordinate-free linear map, is not a first-class citizen. A matrix is a linear mapping plus an arbitrary choice of basis; by making the mapping itself first-class and the coordinate representation derived, operations that appear distinct collapse into a single compositional abstraction over index spaces.

The Mathematics of Arrays (MoA), developed by Lenore Mullin, provides the formal foundation for this approach. MoA treats arrays not as containers of data but as functions from index spaces to values, and array operations as transformations of those index spaces. Libraries such as NumPy already exploit this insight as an internal optimization -- transpose, slice, and reshape can return views backed by the same storage rather than copies. What MoA makes possible, and what array morphisms realize, is promoting these optimizations from private implementation details to first-class abstractions. Index transformation become objects that can be inspected, composed, and extended. Users are no longer dependent on library authors having anticipated a particular combination of operations. Any composition of affine index functions is itself an affine index function, and the index algebra guarantees this closure. MoA's reduction theory (the Psi Calculus) then provides the rules for simplifying chains of such transformations before any data is ever accessed.

Array morphisms are the realization of MoA's index-transformation algebra as first-class Scheme objects. An array morphism pairs a shape, which serves as a description of the index space, with an index function that maps coordinates in that space to coordinates in some underlying storage. Morphisms compose: applying a transpose to a slice to a reshape produces a single new morphism whose index function is the composition of three affine maps, with no intermediate allocation. The representation is stratified: at the bottom lies a flat, contiguous buffer; above it, any number of morphisms may be stacked, each a pure description of how to interpret that buffer. Morphism evaluation, the act of reading or writing data, occurs only when explicitly requested, and only at the coordinates actually needed.

Key Features

Index Functions as First-Class Objects - Inspectable, composable, and3 extensible algebraic data
Unified Abstraction - No split between "generalized" and "specialized" arrays; all operations compose uniformly
Zero-Copy Views - Structural operations (reshape, transpose, slice) via affine index functions
Memory Reuse - Two-phase trace/replay execution with buffer planning via graph coloring
BLAS Integration - Transparent dispatch to optimized linear algebra kernels
Type Safety - Multiple element types (f64, f32, s64, s32, u32, u64) with automatic promotion
Category-Theoretic Foundation - Morphisms as structure-preserving transformations between arrays

Array Morphisms vs. SRFI-231

Index Functions as First-Class Objects: In SRFI-231, the indexer of a specialized array is exposed as an opaque procedure. Array morphisms make index functions inspectable algebraic data -- tagged records with known structure (affine, window, compute, composed). This makes the composition law programmable.
Unified Abstraction Rather Than Two-Tier Hierarchy: SRFI-231 maintains a hard split between generalized arrays (domain plus opaque getter/setter) and specialized arrays (domain plus flat buffer plus affine indexer). Array morphisms replace this with a single stratified structure where structural and computational morphisms compose uniformly.
Algebraic Simplification: SRFI-231 has no mechanism to simplify chains of operations. Array morphisms include a simplification layer: `reshape(reshape(x, s1), s2)` rewrites to `reshape(x, s2)`; `transpose(transpose(x, p1), p2)` rewrites to the composed permutation, and if the result is identity, the operations cancel entirely.
Zero-Copy Structural Operations: SRFI-231's `specialized-array-reshape` requires contiguous arrays and fails otherwise. Array morphisms treat reshape as a pure index transformation always valid, with contiguity being a property of execution strategy.
Deferred Evaluation with Buffer Reuse: SRFI-231 separates specification from execution but without any systematic approach. Array morphisms provide a two-phase execution model with trace/replay and automatic buffer planning via liveness analysis.
Batch Awareness: SRFI-231 has no notion of batch dimension. Array morphisms track the batch axis as a first-class attribute, enabling batch-aware broadcasting, convolution, and reduction.

Type System

Array morphisms support multiple element types with automatic promotion:

'f64 - Double precision floating point (64-bit)
'f32 - Single precision floating point (32-bit)
's64 - 64-bit signed integer
's32 - 32-bit signed integer
'u64 - 64-bit unsigned integer
'u32 - 32-bit unsigned integer

Type promotion rules:

Mixed operations promote to higher precision type
Transcendental functions promote integers to floating point
Reductions preserve dtype (mean promotes to float)

Core Concepts

Morphisms vs Arrays

Array morphisms represent computation as structure-preserving transformations. There are two types:

Concrete Arrays - Materialized data with shape, dtype, and strides
Abstract Morphisms - Deferred computations represented as expression trees

;; Concrete array - data is stored
(define concrete (morph-from-list '(1.0 2.0 3.0) '(3) 'f64))

;; Abstract morphism - represents computation
(define abstract (morph+ concrete concrete))

;; Realization materializes the morphism
(define result (realize abstract))  ; Now concrete

Index Functions

Index functions describe transformations algebraically:

Affine Index Functions - Pure transformations (reshape, transpose, slice)
Compute Index Functions - Element-wise arithmetic operations
Window Index Functions - Convolution and pooling operations
Reduction Index Functions - Aggregate operations (sum, mean, max)

Zero-Copy Views

Structural operations manipulate array views without copying:

(define x (morph-from-list '(0.0 1.0 2.0 3.0 4.0 5.0 6.0 7.0) '(8) 'f64))

;; Non-contiguous slice (stride 2)
(define strided (morph-slice x '(0) '(8) 2))  ; (0.0 2.0 4.0 6.0)

;; Reshape works even on non-contiguous views
(define as-2x2 (morph-reshape strided #(2 2)))  ; ((0.0 2.0) (4.0 6.0))

Installation

chicken-install array-morphisms

Module Imports

(import array-morphisms-core)       ; Core data types and utilities
(import array-morphisms-basic-ops)  ; Arithmetic operations
(import array-morphisms-structural-ops)  ; Reshape, transpose, slice
(import array-morphisms-realization)     ; Materialization
(import array-morphisms-context)         ; Memory reuse contexts
(import array-morphisms-batch-ops)       ; Batch operations

Core API

Morphism Creation

[procedure] (morph-from-list lst shape dtype)

Create a concrete morphism from a nested list.

[procedure] (make-morphism data shape dtype #!key batch-axis)

Create a concrete morphism from a typed vector.

Morphism Properties

[procedure] (morph-shape m)

Get the shape of morphism m.

[procedure] (morph-dtype m)

Get the data type of morphism m.

[procedure] (morph-size m)

Get the total number of elements in morphism m.

[procedure] (morph-rank m)

Get the number of dimensions in morphism m.

[procedure] (concrete-array? m)

Check if morphism is concrete (materialized).

[procedure] (abstract-morphism? m)

Check if morphism is abstract (deferred).

[procedure] (batched? m)

Check if morphism carries a batch dimension.

[procedure] (batch-size m)

Get the batch size (error if not batched).

Realization

[procedure] (realize m)

Materialize morphism to concrete array. Uses context if active.

[procedure] (realize! m)

Alias for realize.

[procedure] (morph->list m)

Convert concrete morphism to nested Scheme list.

Arithmetic Operations

All operations are lazy by default and create morphism expressions.

[procedure] (morph+ m1 m2)
[procedure] (morph- m1 m2)
[procedure] (morph* m1 m2)
[procedure] (morph/ m1 m2)
[procedure] (morph-pow m1 m2)

Element-wise binary operations with broadcasting.

[procedure] (morph-negate m)
[procedure] (morph-abs m)
[procedure] (morph-sqrt m)
[procedure] (morph-exp m)
[procedure] (morph-log m)
[procedure] (morph-sin m)
[procedure] (morph-cos m)
[procedure] (morph-tan m)
[procedure] (morph-floor m)
[procedure] (morph-ceiling m)

Element-wise unary operations.

Comparison Operations

[procedure] (morph> m1 m2)
[procedure] (morph< m1 m2)
[procedure] (morph= m1 m2)
[procedure] (morph>= m1 m2)
[procedure] (morph<= m1 m2)

Element-wise comparison (returns 1.0 for true, 0.0 for false).

Functional Operations

[procedure] (morph-map func m)

Apply function element-wise to morphism.

[procedure] (morph-reduce op m #!optional axes keepdims?)

Reduce morphism along specified axes. Operations: 'sum, 'mean, 'max, 'min, 'prod.

Structural Operations

[procedure] (morph-reshape m new-shape)

Reshape morphism to new shape (zero-copy). Supports -1 for automatic inference.

[procedure] (morph-transpose m #!optional permutation)

Transpose morphism by permuting dimensions. Default reverses all axes.

[procedure] (morph-slice m start end #!optional step)

Extract subarray via affine transformation. Supports negative indices and step.

[procedure] (morph-stack morphisms axis)

Stack morphisms along new dimension.

[procedure] (morph-concat morphisms axis)

Concatenate morphisms along existing dimension.

[procedure] (morph-squeeze m #!optional axes)

Remove dimensions of size 1.

[procedure] (morph-unsqueeze m axis)

Add dimension of size 1 at specified axis.

[procedure] (morph-pad m padding #!optional mode value)

Pad morphism with constant, edge, or reflect mode.

Matrix Operations (BLAS)

[procedure] (morph-matmul A B)

Matrix-matrix multiply morphism: (M,K) x (K,N) -> (M,N).

[procedure] (morph-matvec A v)

Matrix-vector multiply morphism: (M,N) x (N,) -> (M,).

[procedure] (morph-dot v1 v2)

Vector dot product morphism: (N,) . (N,) -> scalar.

[procedure] (morph-axpy alpha x y)

AXPY morphism: result = alpha * x + y.

Convolution Operations

[procedure] (im2col-morph input kernel-size #!optional stride padding)

Image-to-column morphism for convolution.

[procedure] (col2im-morph col output-shape kernel-size #!optional stride)

Column-to-image morphism (adjoint of im2col).

Batch Operations

[procedure] (add-batch-dimension m #!optional axis)

Add a size-1 batch dimension to non-batched morphism.

[procedure] (remove-batch-dimension m)

Remove size-1 batch dimension.

[procedure] (extract-batch-element m i #!optional axis)

Extract a single element from the batch dimension (zero-copy).

[procedure] (stack-into-batch morphisms axis)

Stack list of morphisms into a batch.

[procedure] (batch-map func m #!optional axis)

Apply function to each element along batch dimension.

[procedure] (batch-reduce op m #!optional keepdims?)

Reduce across batch dimension.

[procedure] (batch-zip func m1 m2)

Combine two batched morphisms element-wise.

[procedure] (morph-batch-matmul A B)

Batched matrix multiply with automatic shape handling.

Memory Reuse Context

[procedure] (make-morphism-context)

Create execution context for buffer reuse.

[procedure] (realize/ctx ctx m)

Realize morphism using context for buffer planning.

[procedure] (finalize-context! ctx)

Complete trace phase and compute buffer assignment.

[procedure] (reset-context! ctx)

Reset context for replay phase.

[procedure] (context-stats ctx)

Get statistics about allocations in context.

[procedure] (print-context-plan ctx)

Display buffer assignment plan.

BLAS Backend Configuration

[procedure] (register-blas-backend! backend)

[procedure] (blas-enabled?)

Check if BLAS optimization is enabled.

[procedure] (enable-blas!)

Enable BLAS optimization.

[procedure] (disable-blas!)

Disable BLAS optimization (use pure Scheme).

[procedure] (blas-available?)

Check if a BLAS backend is registered.

Examples

Basic Operations

(import array-morphisms-core
        array-morphisms-basic-ops
        array-morphisms-structural-ops
        array-morphisms-realization)

;; Create arrays
(define a (morph-from-list '(1.0 2.0 3.0 4.0 5.0) '(5) 'f64))
(define b (morph-from-list '(2.0 4.0 6.0 8.0 10.0) '(5) 'f64))

;; Build lazy computation
(define squared (morph-map (lambda (x) (* x x)) a))
(define complex (morph+ squared (morph-map (lambda (x) (* x 2)) b)))

;; Force computation
(morph->list (realize complex))

Layer Normalization

(define (layer-norm x eps)
  (let* ((mean (morph-reduce 'mean x '(0)))
         (centered (morph- x mean))
         (variance (morph-reduce 'mean 
                               (morph* centered centered) 
                               '(0)))
         (std (morph-sqrt (morph+ variance eps))))
    (morph/ centered std)))

;; Create 4x6 activation matrix
(define activations
  (morph-from-list (map (lambda (i) (+ i 1.0)) (iota 24))
                   #(4 6) 'f64))

;; Apply layer normalization
(define normalized (layer-norm activations 1e-5))
(morph->list (realize normalized))

Memory Reuse with Context

(import array-morphisms-context)

;; Create context for repeated computation
(define ctx (make-morphism-context))

;; Trace phase
(define result1 (realize/ctx ctx expensive-morphism))
(finalize-context! ctx)

;; Replay phase - reuses buffers
(reset-context! ctx)
(define result2 (realize/ctx ctx expensive-morphism))

Batched Operations

(import array-morphisms-batch-ops)

;; Create batch of matrices
(define m1 (morph-from-list '(1.0 2.0 3.0 4.0) '(2 2) 'f64))
(define m2 (morph-from-list '(2.0 4.0 6.0 8.0) '(2 2) 'f64))
(define batch (morph-stack (list m1 m2) 0))  ; shape (2 2 2)

;; Batch map: apply function to each matrix
(define doubled (batch-map 
                  (lambda (m) (morph-map (lambda (x) (* x 2)) m))
                  batch))

;; Batch reduce: sum across batch dimension
(define summed (batch-reduce 'sum doubled))

Zero-Copy Structural Chain

(define x (morph-from-list '(0.0 1.0 2.0 3.0 4.0 5.0 6.0 7.0 
                             8.0 9.0 10.0 11.0 12.0 13.0 14.0 15.0)
                           '(16) 'f64))

;; Stride-2 downsample
(define x-even (morph-slice x '(0) '(16) 2))      ; shape (8)

;; Stride-2 again - composed becomes stride-4
(define x-quarter (morph-slice x-even '(0) '(8) 2)) ; shape (4)

;; Reshape and transpose - all zero-copy views
(define as-2x2 (morph-reshape x-quarter #(2 2)))
(define transposed (morph-transpose as-2x2 '(1 0)))

;; Only when we realize is data accessed
(morph->list (realize transposed))

BLAS Integration

(import array-morphisms-blas-exec)

;; Matrix multiply automatically uses BLAS when available
(define A (morph-from-list '(1.0 2.0 3.0 4.0 5.0 6.0) '(2 3) 'f64))
(define B (morph-from-list '(7.0 8.0 9.0 10.0 11.0 12.0) '(3 2) 'f64))

;; This will dispatch to BLAS GEMM if a backend is registered
(define C (morph-matmul A B))

;; Register a BLAS backend (e.g., CRUNCH-compiled kernels)
(cond-expand
  ((library chicken-crunch-blas)
   (import chicken-crunch-blas)
   (register-blas-backend! (make-crunch-blas-backend))))

(morph->list (realize C))

Performance Notes

Lazy Evaluation: Operations build expression trees, materialize once
Zero-Copy Views: Structural operations (reshape, transpose, slice) are free
Memory Reuse: Use contexts for repeated computations to enable buffer planning
Batch Operations: Process multiple arrays together for efficiency
BLAS Dispatch: Automatically uses optimized kernels for large operations

datatype
matchable
srfi-1
srfi-4
srfi-69

License

GPL-3

References

Mullin, Lenore M.R. A Mathematics of Arrays. Ph.D. Thesis, Syracuse University, 1988.

The MoA formalism provides the theoretical foundation for index transformations as compositional abstractions.