Real robot controllers run on digital hardware: encoder measurements are sampled at discrete instants, control laws are evaluated numerically, and actuator commands are held constant between updates by power electronics (PWM, current loops, amplifiers). At the modeling level, we must relate a continuous-time robot model to a discrete-time model consistent with:

• A sampling period \( T \) (e.g. \( T = 1\text{ ms} \)).
• A sample-and-hold structure: the input is piecewise constant on each interval.
• Digital controller updates \( \mathbf{u}[k] \) at times \( t_k = kT \).

For a linearized robot joint or manipulator model we use the continuous-time state-space form

\[ \dot{\mathbf{x}}(t) = A \mathbf{x}(t) + B \mathbf{u}(t), \quad \mathbf{y}(t) = C \mathbf{x}(t) + D \mathbf{u}(t), \]

where \( \mathbf{x}(t) \) typically stacks joint positions and velocities, \( \mathbf{u}(t) \) are torques or motor voltages, and \( \mathbf{y}(t) \) are measured outputs. Under sampling and ZOH, we will derive a discrete-time model

\[ \mathbf{x}[k+1] = \Phi \, \mathbf{x}[k] + \Gamma \, \mathbf{u}[k], \quad \mathbf{y}[k] = C \mathbf{x}[k] + D \mathbf{u}[k], \]

where \( \Phi \) and \( \Gamma \) depend on \( A,B \) and the sampling period \( T \).

flowchart TD
  R["Reference r[k]"] --> E["Error computation"]
  Y["Measured output y[k]"] --> E
  E --> Ctr["Digital controller (software)"]
  Ctr --> Uk["Control command u[k]"]
  Uk --> ZOH["Zero order hold (power electronics)"]
  ZOH --> Plant["Robot joint/arm (continuous plant)"]
  Plant --> Ys["Continuous output y(t)"]
  Ys --> Samp["Sampler (sensors)"]
  Samp --> Y
        

The block diagram emphasizes that only the plant and actuators evolve continuously. Everything else is inherently discrete and must be modeled with sampled-data mathematics.

Let \( T > 0 \) be the sampling period. Sampling instants are

\[ t_k = k T, \quad k \in \mathbb{Z}_{\ge 0}. \]

For any continuous-time signal \( \mathbf{v}(t) \) (e.g. joint position, velocity, torque), its sampled sequence is

\[ \mathbf{v}[k] \triangleq \mathbf{v}(t_k) = \mathbf{v}(kT). \]

In a digital controller, the state, input, and output sequences are defined as

\[ \mathbf{x}[k] \triangleq \mathbf{x}(kT), \quad \mathbf{u}[k] \triangleq \mathbf{u}(kT^-), \quad \mathbf{y}[k] \triangleq \mathbf{y}(kT), \]

where \( kT^- \) denotes the left limit at \( t = kT \) (relevant when a ZOH changes the input exactly at \( t_k \)).

The choice of \( T \) is constrained by:

• Sensing and actuation hardware (maximum feasible rate).
• Computation time of the control law on the CPU.
• Bandwidth of the robot dynamics: typically \( T \) is chosen such that the sampling frequency \( \omega_s = 2\pi / T \) is much larger than the dominant closed-loop bandwidth (rule-of-thumb factors between 10 and 30).

In later lessons we will quantify how \( T \) affects stability and performance; here we focus on deriving the correct discrete-time model for any given \( T \).

A zero-order hold (ZOH) reconstructs a continuous-time signal from a sequence \( \mathbf{u}[k] \) by holding each value constant over one sampling interval:

\[ \mathbf{u}(t) = \mathbf{u}[k], \quad kT \le t < (k+1)T, \quad k = 0,1,2,\dots \]

This is an excellent model of many actuator interfaces in robotics: a DAC or PWM module keeps the commanded value fixed until the next update. The ZOH has a continuous-time transfer function (from the discrete impulse sequence to the held signal) given, in Laplace domain, by

\[ G_{\text{ZOH}}(s) = \frac{1 - e^{-sT}}{s}. \]

