Consider an \( n \)-DOF serial robot manipulator with joint torque vector \( \boldsymbol{\tau} \in \mathbb{R}^n \) and end-effector wrench \( \mathbf{w} \in \mathbb{R}^m \) (typically \( m=6 \) in 3D, or \( m=3 \) in planar cases). From virtual work and manipulator statics (previous lessons), the static equilibrium relation in configuration \( \mathbf{q} \) is

\[ \boldsymbol{\tau} = \mathbf{J}^T(\mathbf{q})\,\mathbf{w}, \]

where \( \mathbf{J}(\mathbf{q}) \) is the (geometric) Jacobian. Real actuators are limited: each joint has saturation bounds

\[ -\tau_i^{\max} \le \tau_i \le \tau_i^{\max},\quad i=1,\dots,n, \]

collected in the vector of limits \( \boldsymbol{\tau}^{\max} \in \mathbb{R}^n_{>0} \). The end-effector load capacity at configuration \( \mathbf{q} \) is the set of all wrenches \( \mathbf{w} \) for which static equilibrium is possible without violating actuator bounds.

Formally, the admissible joint-torque set is a hyper-rectangle (box)

\[ \mathcal{T} \;=\; \left\{ \boldsymbol{\tau} \in \mathbb{R}^n \,\middle|\, -\boldsymbol{\tau}^{\max} \le \boldsymbol{\tau} \le \boldsymbol{\tau}^{\max} \right\}, \]

and, using \( \boldsymbol{\tau} = \mathbf{J}^T \mathbf{w} \), the wrench capacity set is the pre-image

\[ \mathcal{W}(\mathbf{q}) \;=\; \left\{ \mathbf{w} \in \mathbb{R}^m \,\middle|\, \mathbf{J}^T(\mathbf{q})\,\mathbf{w} \in \mathcal{T} \right\}. \]

This set completely describes which external loads the manipulator can statically sustain at that pose. Our goal is to characterize \( \mathcal{W}(\mathbf{q}) \) analytically and algorithmically.

Substituting \( \boldsymbol{\tau} = \mathbf{J}^T\mathbf{w} \) into the box constraints, we obtain a system of linear inequalities in \( \mathbf{w} \):

\[ -\boldsymbol{\tau}^{\max} \le \mathbf{J}^T(\mathbf{q})\,\mathbf{w} \le \boldsymbol{\tau}^{\max}, \]

interpreted component-wise. Equivalently, for each joint \( i \) the feasible wrenches satisfy

\[ -\tau_i^{\max} \le \mathbf{j}_i^T(\mathbf{q})\,\mathbf{w} \le \tau_i^{\max}, \quad \mathbf{j}_i(\mathbf{q}) \in \mathbb{R}^m \]

where \( \mathbf{j}_i \) is the \( i \)-th column of the Jacobian \( \mathbf{J} \). Thus

\[ \mathcal{W}(\mathbf{q}) = \bigcap_{i=1}^n \left\{ \mathbf{w} \in \mathbb{R}^m \;\middle|\; -\tau_i^{\max} \le \mathbf{j}_i^T \mathbf{w} \le \tau_i^{\max} \right\}. \]

Each pair of inequalities for joint \( i \) defines a pair of parallel half-spaces in wrench space. The intersection of finitely many half-spaces is a convex polyhedron; because all inequalities are bounded in both directions, we obtain a convex, centrally symmetric polytope.

Proposition (Convexity and Symmetry). The wrench capacity set \( \mathcal{W}(\mathbf{q}) \) is convex, compact, and satisfies \( \mathbf{w} \in \mathcal{W}(\mathbf{q}) \Rightarrow -\mathbf{w} \in \mathcal{W}(\mathbf{q}) \).

Sketch of proof.

