Throughout this course, the configuration, velocities, and forces of a robot will be encoded as vectors and matrices. Even when the true relationships are nonlinear, we often work with local linear approximations. For example, for generalized coordinates \( \mathbf{q} \in \mathbb{R}^n \) and some task-space variable \( \mathbf{x} \in \mathbb{R}^m \), later lessons will rely heavily on linear relations of the form

\[ \mathbf{x} \approx \mathbf{x}_0 + \mathbf{J}(\mathbf{q}_0)\,\Delta\mathbf{q}, \qquad \dot{\mathbf{x}} = \mathbf{J}(\mathbf{q})\,\dot{\mathbf{q}}, \]

where \( \mathbf{J}(\mathbf{q}) \) is a matrix, and solving for unknown vectors (e.g. \( \Delta\mathbf{q} \), \( \dot{\mathbf{q}} \)) requires understanding rank and null spaces. In this lesson we focus purely on the linear algebraic objects; later chapters interpret them as geometric and physical quantities.

flowchart TD
  A["Robot geometry (links, joints)"] --> B["Generalized coordinates q"]
  B --> C["Linear algebra objects: vectors & matrices"]
  C --> D["Equations of form A x = b"]
  D --> E["Solve / analyze: rank, Null(A), Range(A)"]
  E --> F["Use solutions in kinematics & dynamics models"]
        

A (column) vector in \( \mathbb{R}^n \) is written as

\[ \mathbf{v} = \begin{bmatrix} v_1 \\ v_2 \\ \vdots \\ v_n \end{bmatrix} \in \mathbb{R}^n. \]

In robotics, we will use vectors for joint positions, velocities, wrenches, etc. The standard inner product (dot product) of vectors \( \mathbf{u}, \mathbf{v} \in \mathbb{R}^n \) is

\[ \mathbf{u}^\top \mathbf{v} = \sum_{i=1}^n u_i v_i. \]

This induces the Euclidean norm

\[ \| \mathbf{v} \|_2 = \sqrt{\mathbf{v}^\top \mathbf{v}} = \sqrt{v_1^2 + \cdots + v_n^2}. \]

Two vectors are orthogonal if \( \mathbf{u}^\top \mathbf{v} = 0 \). Orthogonality plays a key role when decomposing motion or forces into components along different directions.

The projection of \( \mathbf{v} \) onto a nonzero vector \( \mathbf{u} \) is

\[ \operatorname{proj}_{\mathbf{u}}(\mathbf{v}) = \frac{\mathbf{u}^\top \mathbf{v}}{\mathbf{u}^\top \mathbf{u}}\,\mathbf{u}. \]

Projections are central in least-squares problems, where we project a desired task velocity onto the subspace of velocities that can be generated by the robot.

A matrix \( A \in \mathbb{R}^{m \times n} \) defines a linear map \( A : \mathbb{R}^n \to \mathbb{R}^m \) via \( \mathbf{y} = A\mathbf{x} \). Linearity means:

\[ A(\alpha \mathbf{x}_1 + \beta \mathbf{x}_2) = \alpha A\mathbf{x}_1 + \beta A\mathbf{x}_2 \quad \text{for all scalars } \alpha,\beta \text{ and vectors } \mathbf{x}_1, \mathbf{x}_2. \]

Writing \( A = [\,\mathbf{a}_1\ \mathbf{a}_2\ \cdots\ \mathbf{a}_n\,] \) in terms of its columns and \( \mathbf{x} = [\,x_1,\dots,x_n\,]^\top \), we obtain

\[ A\mathbf{x} = x_1 \mathbf{a}_1 + x_2 \mathbf{a}_2 + \cdots + x_n \mathbf{a}_n, \]

so \( A\mathbf{x} \) is a linear combination of the columns. In robotics, matrices will encode coordinate transforms, velocity mappings, and dynamic couplings.