More important for state-space analysis is how the ZOH affects the plant dynamics. Over each interval \( [kT,(k+1)T) \), the input is constant, so the continuous-time system

\[ \dot{\mathbf{x}}(t) = A \mathbf{x}(t) + B \mathbf{u}[k], \quad kT \le t < (k+1)T, \]

is a linear time-invariant system with constant forcing. This allows closed-form solutions and leads to an exact discrete-time model, without numerical approximations.

Consider again the continuous-time linear system

\[ \dot{\mathbf{x}}(t) = A \mathbf{x}(t) + B \mathbf{u}(t). \]

The general solution between times \( t_0 \) and \( t \) is

\[ \mathbf{x}(t) = e^{A(t - t_0)} \mathbf{x}(t_0) + \int_{t_0}^{t} e^{A(t - \sigma)} B \mathbf{u}(\sigma)\, d\sigma. \]

Set \( t_0 = kT \) and \( t = (k+1)T \). Under ZOH, \( \mathbf{u}(\sigma) = \mathbf{u}[k] \) for \( \sigma \in [kT,(k+1)T) \). Then

\[ \mathbf{x}((k+1)T) = e^{AT} \mathbf{x}(kT) + \left( \int_{kT}^{(k+1)T} e^{A((k+1)T - \sigma)} B \, d\sigma \right) \mathbf{u}[k]. \]

Substitute \( \sigma = kT + \tau \) with \( \tau \in [0,T] \), so \( (k+1)T - \sigma = T - \tau \):

\[ \int_{kT}^{(k+1)T} e^{A((k+1)T - \sigma)} B \, d\sigma = \int_{0}^{T} e^{A(T - \tau)} B \, d\tau. \]

Because \( e^{A(T - \tau)} = e^{AT} e^{-A\tau} \), we can factor out \( e^{AT} \):

\[ \int_{0}^{T} e^{A(T - \tau)} B \, d\tau = e^{AT} \int_{0}^{T} e^{-A\tau} B \, d\tau. \]

Defining the discrete-time state as \( \mathbf{x}[k] = \mathbf{x}(kT) \), we obtain the exact discrete-time model

\[ \mathbf{x}[k+1] = \Phi \, \mathbf{x}[k] + \Gamma \, \mathbf{u}[k], \]

with

\[ \Phi = e^{AT}, \qquad \Gamma = \int_{0}^{T} e^{A\tau} B \, d\tau. \]

When \( A \) is nonsingular, an equivalent closed form for \( \Gamma \) is obtained by integrating by parts:

\[ \Gamma = \left( \int_{0}^{T} e^{A\tau} \, d\tau \right) B = A^{-1}(e^{AT} - I) B. \]

4.1 Scalar example: first-order actuator model

Let a scalar plant model a first-order actuator:

\[ \dot{x}(t) = a x(t) + b u(t), \quad x(0) = x_0. \]

The exact discrete-time model under ZOH with period \( T \) is:

• State transition:

  \[ \Phi = e^{aT}. \]

• Input gain (for \( a \neq 0 \)):

  \[ \Gamma = \int_{0}^{T} e^{a\tau} b \, d\tau = \frac{b}{a}(e^{aT} - 1). \]

Hence

\[ x[k+1] = e^{aT} x[k] + \frac{b}{a}(e^{aT} - 1) u[k]. \]

For \( a = 0 \), the continuous-time system is \( \dot{x} = b u \), so \( x((k+1)T) - x(kT) = b T u[k] \) and

\[ x[k+1] = x[k] + b T u[k], \quad \Phi = 1, \; \Gamma = bT. \]

4.2 Linearized robot joint model

A common linear model for a single revolute joint around an equilibrium is

\[ \mathbf{x}(t) = \begin{bmatrix} q(t) \\ \dot{q}(t) \end{bmatrix}, \quad \dot{\mathbf{x}}(t) = \underbrace{ \begin{bmatrix} 0 & 1 \\ 0 & -\frac{b}{J} \end{bmatrix} }_{A} \mathbf{x}(t) + \underbrace{ \begin{bmatrix} 0 \\ \frac{1}{J} \end{bmatrix} }_{B} u(t), \]