• Convexity: each inequality of the form \( a^T \mathbf{w} \le b \) defines a convex set. The intersection of convex sets is convex, so \( \mathcal{W}(\mathbf{q}) \) is convex.
• Compactness: the pre-image of the bounded set \( \mathcal{T} \) under linear map \( \mathbf{J}^T \) is bounded whenever \( \mathbf{J} \) has full row rank in the wrench subspace of interest; closedness follows from continuity.
• Symmetry: if \( \mathbf{w} \) satisfies \( -\boldsymbol{\tau}^{\max} \le \mathbf{J}^T\mathbf{w} \le \boldsymbol{\tau}^{\max} \), then multiplying by \( -1 \) yields \( -\boldsymbol{\tau}^{\max} \le \mathbf{J}^T(-\mathbf{w}) \le \boldsymbol{\tau}^{\max} \), so \( -\mathbf{w} \in \mathcal{W}(\mathbf{q}) \).

The polytope \( \mathcal{W}(\mathbf{q}) \) is sometimes called the wrench polytope or static load capacity polytope. Its shape reflects:

• the geometry (columns) of \( \mathbf{J}(\mathbf{q}) \), and
• the anisotropy of actuator limits \( \boldsymbol{\tau}^{\max} \).

The high-level steps to construct a static load capacity model at a given configuration are summarized below.

flowchart TD
  Q["Manipulator model (q, J(q))"] --> L["Joint torque limits 'tau_max'"]
  Q --> S["Static relation 'tau = J^T w'"]
  L --> C["Feasible torque box T"]
  S --> W["Feasible wrench set W"]
  C --> W
  W --> I["Load indices, maps, and safety checks"]
        

In practice, analytic descriptions (half-space form) can be complemented with numerical representations of \( \mathcal{W}(\mathbf{q}) \) (vertex enumeration, sampling, or support functions) depending on downstream use (e.g., planning, design, or safety verification).

While the full polytope \( \mathcal{W}(\mathbf{q}) \) is useful conceptually, many design questions concern the maximum load in a given direction. Let \( \hat{\mathbf{w}} \in \mathbb{S}^{m-1} \) be a unit wrench direction (e.g., vertical downwards force). We define the directional load capacity as

\[ \lambda^{\ast}(\hat{\mathbf{w}}) \;=\; \sup \left\{ \lambda \ge 0 \,\middle|\, \lambda \hat{\mathbf{w}} \in \mathcal{W}(\mathbf{q}) \right\}. \]

Using \( \boldsymbol{\tau} = \mathbf{J}^T(\mathbf{q}) \,\lambda \hat{\mathbf{w}} \), the actuator bounds become

\[ -\tau_i^{\max} \le \lambda\, \mathbf{j}_i^T \hat{\mathbf{w}} \le \tau_i^{\max}, \quad i=1,\dots,n. \]

Taking absolute values and exploiting symmetry yields the scalar inequalities

\[ \lambda\,\left| \mathbf{j}_i^T \hat{\mathbf{w}} \right| \le \tau_i^{\max}, \quad i=1,\dots,n. \]

For any joint \( i \) such that \( \left| \mathbf{j}_i^T \hat{\mathbf{w}} \right| > 0 \), we obtain the upper bound \( \lambda \le \tau_i^{\max} / \left| \mathbf{j}_i^T \hat{\mathbf{w}} \right| \). Joints with \( \mathbf{j}_i^T \hat{\mathbf{w}} = 0 \) do not restrict this particular direction.

Proposition (Directional capacity under box torque limits). Assume at least one joint has \( \mathbf{j}_i^T \hat{\mathbf{w}} \neq 0 \). Then

\[ \lambda^{\ast}(\hat{\mathbf{w}}) \;=\; \min_{i: \,\mathbf{j}_i^T \hat{\mathbf{w}} \neq 0} \frac{\tau_i^{\max}}{\left| \mathbf{j}_i^T \hat{\mathbf{w}} \right|}. \]

Proof.

• Feasibility requires that all inequalities \( \lambda\,\left| \mathbf{j}_i^T \hat{\mathbf{w}} \right| \le \tau_i^{\max} \) hold, so \( \lambda \) cannot exceed the minimum of the right-hand sides.
• For any \( \lambda \) below this minimum, all inequalities are strictly satisfied, hence \( \lambda \hat{\mathbf{w}} \in \mathcal{W}(\mathbf{q}) \). Thus the supremum equals the minimum bound.