Composition of linear maps corresponds to matrix multiplication: if \( A \in \mathbb{R}^{p \times m} \) and \( B \in \mathbb{R}^{m \times n} \), then \( A(B\mathbf{x}) = (AB)\mathbf{x} \) for all \( \mathbf{x} \in \mathbb{R}^n \).

This is exactly how we will later chain multiple joint and link transforms to get a full kinematic chain: product of homogeneous transform matrices corresponds to the composition of the underlying rigid-body motions.

Vectors \( \mathbf{v}_1,\dots,\mathbf{v}_k \in \mathbb{R}^n \) are linearly independent if the only coefficients satisfying

\[ \alpha_1 \mathbf{v}_1 + \cdots + \alpha_k \mathbf{v}_k = \mathbf{0} \]

are \( \alpha_1 = \cdots = \alpha_k = 0 \). Their span is the set of all linear combinations \( \operatorname{span}\{\mathbf{v}_1,\dots,\mathbf{v}_k\} \subseteq \mathbb{R}^n \).

For a matrix \( A \in \mathbb{R}^{m \times n} \), the column space (or range) is

\[ \mathcal{R}(A) = \{ A\mathbf{x} : \mathbf{x} \in \mathbb{R}^n \} = \operatorname{span}\{\mathbf{a}_1,\dots,\mathbf{a}_n\}, \]

where \( \mathbf{a}_j \) are the columns of \( A \). The rank of \( A \), denoted \( \operatorname{rank}(A) \), is defined as

\[ \operatorname{rank}(A) = \dim \mathcal{R}(A), \]

i.e., the maximal number of linearly independent columns (equivalently, rows).

One can show that \( \operatorname{rank}(A) \leq \min\{m,n\} \). Sketch of proof:

• Every set of more than \( m \) vectors in \( \mathbb{R}^m \) is linearly dependent.
• The columns of \( A \) lie in \( \mathbb{R}^m \), so at most \( m \) can be independent.
• Similarly, rows lie in \( \mathbb{R}^n \), so at most \( n \) can be independent.

Therefore at most \( \min\{m,n\} \) independent columns (or rows) exist.

In robotics, rank will later determine whether a manipulator can instantaneously move in all task-space directions (full rank) or is singular (rank deficiency).

The null space (or kernel) of \( A \in \mathbb{R}^{m \times n} \) is

\[ \mathcal{N}(A) = \{\,\mathbf{x} \in \mathbb{R}^n : A\mathbf{x} = \mathbf{0}\,\}. \]

It consists of all inputs mapping to zero under the linear transformation. The dimension of the null space is called the nullity of \( A \) and is denoted \( \dim \mathcal{N}(A) \).

A fundamental theorem in linear algebra is the rank–nullity theorem:

\[ \operatorname{rank}(A) + \dim \mathcal{N}(A) = n \qquad \text{for } A \in \mathbb{R}^{m \times n}. \]

Proof sketch. Perform Gaussian elimination to reduce \( A \) to a row-echelon form \( R \). Elementary row operations do not change the solution set \( \mathcal{N}(A) \) and do not change the number of pivot columns.

• The number of pivot columns in \( R \) equals \( \operatorname{rank}(A) \).
• The number of free variables equals \( n - \operatorname{rank}(A) \).
• Each free variable corresponds to one basis vector of the null space, so \( \dim \mathcal{N}(A) = n - \operatorname{rank}(A) \).

Rearranging gives the claimed identity.

Later, null spaces will describe:

• Self motions of redundant manipulators (joint motions that leave the end-effector pose unchanged).
• Internal forces in constrained dynamics that do not accelerate the robot.

Many modeling problems reduce to solving a linear system \( A\mathbf{x} = \mathbf{b} \), where \( A \in \mathbb{R}^{m \times n} \) is known, \( \mathbf{b} \in \mathbb{R}^m \) is given, and \( \mathbf{x} \in \mathbb{R}^n \) is unknown.

Basic classification:

• Underdetermined (\( m < n \)): often infinitely many solutions (if consistent).
• Square (\( m = n \)): a unique solution if \( \det(A) \neq 0 \).
• Overdetermined (\( m > n \)): may have no exact solution; we then seek a best approximate solution in the least-squares sense.

When the system is consistent (there exists at least one solution), the general solution has the form

\[ \mathbf{x} = \mathbf{x}_p + \mathbf{z}, \qquad \mathbf{z} \in \mathcal{N}(A), \]

where \( \mathbf{x}_p \) is any particular solution. The null space describes all directions along which we can move without changing \( A\mathbf{x} \).

For inconsistent systems, the least-squares solution minimizes \( \|A\mathbf{x} - \mathbf{b}\|_2^2 \). Differentiating with respect to \( \mathbf{x} \) and setting the gradient to zero yields the normal equations:

\[ A^\top A \mathbf{x} = A^\top \mathbf{b}. \]

If \( A^\top A \) is invertible (which requires \( \operatorname{rank}(A) = n \)), the unique least-squares solution is

\[ \mathbf{x} = (A^\top A)^{-1} A^\top \mathbf{b}. \]

This relation will reappear when we discuss numerical inverse kinematics based on locally linear approximations.

Consider a configuration vector \( \mathbf{q} \in \mathbb{R}^n \) of joint variables. If we linearize the relationship between configuration and some task quantity \( \mathbf{x} \in \mathbb{R}^m \) around \( \mathbf{q}_0 \), we obtain an equation of the form

\[ \Delta \mathbf{x} \approx A\,\Delta \mathbf{q}, \]

where \( A \in \mathbb{R}^{m \times n} \) is the local linear mapping (in later lessons this will become a Jacobian matrix).

• If \( \operatorname{rank}(A) = m \) and \( n \geq m \), we can realize any small task displacement \( \Delta \mathbf{x} \) by choosing \( \Delta \mathbf{q} \), but the solution may not be unique.
• If \( \operatorname{rank}(A) < m \), some directions in task space are unreachable at that configuration. The unreachable directions form a subspace orthogonal to the row space of \( A \).
• If \( \dim \mathcal{N}(A) > 0 \), then there exist nonzero joint motions \( \Delta \mathbf{q} \in \mathcal{N}(A) \) that do not change \( \mathbf{x} \) at first order.

These geometric interpretations rely entirely on the linear algebra developed in this lesson; later chapters will specialize \( A \) to concrete kinematic and dynamic matrices.

In Python, linear algebra for robotics is typically handled with NumPy and optionally SciPy, while higher-level robotics libraries (e.g. Python versions of the Robotics Toolbox) build on these. Here we demonstrate:

• Constructing matrices representing simple kinematic constraints.
• Computing rank and an approximate basis for the null space.
• A pure Python implementation of Gaussian elimination for rank and null space.

import numpy as np

# Example: simple planar constraint relating joint velocities to an end-effector velocity
# Suppose A * qdot = vx, where qdot in R^3 and vx in R^2
A = np.array([[1.0,  1.0, 0.0],
              [0.0,  1.0, 1.0]])

print("A =\n", A)

# Numerical rank using NumPy (SVD-based)
rank = np.linalg.matrix_rank(A)
print("rank(A) =", rank)

# Approximate null space using SVD
def null_space(A, tol=1e-10):
    U, S, Vt = np.linalg.svd(A)
    # Singular values close to zero correspond to null space directions
    mask = S < tol
    # Rows of Vt corresponding to small singular values span the null space
    # (Note: Vt has shape (n, n) for full SVD)
    return Vt[mask].T

N = null_space(A)
print("Null space basis (columns):\n", N)

# Check that A * z is numerically close to zero for each basis vector z
for i in range(N.shape[1]):
    z = N[:, i]
    print(f"Check A @ z_{i} =", A @ z)
      