where \( J \) is the equivalent inertia and \( b \) is viscous damping. For a sampling period \( T \), the exact discrete matrices are

\[ \Phi = e^{AT}, \qquad \Gamma = \int_{0}^{T} e^{A\tau} B \, d\tau, \]

which can be evaluated analytically for this \( 2 \times 2 \) system or numerically using matrix exponentials. These discrete matrices are then used in all subsequent discrete-time control-design steps.

Exact discretization via \( e^{AT} \) and the integral for \( \Gamma \) is theoretically clean but sometimes computationally heavy on embedded hardware, or difficult to derive by hand for complex models. Several approximate schemes are widely used. They can be interpreted as numerical integration methods applied to \( \dot{\mathbf{x}} = A \mathbf{x} + B \mathbf{u} \).

5.1 Forward Euler (explicit Euler)

Approximate the derivative at the left endpoint:

\[ \dot{\mathbf{x}}(t) \approx \frac{\mathbf{x}((k+1)T) - \mathbf{x}(kT)}{T}, \quad t \approx kT. \]

Substituting into the state equation with constant input \( \mathbf{u}(t) = \mathbf{u}[k] \) yields

\[ \frac{\mathbf{x}((k+1)T) - \mathbf{x}(kT)}{T} \approx A \mathbf{x}(kT) + B \mathbf{u}[k]. \]

Therefore

\[ \mathbf{x}[k+1] \approx (I + T A) \mathbf{x}[k] + T B \mathbf{u}[k], \]

i.e.

\[ \Phi_{\text{FE}} = I + T A, \quad \Gamma_{\text{FE}} = T B. \]

Forward Euler is first-order accurate in \( T \), i.e. the local truncation error is \( \mathcal{O}(T^2) \).

5.2 Backward Euler (implicit Euler)

Approximate the derivative at the right endpoint:

\[ \dot{\mathbf{x}}(t) \approx \frac{\mathbf{x}((k+1)T) - \mathbf{x}(kT)}{T}, \quad t \approx (k+1)T. \]

Now plug in \( \mathbf{x}((k+1)T) \) into the right-hand side:

\[ \frac{\mathbf{x}((k+1)T) - \mathbf{x}(kT)}{T} \approx A \mathbf{x}((k+1)T) + B \mathbf{u}[k]. \]

Rearranging,

\[ (I - T A)\mathbf{x}[k+1] \approx \mathbf{x}[k] + T B \mathbf{u}[k], \]

and assuming \( I - T A \) is invertible,

\[ \mathbf{x}[k+1] \approx (I - T A)^{-1} \mathbf{x}[k] + (I - T A)^{-1} T B \mathbf{u}[k]. \]

Thus

\[ \Phi_{\text{BE}} = (I - T A)^{-1}, \quad \Gamma_{\text{BE}} = (I - T A)^{-1} T B. \]

Backward Euler is also first-order accurate, but has better stability properties for stiff systems (e.g. highly damped robot joints).

5.3 Bilinear (Tustin) approximation (preview)

Another common approach uses the bilinear transform in the transfer-function domain, replacing \( s \) by \( \tfrac{2}{T} \tfrac{z-1}{z+1} \). This maps the entire left-half \( s \)-plane to the interior of the unit circle in the \( z \)-plane and gives second-order accuracy in \( T \). Detailed frequency-domain analysis will be left for later digital control theory; for now, remember that Tustin is a more accurate but also more complex alternative to Euler methods.

flowchart TD
  Ct["Continuous model: 'x_dot = A x + B u'"] --> Ts["Choose sampling period T"]
  Ts --> Hold["Assume input hold type: 'zero order'"]
  Hold --> PhiStep["Compute Phi = exp(A*T) (or approximate)"]
  PhiStep --> GammaStep["Compute Gamma = integral_0^T exp(A*tau) B dtau (or approximation)"]
  GammaStep --> Disc["Form discrete model: 'x[k+1] = Phi x[k] + Gamma u[k]'"]
  Disc --> Impl["Implement in digital controller loop"]
        