This formula provides an efficient way to compute directional capacity: only the Jacobian column projections on \( \hat{\mathbf{w}} \) and the torque limits enter.

The box model reflects independent joint limits. Another common model uses a quadratic bound on the joint torque norm (e.g., when amplifiers or thermal constraints couple joints):

\[ \left\| \boldsymbol{\tau} \right\|_2^2 \;=\; \boldsymbol{\tau}^T \boldsymbol{\tau} \le (\tau_{\max})^2. \]

Substituting \( \boldsymbol{\tau} = \mathbf{J}^T \mathbf{w} \) yields the wrench constraint

\[ \left\| \mathbf{J}^T(\mathbf{q})\,\mathbf{w} \right\|_2^2 \le (\tau_{\max})^2, \]

or, equivalently,

\[ \mathbf{w}^T \,\mathbf{J}(\mathbf{q})\,\mathbf{J}^T(\mathbf{q})\,\mathbf{w} \le (\tau_{\max})^2. \]

Provided \( \mathbf{J}(\mathbf{q})\,\mathbf{J}^T(\mathbf{q}) \) is positive definite in the wrench subspace, the set

\[ \mathcal{W}_{\mathrm{ell}}(\mathbf{q}) \;=\; \left\{ \mathbf{w} \in \mathbb{R}^m \,\middle|\, \mathbf{w}^T \mathbf{J}\mathbf{J}^T \mathbf{w} \le (\tau_{\max})^2 \right\} \]

is an origin-centered ellipsoid in wrench space, sometimes called the force manipulability ellipsoid.

The corresponding directional capacity is

\[ \lambda^{\ast}_{\mathrm{ell}}(\hat{\mathbf{w}}) = \sup \left\{ \lambda \ge 0 \,\middle|\, \lambda^2 \hat{\mathbf{w}}^T \mathbf{J}\mathbf{J}^T \hat{\mathbf{w}} \le (\tau_{\max})^2 \right\} = \frac{\tau_{\max}}{\left\| \mathbf{J}^T \hat{\mathbf{w}} \right\|_2}. \]

If \( \sigma_1 \ge \dots \ge \sigma_m > 0 \) are the singular values of \( \mathbf{J} \), then the principal semi-axes of \( \mathcal{W}_{\mathrm{ell}} \) have lengths \( \tau_{\max} / \sigma_i \). Thus near kinematic singularities (small \( \sigma_i \)) the ellipsoid is extremely elongated and directional load capacity degrades sharply in certain directions.

In summary:

• The box model yields a polytope \( \mathcal{W}(\mathbf{q}) \).
• The norm model yields an ellipsoid \( \mathcal{W}_{\mathrm{ell}}(\mathbf{q}) \).
• Both are consistent with the same static mapping \( \boldsymbol{\tau} = \mathbf{J}^T\mathbf{w} \); the difference lies in how actuator constraints are modeled.

Consider a planar 2R robot with link lengths \( \ell_1, \ell_2 \) and joint angles \( q_1, q_2 \). The end-effector position is

\[ x = \ell_1 \cos q_1 + \ell_2 \cos(q_1+q_2), \quad y = \ell_1 \sin q_1 + \ell_2 \sin(q_1+q_2). \]

The planar translational Jacobian (mapping joint velocities to end-effector translational velocity) is

\[ \mathbf{J}(q_1,q_2) = \begin{bmatrix} \dfrac{\partial x}{\partial q_1} & \dfrac{\partial x}{\partial q_2} \\ \dfrac{\partial y}{\partial q_1} & \dfrac{\partial y}{\partial q_2} \end{bmatrix} = \begin{bmatrix} -\ell_1 \sin q_1 - \ell_2 \sin(q_1+q_2) & -\ell_2 \sin(q_1+q_2) \\ \ell_1 \cos q_1 + \ell_2 \cos(q_1+q_2) & \ell_2 \cos(q_1+q_2) \end{bmatrix}. \]