To understand the structure better, we implement a basic row-reduction algorithm (without pivoting) to compute rank and one null-space basis in pure Python.

def rref(A, tol=1e-12):
    """
    Compute a reduced row echelon form of A (copy), returning R and pivot column indices.
    No pivoting (for simplicity); not numerically robust but fine for small examples.
    """
    A = A.astype(float).copy()
    m, n = A.shape
    pivots = []
    row = 0
    for col in range(n):
        if row >= m:
            break
        # Find a row with nonzero entry in this column
        pivot_row = None
        for r in range(row, m):
            if abs(A[r, col]) > tol:
                pivot_row = r
                break
        if pivot_row is None:
            continue  # no pivot in this column
        # Swap into position
        if pivot_row != row:
            A[[row, pivot_row], :] = A[[pivot_row, row], :]
        # Normalize pivot row
        A[row, :] /= A[row, col]
        # Eliminate other entries in this column
        for r in range(m):
            if r != row:
                A[r, :] -= A[r, col] * A[row, :]
        pivots.append(col)
        row += 1
    return A, pivots

A = np.array([[1.0, 1.0, 0.0],
              [0.0, 1.0, 1.0]])
R, pivots = rref(A)
print("RREF(A) =\n", R)
print("Pivot columns =", pivots)
rank = len(pivots)

# Construct a null-space basis by setting free variables to parameters
m, n = A.shape
free_cols = [j for j in range(n) if j not in pivots]
basis = []
for free in free_cols:
    z = np.zeros(n)
    z[free] = 1.0
    # Solve for pivot variables from RREF: pivot columns are leading 1's
    for i, col in enumerate(pivots):
        z[col] = -R[i, free]
    basis.append(z)

N_manual = np.column_stack(basis) if basis else np.zeros((n, 0))
print("Manual null-space basis:\n", N_manual)
      

These routines will be useful when verifying analytic rank and null-space computations for small Jacobian- or constraint-like matrices later in the course.

In C++ robotics code (for example in ROS-based stacks), the Eigen library is a de-facto standard for vector and matrix computations. Below we show how to compute rank and null space for a small matrix using Eigen.

#include <iostream>
#include <Eigen/Dense>

using Eigen::MatrixXd;
using Eigen::VectorXd;

int main() {
    MatrixXd A(2, 3);
    A << 1.0, 1.0, 0.0,
          0.0, 1.0, 1.0;

    std::cout << "A =\n" << A << std::endl;

    // Compute SVD for rank and null space
    Eigen::JacobiSVD<MatrixXd> svd(
        A, Eigen::ComputeFullU | Eigen::ComputeFullV);
    auto S = svd.singularValues();
    auto V = svd.matrixV();

    double tol = 1e-10;
    int n = A.cols();
    int rank = 0;
    for (int i = 0; i < S.size(); ++i) {
        if (S(i) > tol) {
            rank++;
        }
    }
    std::cout << "rank(A) = " << rank << std::endl;

    // Columns of V corresponding to zero singular values span the null space
    std::cout << "Approximate null-space basis:\n";
    for (int i = 0; i < n; ++i) {
        if (S(i) <= tol) {
            std::cout << "basis vector " << i << ":\n"
                      << V.col(i) << "\n";
        }
    }
    return 0;
}
      

These primitives will be directly reusable when implementing kinematic and dynamic algorithms in C++ later in the course, where Eigen matrices will represent transforms, Jacobians, and inertia tensors.

For Java-based robotics or simulation tools, popular linear algebra libraries include EJML and Apache Commons Math. The following example uses EJML-style pseudo-code to demonstrate rank and null space computations.

// Using EJML (http://ejml.org) style
import org.ejml.simple.SimpleMatrix;