We now implement discretization for the linearized single-joint robot model from Section 4.2 in several languages. We assume numerical values \( J = 0.01 \), \( b = 0.1 \), sampling period \( T = 0.001 \) s (1 kHz), and torque input \( u \).

\[ A = \begin{bmatrix} 0 & 1 \\ 0 & -10 \end{bmatrix}, \quad B = \begin{bmatrix} 0 \\ 100 \end{bmatrix}. \]

We show both exact methods (using high-level libraries) and simple Euler approximations suitable for embedded code.

6.1 Python (NumPy, SciPy, python-control, ROS)

In Python, it is natural to use numpy for linear algebra, scipy.linalg.expm for matrix exponentials, and the python-control toolbox for system objects. In practical robotics, such a controller would run inside a ROS node (rclpy) at a fixed rate, using the discrete model only conceptually (the actual implementation often uses continuous-time terms evaluated in discrete time).

import numpy as np
from scipy.linalg import expm
from control import ss, c2d  # python-control package

# Continuous-time joint model
J = 0.01
b = 0.1
A = np.array([[0.0,   1.0],
              [0.0, -b/J]])     # [[0, 1], [0, -10]]
B = np.array([[0.0],
              [1.0/J]])         # [[0], [100]]
C = np.eye(2)
D = np.zeros((2, 1))

Ts = 1e-3  # 1 kHz sampling

# Exact discrete-time matrices via matrix exponential
Phi = expm(A * Ts)
# Use Gamma = A^{-1}(Phi - I) B when A is nonsingular
Gamma = np.linalg.solve(A, (Phi - np.eye(2))) @ B

print("Phi =\n", Phi)
print("Gamma =\n", Gamma)

# Using python-control for a ZOH discretization
sys_c = ss(A, B, C, D)
sys_d = c2d(sys_c, Ts, method='zoh')

print("Phi (control) =\n", sys_d.A)
print("Gamma (control) =\n", sys_d.B)

# Discrete-time step in a control loop (conceptual)
def joint_step(x_k, u_k):
    """
    x_k: 2x1 numpy array [q, dq]
    u_k: scalar torque at step k
    """
    return Phi @ x_k + Gamma * u_k

# In a ROS controller, joint_step would be called at 1 kHz
# using the current torque command u_k and state estimate x_k.
      

6.2 C++ (Eigen, ros_control)

In C++, the Eigen library is standard for linear algebra. Embedded robot controllers inside ros_control often use simple Euler updates, especially when \( T \) is small and dynamics are moderately fast. Below we use forward Euler:

#include <Eigen/Dense>

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

struct JointDiscModel {
    Matrix2d Phi;
    Vector2d Gamma;

    JointDiscModel(double J, double b, double Ts) {
        Matrix2d A;
        Vector2d B;
        A << 0.0, 1.0,
               0.0, -b / J;
        B << 0.0,
               1.0 / J;

        // Forward Euler discretization: Phi = I + Ts*A, Gamma = Ts*B
        Matrix2d I = Matrix2d::Identity();
        Phi = I + Ts * A;
        Gamma = Ts * B;
    }

    void step(Vector2d &x, double u) const {
        // x = [q, dq]^T, u = torque
        x = Phi * x + Gamma * u;
    }
};

int main() {
    double J = 0.01;
    double b = 0.1;
    double Ts = 1e-3;

    JointDiscModel model(J, b, Ts);
    Vector2d x;
    x << 0.0, 0.0;  // start at rest
    double u = 1.0;   // constant torque

    for (int k = 0; k < 1000; ++k) { // simulate 1 second
        model.step(x, u);
        // send x(0) as simulated joint angle, etc.
    }
}
      