Let the external planar force be \( \mathbf{f} = [f_x, f_y]^T \). The static relation is

\[ \begin{bmatrix} \tau_1 \\ \tau_2 \end{bmatrix} = \mathbf{J}^T(q_1,q_2) \begin{bmatrix} f_x \\ f_y \end{bmatrix}. \]

Suppose we are interested in a pure vertical downward load \( \hat{\mathbf{f}} = [0, -1]^T \), and we consider loads of the form \( \mathbf{f} = \lambda \hat{\mathbf{f}} \). Then

\[ \begin{bmatrix} \tau_1 \\ \tau_2 \end{bmatrix} = \lambda \mathbf{J}^T(q_1,q_2) \begin{bmatrix} 0 \\ -1 \end{bmatrix} = \lambda \begin{bmatrix} -\dfrac{\partial y}{\partial q_1} \\ -\dfrac{\partial y}{\partial q_2} \end{bmatrix} = -\lambda \begin{bmatrix} \ell_1 \cos q_1 + \ell_2 \cos(q_1+q_2) \\ \ell_2 \cos(q_1+q_2) \end{bmatrix}. \]

With symmetric torque limits \( |\tau_1| \le \tau_1^{\max} \), \( |\tau_2| \le \tau_2^{\max} \), we obtain

\[ \lambda \le \frac{\tau_1^{\max}}{\left| \ell_1 \cos q_1 + \ell_2 \cos(q_1+q_2) \right|},\quad \lambda \le \frac{\tau_2^{\max}}{\left| \ell_2 \cos(q_1+q_2) \right|}. \]

Hence the directional capacity in the vertical downward direction is

\[ \lambda^{\ast}(\hat{\mathbf{f}}) = \min\!\left( \frac{\tau_1^{\max}}{\left| \ell_1 \cos q_1 + \ell_2 \cos(q_1+q_2) \right|}, \frac{\tau_2^{\max}}{\left| \ell_2 \cos(q_1+q_2) \right|} \right). \]

As a concrete numeric example, take \( \ell_1 = \ell_2 = 1 \) and \( q_1 = 0, q_2 = 0 \) (arm fully stretched horizontally). Then \( \cos q_1 = 1 \), \( \cos(q_1+q_2) = 1 \), so

\[ \ell_1 \cos q_1 + \ell_2 \cos(q_1+q_2) = 2, \quad \ell_2 \cos(q_1+q_2) = 1. \]

If \( \tau_1^{\max} = \tau_2^{\max} = \tau_{\max} \), we obtain

\[ \lambda^{\ast}(\hat{\mathbf{f}}) = \min\!\left( \frac{\tau_{\max}}{2}, \tau_{\max} \right) = \frac{\tau_{\max}}{2}. \]

Thus the maximum sustainable vertical load at this configuration is \( \tau_{\max}/2 \), revealing a strong lever-arm effect on joint 1. Closer to a folded configuration, the cosine terms change and the capacity in the same direction can increase or decrease significantly.

The formula from Section 4 leads to a simple algorithm for evaluating directional capacity under box torque limits. This is illustrated below.

flowchart TD
  D["Input: J(q), tau_max, direction d"] --> N["Normalize d to unit 'hat_d'"]
  N --> V["Compute v = J^T * hat_d"]
  V --> L1["For each joint i: if v_i != 0, lambda_i = tau_max_i / |v_i|"]
  L1 --> M["lambda* = min over i of lambda_i"]
  M --> O["Return lambda* and active limiting joints"]
        

This procedure can be repeated over a set of directions \( \hat{\mathbf{w}} \) (e.g., on a sphere grid) to approximate the entire wrench polytope by its support function, which is often sufficient for design and safety analyses.

We now implement the directional load capacity for the planar 2R arm in Python using numpy. We restrict to planar forces \( \mathbf{f} \in \mathbb{R}^2 \) and joint torque limits \( \tau_1^{\max}, \tau_2^{\max} \).