public class LinearAlgebraDemo {
    public static void main(String[] args) {
        double[][] dataA = {
            {1.0, 1.0, 0.0},
            {0.0, 1.0, 1.0}
        };
        SimpleMatrix A = new SimpleMatrix(dataA);

        System.out.println("A =");
        A.print();

        // SVD for rank and null space
        var svd = A.svd();
        SimpleMatrix W = svd.getW(); // singular values on diagonal
        SimpleMatrix V = svd.getV();

        double tol = 1e-10;
        int rank = 0;
        int n = A.numCols();
        for (int i = 0; i < n; ++i) {
            double sigma = W.get(i, i);
            if (sigma > tol) {
                rank++;
            }
        }
        System.out.println("rank(A) = " + rank);

        System.out.println("Approximate null-space basis:");
        for (int i = 0; i < n; ++i) {
            double sigma = W.get(i, i);
            if (sigma <= tol) {
                System.out.println("basis vector " + i + ":");
                V.extractVector(false, i).print();
            }
        }
    }
}
      

Even though the syntax differs from C++, the underlying concepts (SVD, rank, null space) are identical, and the same reasoning will carry over to Java-based robotics software.

MATLAB is widely used for robotics teaching and research. Its syntax is designed around matrices and vectors, and Simulink is naturally suited for block-diagram representations of linear systems.

% Define the same matrix A
A = [1 1 0;
     0 1 1];

disp('A =');
disp(A);

% Rank and null space
rA = rank(A);
disp(['rank(A) = ' num2str(rA)]);

N = null(A, 'r');  % 'r' requests a rational (exact) basis when possible
disp('Null-space basis (columns):');
disp(N);

% Verify that A*N is (numerically) zero
disp('A * N =');
disp(A * N);

% Example: solve a least-squares problem A*x ≈ b
b = [1; 2];
x_ls = A\b;   % MATLAB backslash chooses least-squares for overdetermined systems
disp('Least-squares solution x_ls =');
disp(x_ls);

% Simulink interpretation:
% A simple block diagram with a Constant block for x, a Gain block for A,
% and a Sum block for subtracting b represents the error A*x - b.
      

Later, MATLAB/Simulink models will integrate these blocks with kinematic and dynamic subsystems, but at their core they remain linear algebra operations on vectors and matrices.

Wolfram Mathematica is powerful for symbolic derivations of kinematics and dynamics. Here we focus on basic linear algebra functions that will be reused later.

(* Define a symbolic 2x3 matrix and compute its rank and null space *)
A = { {1, 1, 0},
     {0, 1, 1} };

MatrixForm[A]

rankA = MatrixRank[A]
nullA  = NullSpace[A]

(* Verify that A . z == {0,0} for each basis vector z in nullA *)
Table[A . z, {z, nullA}]

(* Least-squares solution of A.x == b *)
b = {1, 2};
xLS = LinearSolve[Transpose[A].A, Transpose[A].b]

(* Alternatively, use the symbolic pseudo-inverse *)
xLS2 = PseudoInverse[A].b
      

Symbolic manipulation of these matrices will later support deriving closed-form expressions for simple robot kinematics and dynamics, which can then be validated numerically.

In this lesson we reviewed the linear algebra objects that underlie all subsequent kinematics and dynamics modeling:

• Vectors and inner products define geometric notions of angle, length, and projection.
• Matrices represent linear maps; their columns span the range, and their rank measures the dimension of that range.
• The null space captures directions that map to zero, and the rank–nullity theorem organizes the dimensions of domain, range, and null space.
• Solving \( A\mathbf{x} = \mathbf{b} \) relies on rank conditions; least-squares solutions arise from normal equations.
• Robot kinematic and dynamic equations will be expressed in exactly this language, with matrices encoding coupling between coordinates, velocities, and forces.
• Practical implementations in Python (NumPy/SciPy), C++ (Eigen), Java (EJML), MATLAB/Simulink, and Mathematica all share the same underlying linear algebra.

The next lesson builds on this foundation to introduce eigenvalues and eigenvectors, which will provide deeper insight into the structure of matrices that arise in robotics (e.g. inertia matrices and linearized dynamics).