In a ros_control position or effort controller, the step method would be called in the real-time control loop, using measured states instead of simulated ones.

6.3 Java (EJML, ROSJava)

Java is less common for low-level robot control, but libraries such as EJML provide efficient linear algebra, and ROSJava enables integration with ROS-based systems. Here we implement the same forward Euler discretization with simple arrays (for clarity):

public class JointDiscModel {
    private final double[][] Phi;  // 2x2
    private final double[] Gamma;  // 2x1

    public JointDiscModel(double J, double b, double Ts) {
        double[][] A = {
            {0.0,        1.0},
            {0.0, -b / J}
        };
        double[] B = {0.0, 1.0 / J};

        Phi = new double[2][2];
        double[][] I = {
            {1.0, 0.0},
            {0.0, 1.0}
        };

        // Phi = I + Ts*A
        for (int i = 0; i < 2; ++i) {
            for (int j = 0; j < 2; ++j) {
                Phi[i][j] = I[i][j] + Ts * A[i][j];
            }
        }

        Gamma = new double[2];
        for (int i = 0; i < 2; ++i) {
            Gamma[i] = Ts * B[i];
        }
    }

    public void step(double[] x, double u) {
        double x0 = Phi[0][0] * x[0] + Phi[0][1] * x[1] + Gamma[0] * u;
        double x1 = Phi[1][0] * x[0] + Phi[1][1] * x[1] + Gamma[1] * u;
        x[0] = x0;
        x[1] = x1;
    }

    public static void main(String[] args) {
        double J = 0.01, b = 0.1, Ts = 1e-3;
        JointDiscModel model = new JointDiscModel(J, b, Ts);
        double[] x = {0.0, 0.0};
        double u = 1.0;

        for (int k = 0; k < 1000; ++k) {
            model.step(x, u);
            // publish x to ROSJava topics, log, etc.
        }
    }
}
      

6.4 MATLAB / Simulink (Control System Toolbox, Robotics Toolbox)

MATLAB has built-in support for continuous and discrete linear systems via ss and c2d. In robotics workflows, this often combines with the Robotics System Toolbox and Simulink models.

J = 0.01;
b = 0.1;
A = [0 1;
     0 -b/J];
B = [0;
     1/J];
C = eye(2);
D = zeros(2,1);

Ts = 1e-3;

sysc = ss(A, B, C, D);
% Exact ZOH discretization
sysd = c2d(sysc, Ts, 'zoh');

Phi   = sysd.A;
Gamma = sysd.B;

disp('Phi =');   disp(Phi);
disp('Gamma ='); disp(Gamma);

% Simulink implementation:
%  - Use a "State-Space" block with A=Phi, B=Gamma, C=eye(2), D=zeros(2,1).
%  - Drive the block with discrete torque command u[k] from a "Zero-Order Hold".
%  - Set the Simulink solver to "discrete (no continuous states)" with step Ts.
      

6.5 Wolfram Mathematica (StateSpaceModel, ToDiscreteTimeModel)

Mathematica offers symbolic and numeric tools for state-space models, including exact discretization under ZOH.

J = 0.01;
b = 0.1;

A = { {0, 1},
     {0, -b/J} };
B = { {0},
     {1/J} };
C = IdentityMatrix[2];
D = { {0}, {0} };

sysc = StateSpaceModel[{A, B, C, D}];

Ts = 0.001;

(* Exact zero-order-hold discretization *)
sysd = ToDiscreteTimeModel[sysc, Ts, Method -> "ZeroOrderHold"];

Phi   = StateSpaceModelData[sysd, "A"];
Gamma = StateSpaceModelData[sysd, "B"];

Phi // MatrixForm
Gamma // MatrixForm

(* One-step update as a pure function *)
jointStep[x_, u_] := Phi . x + Gamma . {u};
      

Mathematica can also symbolically compute \( e^{AT} \) and \( \Gamma \) for small systems, which is useful to derive closed-form formulas before implementing approximate versions on an embedded target.