import numpy as np

def planar_2r_jacobian(q1, q2, l1, l2):
    """
    Return the 2x2 translational Jacobian J(q) for a planar 2R arm.
    """
    s1 = np.sin(q1)
    c1 = np.cos(q1)
    s12 = np.sin(q1 + q2)
    c12 = np.cos(q1 + q2)

    j11 = -l1 * s1 - l2 * s12
    j12 = -l2 * s12
    j21 =  l1 * c1 + l2 * c12
    j22 =  l2 * c12
    J = np.array([[j11, j12],
                  [j21, j22]], dtype=float)
    return J

def directional_capacity_box(J, f_dir, tau_max, eps=1e-9):
    """
    Directional load capacity under box torque limits.

    Parameters
    ----------
    J : (2, 2) ndarray
        Planar Jacobian at configuration q.
    f_dir : (2,) ndarray
        Desired force direction (need not be unit length).
    tau_max : (2,) ndarray
        Positive torque limits [tau1_max, tau2_max].
    eps : float
        Tolerance for treating small Jacobian projections as zero.

    Returns
    -------
    lam_star : float
        Maximum scalar lambda such that lambda * f_hat is feasible.
    f_hat : (2,) ndarray
        Normalized direction vector.
    active_joints : list of int
        Indices of joints that are limiting at lam_star.
    """
    f_dir = np.asarray(f_dir, dtype=float)
    tau_max = np.asarray(tau_max, dtype=float)
    norm = np.linalg.norm(f_dir)
    if norm < eps:
        raise ValueError("Direction vector f_dir must be nonzero.")
    f_hat = f_dir / norm

    v = J.T @ f_hat  # mapping from unit force to joint torques
    lambdas = []
    active_indices = []
    for i, v_i in enumerate(v):
        if abs(v_i) < eps:
            # This joint does not constrain this direction
            continue
        lam_i = tau_max[i] / abs(v_i)
        lambdas.append(lam_i)
        active_indices.append(i)

    if not lambdas:
        # In this (rare) case, the direction is unobservable in joint torques
        # and is not constrained by the box model.
        return np.inf, f_hat, []

    lam_star = min(lambdas)
    # Joints achieving the minimum are limiting
    limiting = [active_indices[i]
                for i, lam_i in enumerate(lambdas)
                if abs(lam_i - lam_star) < 1e-6]

    return lam_star, f_hat, limiting

if __name__ == "__main__":
    # Example: l1 = l2 = 1, q1 = 0, q2 = 0, symmetric torque limits
    l1 = l2 = 1.0
    q1 = 0.0
    q2 = 0.0
    J = planar_2r_jacobian(q1, q2, l1, l2)
    tau_max = np.array([10.0, 10.0])  # arbitrary units
    f_dir = np.array([0.0, -1.0])     # downward

    lam_star, f_hat, limiting = directional_capacity_box(J, f_dir, tau_max)
    print("J(q) =\n", J)
    print("Unit direction f_hat:", f_hat)
    print("Directional capacity lambda*:", lam_star)
    print("Limiting joints:", limiting)
      

This code implements the theoretical formula in Section 4 directly and scales easily to higher-dimensional Jacobians with minor modifications.

In C++, we can use the Eigen library for linear algebra and implement an analogous function for 2-DOF manipulators. The same logic extends to \( n \)-DOF arms using VectorXd and MatrixXd.

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

using Eigen::Matrix2d;
using Eigen::Vector2d;

Matrix2d planar2RJacobian(double q1, double q2, double l1, double l2) {
    double s1  = std::sin(q1);
    double c1  = std::cos(q1);
    double s12 = std::sin(q1 + q2);
    double c12 = std::cos(q1 + q2);

    Matrix2d J;
    J(0,0) = -l1 * s1 - l2 * s12;
    J(0,1) = -l2 * s12;
    J(1,0) =  l1 * c1 + l2 * c12;
    J(1,1) =  l2 * c12;
    return J;
}