Problem 1 (Derivation of the ZOH discrete-time model). Let \( \dot{\mathbf{x}}(t) = A \mathbf{x}(t) + B \mathbf{u}(t) \), with \( \mathbf{u}(t) \) generated by a zero-order hold from \( \mathbf{u}[k] \). Derive the exact discrete-time model \( \mathbf{x}[k+1] = \Phi \mathbf{x}[k] + \Gamma \mathbf{u}[k] \), and show that \( \Phi = e^{AT} \) and \( \Gamma = \int_{0}^{T} e^{A\tau} B\, d\tau \).

Solution.

1. The general solution for \( \dot{\mathbf{x}} = A \mathbf{x} + B \mathbf{u} \) is

   \[ \mathbf{x}(t) = e^{A(t - t_0)} \mathbf{x}(t_0) + \int_{t_0}^{t} e^{A(t - \sigma)} B \mathbf{u}(\sigma)\, d\sigma. \]

2. Choose \( t_0 = kT \), \( t = (k+1)T \), so

   \[ \mathbf{x}((k+1)T) = e^{AT} \mathbf{x}(kT) + \int_{kT}^{(k+1)T} e^{A((k+1)T - \sigma)} B \mathbf{u}(\sigma)\, d\sigma. \]

3. Under ZOH, \( \mathbf{u}(\sigma) = \mathbf{u}[k] \) for \( \sigma \in [kT,(k+1)T) \), so

   \[ \mathbf{x}((k+1)T) = e^{AT} \mathbf{x}(kT) + \left( \int_{kT}^{(k+1)T} e^{A((k+1)T - \sigma)} B \, d\sigma \right) \mathbf{u}[k]. \]

4. Substituting \( \sigma = kT + \tau \) gives

   \[ \int_{kT}^{(k+1)T} e^{A((k+1)T - \sigma)} B \, d\sigma = \int_{0}^{T} e^{A(T - \tau)} B \, d\tau. \]

   Using \( e^{A(T - \tau)} = e^{AT} e^{-A\tau} \), we can equivalently write

   \[ \Gamma = \int_{0}^{T} e^{A\tau} B \, d\tau, \]

   after changing the dummy integration variable.
5. Define \( \mathbf{x}[k] = \mathbf{x}(kT) \). Then

   \[ \mathbf{x}[k+1] = \Phi \mathbf{x}[k] + \Gamma \mathbf{u}[k], \quad \Phi = e^{AT}, \; \Gamma = \int_{0}^{T} e^{A\tau} B \, d\tau. \]

Problem 2 (Scalar ZOH discretization). Consider \( \dot{x}(t) = a x(t) + b u(t) \) with ZOH input \( u(t) = u[k] \) on \( [kT,(k+1)T) \). Derive the discrete-time equation \( x[k+1] = \Phi x[k] + \Gamma u[k] \) and identify \( \Phi \) and \( \Gamma \).

Solution. The scalar solution from time \( kT \) to \( (k+1)T \) is

\[ x((k+1)T) = e^{aT} x(kT) + \int_{kT}^{(k+1)T} e^{a((k+1)T - \sigma)} b u[k] \, d\sigma. \]

With \( u[k] \) constant, factor it out and integrate:

\[ \int_{kT}^{(k+1)T} e^{a((k+1)T - \sigma)} \, d\sigma = \int_{0}^{T} e^{a\tau} \, d\tau = \frac{1}{a}(e^{aT} - 1), \quad a \neq 0. \]

Thus

\[ x[k+1] = e^{aT} x[k] + \frac{b}{a}(e^{aT} - 1) u[k], \]

so \( \Phi = e^{aT} \), \( \Gamma = \tfrac{b}{a}(e^{aT} - 1) \). When \( a = 0 \), the limit gives \( \Gamma = bT \) and \( \Phi = 1 \).

Problem 3 (Eigenvalue mapping under sampling). Let \( A \) be diagonalizable with eigenvalues \( \lambda_i \). Show that the eigenvalues of \( \Phi = e^{AT} \) are \( z_i = e^{\lambda_i T} \).

Solution. If \( A = V \Lambda V^{-1} \) with \( \Lambda = \operatorname{diag}(\lambda_1,\dots,\lambda_n) \), then

\[ \Phi = e^{AT} = V e^{\Lambda T} V^{-1}, \]

where \( e^{\Lambda T} = \operatorname{diag}(e^{\lambda_1 T},\dots, e^{\lambda_n T}) \). Therefore, the eigenvalues of \( \Phi \) are \( z_i = e^{\lambda_i T} \). This mapping from the continuous-time spectrum to the discrete-time spectrum will be fundamental when we discuss discrete-time stability in the next lesson.

Problem 4 (Error of forward Euler discretization). For the scalar system in Problem 2, compare the exact discrete-time state transition \( \Phi = e^{aT} \) with the forward Euler approximation \( \Phi_{\text{FE}} = 1 + aT \). Derive the first non-zero term of the error expansion in powers of \( T \).

Solution. Expand the exponential:

\[ e^{aT} = 1 + aT + \frac{(aT)^2}{2!} + \frac{(aT)^3}{3!} + \dots. \]

The Euler approximation keeps only \( 1 + aT \). Thus

\[ e^{aT} - (1 + aT) = \frac{(aT)^2}{2!} + \mathcal{O}(T^3). \]

The leading error term scales like \( \tfrac{a^2}{2} T^2 \), confirming that forward Euler has local truncation error of order \( \mathcal{O}(T^2) \).

Problem 5 (Augmented matrix method for \( \Phi \) and \( \Gamma \)). Show that \( \Phi \) and \( \Gamma \) can be obtained simultaneously as blocks of a single matrix exponential:

\[ e^{\begin{bmatrix} A & B \\ 0 & 0 \end{bmatrix} T} = \begin{bmatrix} \Phi & \Gamma \\ 0 & I \end{bmatrix}. \]

Solution (sketch). Consider the augmented system

\[ \frac{d}{dt} \begin{bmatrix} \mathbf{x}(t) \\ \mathbf{u}(t) \end{bmatrix} = \begin{bmatrix} A & B \\ 0 & 0 \end{bmatrix} \begin{bmatrix} \mathbf{x}(t) \\ \mathbf{u}(t) \end{bmatrix}. \]

Under ZOH, \( \mathbf{u}(t) \) is constant on \( [kT,(k+1)T) \), so \( \mathbf{u}(t) = \mathbf{u}[k] \) and \( \dot{\mathbf{u}}(t) = 0 \). The solution from \( kT \) to \( (k+1)T \) is

\[ \begin{bmatrix} \mathbf{x}((k+1)T) \\ \mathbf{u}((k+1)T) \end{bmatrix} = e^{\begin{bmatrix} A & B \\ 0 & 0 \end{bmatrix} T} \begin{bmatrix} \mathbf{x}(kT) \\ \mathbf{u}(kT) \end{bmatrix}. \]

Because \( \mathbf{u}((k+1)T) = \mathbf{u}(kT) = \mathbf{u}[k] \), the lower-right block must be \( I \), and the upper blocks are precisely \( \Phi \) and \( \Gamma \), matching the definition from Section 4. This identity is widely used in numerical computation of ZOH-equivalent discrete models.

In this lesson we built the core bridge between continuous-time robot dynamics and their discrete-time implementations. We formalized sampling, defined the zero-order hold model for actuator interfaces, and derived the exact sample-data discrete-time model \( \mathbf{x}[k+1] = \Phi \mathbf{x}[k] + \Gamma \mathbf{u}[k] \) with \( \Phi = e^{AT} \) and \( \Gamma = \int_{0}^{T} e^{A\tau} B \, d\tau \). We also introduced practical approximate methods such as forward and backward Euler and illustrated how to compute discrete models in Python, C++, Java, MATLAB/Simulink, and Wolfram Mathematica, in the context of a simple joint model. In the next lesson we will analyze discrete-time stability and performance using these discrete models.