double directionalCapacityBox(const Matrix2d& J,
                              const Vector2d& f_dir,
                              const Vector2d& tau_max,
                              double eps = 1e-9) {
    double norm = f_dir.norm();
    if (norm < eps) {
        throw std::runtime_error("Direction vector must be nonzero.");
    }
    Vector2d f_hat = f_dir / norm;
    Vector2d v = J.transpose() * f_hat;

    double lambda_star = std::numeric_limits<double>::infinity();
    for (int i = 0; i < 2; ++i) {
        double v_i = v(i);
        if (std::abs(v_i) < eps) {
            continue;  // does not constrain this direction
        }
        double lambda_i = tau_max(i) / std::abs(v_i);
        if (lambda_i < lambda_star) {
            lambda_star = lambda_i;
        }
    }
    return lambda_star;
}

int main() {
    double l1 = 1.0, l2 = 1.0;
    double q1 = 0.0, q2 = 0.0;

    Matrix2d J = planar2RJacobian(q1, q2, l1, l2);
    Vector2d tau_max(10.0, 10.0);
    Vector2d f_dir(0.0, -1.0);

    double lambda_star = directionalCapacityBox(J, f_dir, tau_max);
    std::cout << "Directional capacity lambda* = " << lambda_star << std::endl;
    return 0;
}
      

For larger robots, one would replace Matrix2d and Vector2d by Eigen::MatrixXd and Eigen::VectorXd, and loop over all joint indices.

10.1 Java Implementation (Using EJML or Plain Arrays)

In Java, a convenient linear algebra library is EJML, but for a 2x2 example we can use plain arrays:

public class Planar2RLoadCapacity {

    public static double[][] jacobian(double q1, double q2, double l1, double l2) {
        double s1  = Math.sin(q1);
        double c1  = Math.cos(q1);
        double s12 = Math.sin(q1 + q2);
        double c12 = Math.cos(q1 + q2);

        double[][] J = new double[2][2];
        J[0][0] = -l1 * s1 - l2 * s12;
        J[0][1] = -l2 * s12;
        J[1][0] =  l1 * c1 + l2 * c12;
        J[1][1] =  l2 * c12;
        return J;
    }

    public static double directionalCapacityBox(double[][] J,
                                                double[] fDir,
                                                double[] tauMax,
                                                double eps) {
        double norm = Math.sqrt(fDir[0] * fDir[0] + fDir[1] * fDir[1]);
        if (norm < eps) {
            throw new IllegalArgumentException("Direction vector must be nonzero.");
        }
        double[] fHat = new double[2];
        fHat[0] = fDir[0] / norm;
        fHat[1] = fDir[1] / norm;

        // v = J^T * f_hat
        double[] v = new double[2];
        v[0] = J[0][0] * fHat[0] + J[1][0] * fHat[1];
        v[1] = J[0][1] * fHat[0] + J[1][1] * fHat[1];

        double lambdaStar = Double.POSITIVE_INFINITY;
        for (int i = 0; i < 2; ++i) {
            double v_i = v[i];
            if (Math.abs(v_i) < eps) {
                continue;
            }
            double lambda_i = tauMax[i] / Math.abs(v_i);
            if (lambda_i < lambdaStar) {
                lambdaStar = lambda_i;
            }
        }
        return lambdaStar;
    }

    public static void main(String[] args) {
        double l1 = 1.0, l2 = 1.0;
        double q1 = 0.0, q2 = 0.0;
        double[][] J = jacobian(q1, q2, l1, l2);
        double[] tauMax = {10.0, 10.0};
        double[] fDir = {0.0, -1.0};

        double lambdaStar = directionalCapacityBox(J, fDir, tauMax, 1e-9);
        System.out.println("Directional capacity lambda* = " + lambdaStar);
    }
}
      

10.2 MATLAB/Simulink Implementation

In MATLAB, we can create an m-file function that encapsulates the same computation:

function [lambda_star, f_hat, limiting] = directionalCapacityBox2R(q1, q2, l1, l2, f_dir, tau_max)
%DIRECTIONALCAPACITYBOX2R  Planar 2R directional load capacity (static).
%
% Inputs:
%   q1, q2   Joint angles (rad)
%   l1, l2   Link lengths
%   f_dir    2x1 desired force direction (need not be unit length)
%   tau_max  2x1 positive torque limits [tau1_max; tau2_max]
%
% Outputs:
%   lambda_star  Maximum scalar lambda such that lambda * f_hat is feasible
%   f_hat        2x1 normalized direction vector
%   limiting     Indices of limiting joints at lambda_star

    eps = 1e-9;
    if norm(f_dir) < eps
        error('Direction vector must be nonzero.');
    end
    f_hat = f_dir / norm(f_dir);

    % Jacobian
    s1  = sin(q1);  c1  = cos(q1);
    s12 = sin(q1 + q2);  c12 = cos(q1 + q2);

    J = [ -l1*s1 - l2*s12, -l2*s12;
           l1*c1 + l2*c12,  l2*c12 ];

    v = J.' * f_hat;  % mapping to torques
    lambda_list = [];
    idx_list = [];

    for i = 1:2
        v_i = v(i);
        if abs(v_i) < eps
            continue;
        end
        lam_i = tau_max(i) / abs(v_i);
        lambda_list(end+1) = lam_i; %#ok<AGROW>
        idx_list(end+1) = i;        %#ok<AGROW>
    end

    if isempty(lambda_list)
        lambda_star = inf;
        limiting = [];
    else
        lambda_star = min(lambda_list);
        limiting = idx_list(abs(lambda_list - lambda_star) < 1e-6);
    end
end
      

In Simulink, the static relation tau = J^T * f can be implemented using a Matrix Gain block (with gain J') followed by Saturation blocks modeling the element-wise torque limits. A MATLAB Function block can then compute the directional capacity or serve as a static safety monitor.

10.3 Wolfram Mathematica Implementation

Mathematica is well-suited for symbolic derivations and numeric evaluation of load capacity:

(* Planar 2R Jacobian *)
J[q1_, q2_, l1_, l2_] := Module[
  {s1, c1, s12, c12},
  s1 = Sin[q1]; c1 = Cos[q1];
  s12 = Sin[q1 + q2]; c12 = Cos[q1 + q2];
  {
    {-l1*s1 - l2*s12, -l2*s12},
    { l1*c1 + l2*c12,  l2*c12}
  }
];

directionalCapacityBox[Jmat_, fdir_, tauMax_] := Module[
  {eps = 10.^-9, fhat, v, lambdas},
  If[Norm[fdir] < eps, Return[Infinity]];
  fhat = fdir/Norm[fdir];
  v = Transpose[Jmat].fhat;
  lambdas = Table[
    If[Abs[v[[i]]] < eps, Infinity,
      tauMax[[i]]/Abs[v[[i]]]
    ],
    {i, Length[v]}
  ];
  Min[lambdas]
];

(* Example usage *)
l1 = 1.; l2 = 1.;
q1 = 0.; q2 = 0.;
Jmat = J[q1, q2, l1, l2];
tauMax = {10., 10.};
fdir = {0., -1.};

lambdaStar = directionalCapacityBox[Jmat, fdir, tauMax]
      

Symbolic manipulations (e.g., leaving q1, q2 as symbols) can expose analytic dependence of load capacity on configuration, which is useful for design optimization or workspace analysis.

In this lesson, we modeled the static load capacity of a manipulator as the set of end-effector wrenches compatible with actuator torque limits via the static relation \( \boldsymbol{\tau} = \mathbf{J}^T \mathbf{w} \). For independent joint bounds, this yields a convex wrench polytope; for norm-bounded torques, an ellipsoidal wrench set emerges. We derived closed-form expressions for directional load capacity and illustrated configuration-dependent effects on a 2R planar arm. Finally, we implemented these models in Python, C++, Java, MATLAB/Simulink, and Mathematica, preparing the ground for constrained equilibrium analysis in the next lesson.