Chapter 28: Computer-Aided Analysis and Design for Linear Control

Lesson 5: Reproducible Design Workflows and Reporting

This lesson develops a rigorous view of reproducible computer-aided workflows for linear controller design. We formalize a design experiment as a mapping from plant model and specifications to controller parameters and performance metrics, then show how to implement such workflows in Python, C++, Java, MATLAB/Simulink, and Wolfram Mathematica. The emphasis is on deterministic scripts, parameter tracking, and report-quality figures and tables that can be regenerated exactly.

1. Reproducibility in Linear Control Design

In a classical linear control setting, a design task starts from a plant transfer function \( G(s) \), performance specifications, and a chosen controller structure \( C(s,\boldsymbol{\theta}) \) with parameter vector \( \boldsymbol{\theta}\in\mathbb{R}^p \). A reproducible workflow ensures that another engineer, possibly years later, can regenerate:

  • the controller parameters \( \boldsymbol{\theta} \),
  • the closed-loop transfer functions and stability margins,
  • the time- and frequency-response plots and numerical tables,
  • all from a precisely specified set of scripts and configuration files.

For a single-input–single-output (SISO) continuous-time plant, we typically specify

\[ G(s) = \frac{b_m s^m + \cdots + b_1 s + b_0}{s^n + a_{n-1} s^{n-1} + \cdots + a_1 s + a_0}, \quad a_0,\dots,a_{n-1}, b_0,\dots,b_m \in \mathbb{R}. \]

A controller family (e.g. PID, lead, lag) is parameterized as

\[ C(s,\boldsymbol{\theta}) = \frac{\beta_q(\boldsymbol{\theta}) s^q + \cdots + \beta_0(\boldsymbol{\theta})}{s^r + \gamma_{r-1}(\boldsymbol{\theta}) s^{r-1} + \cdots + \gamma_0(\boldsymbol{\theta})}, \quad \boldsymbol{\theta} \in \Theta \subset \mathbb{R}^p. \]

The open-loop transfer function is \( L(s,\boldsymbol{\theta}) = G(s)\,C(s,\boldsymbol{\theta}) \). With unity feedback, the closed-loop sensitivity and complementary sensitivity are

\[ S(s,\boldsymbol{\theta}) = \frac{1}{1 + L(s,\boldsymbol{\theta})}, \qquad T(s,\boldsymbol{\theta}) = \frac{L(s,\boldsymbol{\theta})}{1 + L(s,\boldsymbol{\theta})}. \]

A reproducible workflow must therefore specify all ingredients of the mapping

\[ \mathcal{F} : \bigl(G(s), \text{specs}, \text{code}, \text{config}\bigr) \mapsto \bigl(\boldsymbol{\theta}, S(\cdot,\boldsymbol{\theta}), T(\cdot,\boldsymbol{\theta}), \text{metrics}, \text{plots}, \text{tables}\bigr) \]

where “config” includes numerical tolerances, simulation horizons, sampling grids, and software versions. If \( \mathcal{F} \) is deterministic and fully documented, the design is reproducible.

2. Formal Design Records and Performance Functionals

We can formalize a control design experiment as a tuple \( \mathcal{R} = (G, \boldsymbol{\theta}, \mathcal{S}, \mathcal{M}) \) where:

  • \( G \) is the plant transfer function (with dimensions and units recorded),
  • \( \boldsymbol{\theta} \) is the controller parameter vector,
  • \( \mathcal{S} \) is a set of simulation settings (time grid, step size, numerical method, tolerances),
  • \( \mathcal{M} \) is a set of performance functionals.

Typical performance functionals are derived from the step response \( y(t;\boldsymbol{\theta}) \) of the closed loop with unit-step input:

\[ \begin{aligned} M_p(\boldsymbol{\theta}) &= \max_{t\ge 0} \bigl(y(t;\boldsymbol{\theta}) - 1\bigr)_+, \\ t_s(\varepsilon;\boldsymbol{\theta}) &= \inf\bigl\{t \ge 0 : |y(\tau;\boldsymbol{\theta}) - 1| \le \varepsilon \ \forall \tau \ge t\bigr\}, \\ e_{\text{ss}}(\boldsymbol{\theta}) &= \lim_{t \to \infty} \bigl(1 - y(t;\boldsymbol{\theta})\bigr), \end{aligned} \]

where \( (x)_+ = \max\{x,0\} \) and \( \varepsilon \) is a small tolerance (e.g. \(0.02\) for a \(2\%\) settling band). In frequency domain, one also records

\[ \begin{aligned} M_T(\boldsymbol{\theta}) &= \max_{\omega \ge 0} |T(j\omega,\boldsymbol{\theta})|, \\ M_S(\boldsymbol{\theta}) &= \max_{\omega \ge 0} |S(j\omega,\boldsymbol{\theta})|, \\ \omega_{\text{BW}}(\boldsymbol{\theta}) &= \sup\{\omega \ge 0 : |T(j\omega,\boldsymbol{\theta})| \ge 1/\sqrt{2}\}. \end{aligned} \]

A design record is then a labeled collection \( \{\mathcal{R}_k\}_{k=1}^N \) produced by one script, where the index \( k \) may correspond to different controllers, different plant parameter samples (robustness study), or different design choices.

Reproducibility demands that for each record \( \mathcal{R}_k \), there exists a script such that rerunning it regenerates the same \( \boldsymbol{\theta}_k \) and the same numerical values of \( M_p, t_s, e_{\text{ss}}, M_T, M_S, \omega_{\text{BW}} \) up to a clearly specified numerical precision.

3. Script-Driven Workflow and Versioning

GUI-based design (dragging poles on a root locus, interactively tuning sliders) is convenient but often not reproducible: it is difficult to reconstruct exactly which sequence of interactions produced the final design. Instead, we use scripted workflows:

  1. Define plant and units in code (no manual entry into GUI tables).
  2. Encode controller structure and design algorithms as deterministic functions.
  3. Store controller parameters and simulation settings in a configuration file.
  4. Generate all plots and tables via scripts, not manual exporting.
  5. Use version control (e.g. git) for scripts, configuration, and reports.

A typical linear control design workflow can be sketched as:

flowchart TD
  A["Specs (overshoot, settling time, ess)"] --> B["Code: model G(s) and controller family C(s, theta)"]
  B --> C["Config: numeric options, time grid, plant params"]
  C --> D["Run design script (all languages)"]
  D --> E["Computed theta, stability margins, time/frequency responses"]
  E --> F["Generate figures and tables for report"]
  F --> G["Commit code, config, and report to repository"]
        

The role of robotics libraries (for example, Robotics Toolbox for Python, ROS control libraries in C++, Java robotics frameworks such as WPILib, and MATLAB Robotics System Toolbox) is to supply plant models for robot joints, links, and actuators. In a reproducible workflow, these libraries are only accessed through scripted APIs with version numbers recorded, so that the underlying kinematic or dynamic maps remain traceable.

4. Python Workflow for Reproducible Linear Controller Design

Consider the plant \( G(s) = \dfrac{1}{(s+1)^2} \), representative of a lightly damped rigid-body mode with actuator dynamics. We design a PI controller \( C(s,\boldsymbol{\theta}) = K_p + \dfrac{K_i}{s} \) and record:

  • Controller gains \( (K_p, K_i) \),
  • Closed-loop step response and derived metrics,
  • Configuration (final time, time step, numerical library versions).

The following Python script (using the control library and standard scientific stack) implements a reproducible experiment. The same structure can be used when the plant comes from a robotics library such as roboticstoolbox, for example a single revolute joint model.


import json
import platform
from datetime import datetime

import numpy as np
import control as ct   # python-control library
import matplotlib.pyplot as plt

# ------------------------------
# 1. Configuration (single source of truth)
# ------------------------------
config = {
    "experiment_name": "pi_design_second_order",
    "plant": {
        "num": [1.0],
        "den": [1.0, 2.0, 1.0]  # (s + 1)^2
    },
    "controller": {
        "structure": "PI",
        "Kp": 4.0,
        "Ki": 3.0
    },
    "simulation": {
        "t_final": 10.0,
        "n_points": 2000
    },
    "tolerances": {
        "settling_band": 0.02
    },
    "environment": {
        "python_version": platform.python_version(),
        "control_version": ct.__version__,
        "numpy_version": np.__version__
    }
}

# ------------------------------
# 2. Construct plant and controller from config
# ------------------------------
G = ct.TransferFunction(config["plant"]["num"], config["plant"]["den"])
Kp = config["controller"]["Kp"]
Ki = config["controller"]["Ki"]
C = ct.TransferFunction([Kp, Ki], [1.0, 0.0])  # Kp + Ki / s

# Closed-loop with unity feedback
L = C * G
T = ct.feedback(L, 1.0)

# ------------------------------
# 3. Simulate step response and compute metrics
# ------------------------------
t = np.linspace(0.0, config["simulation"]["t_final"], config["simulation"]["n_points"])
t_out, y_out = ct.step_response(T, T=t)

# Overshoot and steady-state error
y_final = y_out[-1]
M_p = float(np.max(y_out) - 1.0)
e_ss = float(1.0 - y_final)

# Settling time: last time when |y - 1| is greater than band
band = config["tolerances"]["settling_band"]
mask_outside = np.abs(y_out - 1.0) > band
if np.any(mask_outside):
    last_outside_index = np.where(mask_outside)[0][-1]
    t_s = float(t_out[last_outside_index])
else:
    t_s = 0.0

metrics = {
    "overshoot": M_p,
    "steady_state_error": e_ss,
    "settling_time": t_s
}

# ------------------------------
# 4. Persist results (JSON) for later reuse
# ------------------------------
result_record = {
    "timestamp": datetime.utcnow().isoformat() + "Z",
    "config": config,
    "metrics": metrics
}

with open("pi_design_second_order_result.json", "w", encoding="utf-8") as f:
    json.dump(result_record, f, indent=2)

# ------------------------------
# 5. Generate reproducible plots
# ------------------------------
plt.figure()
plt.plot(t_out, y_out)
plt.axhline(1.0, linestyle="--")
plt.xlabel("Time (s)")
plt.ylabel("Output y(t)")
plt.title("Closed-loop step response (PI-controlled second-order plant)")
plt.grid(True)
plt.tight_layout()
plt.savefig("pi_design_second_order_step_response.png", dpi=150)
plt.close()

# Note: In a robotics context, G could represent the linearized transfer from
# joint torque to joint angle, obtained for instance from roboticstoolbox or
# another robotics library. The same script structure still guarantees that
# given the same config, metrics and plots are reproducible.
      

The combination of the JSON file, the PNG figure, and the Python script constitutes a complete, machine-readable design record. For a robotic joint, only the plant definition section would change.

5. C++ Workflow with Eigen and Robotics Context

In C++, reproducible workflows commonly rely on a small layer of utility code and linear algebra libraries such as Eigen. In robotics settings, these scripts can be integrated with ROS or other middleware while still maintaining determinism: controller gains are loaded from configuration files (YAML, JSON) and logged back after tuning.

The snippet below implements a simple numerical integration of the closed-loop system for the same plant \( G(s) = \dfrac{1}{(s+1)^2} \) with a PI controller. The plant is represented in state-space form (obtained from the transfer function by standard realization techniques) but the underlying mathematics still describe the same closed-loop dynamics.


#include <iostream>
#include <fstream>
#include <vector>
#include <string>
#include <cmath>

// Simple struct for configuration (could be loaded from JSON/YAML)
struct Config {
    double Kp;
    double Ki;
    double dt;
    double t_final;
    double settling_band;
};

// Closed-loop simulation for G(s) = 1 / (s + 1)^2 with PI controller
// State-space realization (minimal) of plant with input u and output y:
// x_dot = A x + B u, y = C x
// For G(s) = 1 / (s + 1)^2, one possible realization is:
// A = [ -1  1; 0  -1 ], B = [0; 1], C = [1  0]
struct PlantPI {
    Config cfg;
    double x1, x2;  // plant states
    double xI;      // integral state for PI

    explicit PlantPI(const Config& c)
        : cfg(c), x1(0.0), x2(0.0), xI(0.0) {}

    void reset() {
        x1 = 0.0; x2 = 0.0; xI = 0.0;
    }

    // One integration step using forward Euler for simplicity
    double step(double r) {
        double y = x1;
        double e = r - y;
        xI += e * cfg.dt;

        double u = cfg.Kp * e + cfg.Ki * xI;

        double x1_dot = -x1 + x2;
        double x2_dot = -x2 + u;

        x1 += x1_dot * cfg.dt;
        x2 += x2_dot * cfg.dt;
        return y;
    }
};

int main() {
    Config cfg;
    cfg.Kp = 4.0;
    cfg.Ki = 3.0;
    cfg.dt = 0.001;
    cfg.t_final = 10.0;
    cfg.settling_band = 0.02;

    PlantPI system(cfg);
    const int n_steps = static_cast<int>(cfg.t_final / cfg.dt);
    std::vector<double> t_vec;
    std::vector<double> y_vec;
    t_vec.reserve(n_steps + 1);
    y_vec.reserve(n_steps + 1);

    double t = 0.0;
    for (int k = 0; k <= n_steps; ++k) {
        double y = system.step(1.0);  // unit step
        t_vec.push_back(t);
        y_vec.push_back(y);
        t += cfg.dt;
    }

    // Compute simple metrics
    double y_final = y_vec.back();
    double Mp = 0.0;
    for (double y : y_vec) {
        if (y - 1.0 > Mp) Mp = y - 1.0;
    }
    double ess = 1.0 - y_final;

    double t_s = 0.0;
    for (int k = n_steps; k >= 0; --k) {
        if (std::fabs(y_vec[static_cast<std::size_t>(k)] - 1.0) > cfg.settling_band) {
            t_s = t_vec[static_cast<std::size_t>(k)];
            break;
        }
    }

    // Persist results (CSV for plots, text for metrics)
    std::ofstream data_file("cpp_pi_closed_loop.csv");
    data_file << "t,y\n";
    for (std::size_t k = 0; k < t_vec.size(); ++k) {
        data_file << t_vec[k] << "," << y_vec[k] << "\n";
    }

    std::ofstream metrics_file("cpp_pi_metrics.txt");
    metrics_file << "Kp=" << cfg.Kp << "\n";
    metrics_file << "Ki=" << cfg.Ki << "\n";
    metrics_file << "Mp=" << Mp << "\n";
    metrics_file << "ess=" << ess << "\n";
    metrics_file << "ts=" << t_s << "\n";

    std::cout << "Results written to cpp_pi_closed_loop.csv and cpp_pi_metrics.txt\n";

    // In a robotics stack, the same simulation could be wrapped in a ROS node,
    // with cfg loaded from and written back to the ROS parameter server to
    // guarantee that all experiments are driven by explicit parameters.
    return 0;
}
      

Here, reproducibility is achieved by:

  • centralizing all tunable parameters in Config,
  • logging outputs and metrics to deterministic text files,
  • placing the main program under version control together with configuration.

6. Java Workflow and Structured Reporting

Java is less common in traditional control engineering but widely used in robotics, especially in educational and industrial platforms (e.g. FRC robots with WPILib). The following skeleton demonstrates how a reproducible controller design experiment can be encoded in Java using simple numerical integration and collections.


import java.io.FileWriter;
import java.io.IOException;
import java.util.ArrayList;
import java.util.List;

class Config {
    public double Kp = 4.0;
    public double Ki = 3.0;
    public double dt = 0.001;
    public double tFinal = 10.0;
    public double settlingBand = 0.02;
}

// Same plant as before: G(s) = 1 / (s + 1)^2, PI controller
class PlantPI {
    private Config cfg;
    private double x1 = 0.0;
    private double x2 = 0.0;
    private double xI = 0.0;

    public PlantPI(Config cfg) {
        this.cfg = cfg;
    }

    public void reset() {
        x1 = 0.0;
        x2 = 0.0;
        xI = 0.0;
    }

    public double step(double r) {
        double y = x1;
        double e = r - y;
        xI += e * cfg.dt;
        double u = cfg.Kp * e + cfg.Ki * xI;

        double x1dot = -x1 + x2;
        double x2dot = -x2 + u;

        x1 += x1dot * cfg.dt;
        x2 += x2dot * cfg.dt;
        return y;
    }
}

public class PiExperiment {
    public static void main(String[] args) throws IOException {
        Config cfg = new Config();
        PlantPI system = new PlantPI(cfg);

        int nSteps = (int) Math.round(cfg.tFinal / cfg.dt);
        List<Double> tList = new ArrayList<>(nSteps + 1);
        List<Double> yList = new ArrayList<>(nSteps + 1);

        double t = 0.0;
        for (int k = 0; k <= nSteps; ++k) {
            double y = system.step(1.0);
            tList.add(t);
            yList.add(y);
            t += cfg.dt;
        }

        // Metrics
        double yFinal = yList.get(yList.size() - 1);
        double Mp = 0.0;
        for (double y : yList) {
            double val = y - 1.0;
            if (val > Mp) Mp = val;
        }
        double ess = 1.0 - yFinal;

        double ts = 0.0;
        for (int k = nSteps; k >= 0; --k) {
            double y = yList.get(k);
            if (Math.abs(y - 1.0) > cfg.settlingBand) {
                ts = tList.get(k);
                break;
            }
        }

        // CSV export
        try (FileWriter fw = new FileWriter("java_pi_closed_loop.csv")) {
            fw.write("t,y\n");
            for (int k = 0; k < tList.size(); ++k) {
                fw.write(tList.get(k) + "," + yList.get(k) + "\n");
            }
        }

        // Metrics export
        try (FileWriter fw = new FileWriter("java_pi_metrics.txt")) {
            fw.write("Kp=" + cfg.Kp + "\n");
            fw.write("Ki=" + cfg.Ki + "\n");
            fw.write("Mp=" + Mp + "\n");
            fw.write("ess=" + ess + "\n");
            fw.write("ts=" + ts + "\n");
        }

        // In robot control libraries such as WPILib, the same numerical
        // experiment can be run in a test or simulation mode, reading cfg
        // from a configuration file and writing back metrics for documentation.
    }
}
      

Once again, the reproducible unit is the combination of: configuration fields Config, the integrated dynamics in PlantPI, and the exported CSV and text files. In a robotics setting, the step input could correspond, for instance, to a desired joint angle command, and the same logging structure would be used.

7. MATLAB/Simulink Workflow and Robotics Toolboxes

MATLAB and Simulink are central tools in industrial linear control design. A reproducible workflow uses scripts (M-files) and model files (SLX) under version control. Plant models may come from tf, from linearization of Simscape Multibody models, or from Robotics System Toolbox models of manipulator dynamics.

The following script shows a fully scripted experiment for the same plant and PI controller, including generation of Bode and step response plots, and a table of metrics.


% config.m (reproducible configuration)
cfg.experiment_name = 'pi_design_second_order';
cfg.plant.num = 1;
cfg.plant.den = [1 2 1];  % (s + 1)^2
cfg.controller.Kp = 4.0;
cfg.controller.Ki = 3.0;
cfg.sim.t_final = 10;
cfg.sim.n_points = 2000;
cfg.tolerances.settling_band = 0.02;

% plant and controller
s = tf('s');
G = tf(cfg.plant.num, cfg.plant.den);
C = cfg.controller.Kp + cfg.controller.Ki / s;

L = C * G;
T = feedback(L, 1);  % unity feedback

% Bode and step plots (saved to files, no manual export)
fig1 = figure('Visible', 'off');
bode(L);
grid on;
title('Open-loop Bode plot');
saveas(fig1, 'matlab_pi_openloop_bode.png');

fig2 = figure('Visible', 'off');
[y, t] = step(T, linspace(0, cfg.sim.t_final, cfg.sim.n_points));
plot(t, y); hold on;
yline(1, '--');
xlabel('Time (s)');
ylabel('Output y(t)');
title('Closed-loop step response (PI)');
grid on;
saveas(fig2, 'matlab_pi_step_response.png');

close(fig1); close(fig2);

% Compute metrics
y_final = y(end);
Mp = max(y) - 1;
ess = 1 - y_final;

idx_last_outside = find(abs(y - 1) > cfg.tolerances.settling_band, 1, 'last');
if isempty(idx_last_outside)
    ts = 0;
else
    ts = t(idx_last_outside);
end

metrics = table(Mp, ess, ts, ...
    'VariableNames', {'Overshoot', 'SteadyStateError', 'SettlingTime'});

writetable(metrics, 'matlab_pi_metrics.csv');

% Example robotics extension:
% If G is derived from a linearized robot joint model via
% linmod('robot_joint_model') or linearize from Robotics System Toolbox,
% the exact same structure applies and metrics remain reproducible.
      

For Simulink-driven experiments, the controller and plant can be implemented as blocks in a model, e.g. pid_reproducible.slx, with parameters driven from cfg by calls to set_param, and the entire simulation controlled by a single script. Logging to MAT-files and CSV tables completes the reproducible record.

8. Wolfram Mathematica Workflow

Wolfram Mathematica provides symbolic and numerical tools for linear systems through transfer function and state-space objects. A reproducible workflow is naturally expressed as a notebook that runs top-to-bottom without manual intervention. The following code illustrates the same plant and PI controller in Mathematica, along with automatic computation of basic metrics.


(* Configuration *)
cfg = <|
  "Kp" -> 4.0,
  "Ki" -> 3.0,
  "tFinal" -> 10.0,
  "nPoints" -> 2000,
  "SettlingBand" -> 0.02
|>;

s = LaplaceTransformVariable["s"];

(* Plant and controller *)
G = TransferFunctionModel[1/(s + 1)^2, s];
C = TransferFunctionModel[cfg["Kp"] + cfg["Ki"]/s, s];

L = SeriesConnect[C, G];
T = FeedbackConnect[L, 1];

(* Time grid and step response *)
tGrid = Subdivide[0.0, cfg["tFinal"], cfg["nPoints"]];
y = OutputResponse[T, UnitStep[t], {t, 0, cfg["tFinal"]}];
yVals = y /. t -> # & /@ tGrid;

(* Metrics *)
yFinal = Last[yVals];
Mp = Max[yVals] - 1.0;
ess = 1.0 - yFinal;

idxLastOutside = Last@
  Flatten@Position[Abs[yVals - 1.0], _?(# > cfg["SettlingBand"] &)];
ts = If[idxLastOutside === Missing["NotAvailable"], 0.0, tGrid[[idxLastOutside]]];

metrics = <|
  "Overshoot" -> Mp,
  "SteadyStateError" -> ess,
  "SettlingTime" -> ts
|>;

Export["mma_pi_metrics.json", metrics, "JSON"];
Export["mma_pi_step_response.dat", Transpose[{tGrid, yVals}]];

(* A ListLinePlot of Transpose[{tGrid, yVals}] can be generated in the notebook
   and exported as a figure for reports. In robotics contexts, G can instead
   represent a transfer function obtained from a linearized manipulator model. *)
      

Because the configuration is stored in a single association cfg and all exports are performed programmatically, simply re-evaluating the notebook regenerates the same controller, plots, and metrics.

9. Automated Parameter Sweeps, Robustness Studies, and Reporting

Lesson 4 of this chapter introduced parameter sweeps and robustness studies. Reproducibility adds a further requirement: for each sampled parameter set (e.g. plant gains, controller gains), we must be able to regenerate the closed-loop performance metrics and plots. Mathematically, suppose we have a family of plants \( G(s,\boldsymbol{\rho}) \) indexed by parameter vector \( \boldsymbol{\rho} \), such as uncertain masses or stiffnesses. For a fixed controller structure \( C(s,\boldsymbol{\theta}) \), define

\[ J(\boldsymbol{\theta},\boldsymbol{\rho}) = \int_0^{\infty} w(t)\,\bigl(y(t;\boldsymbol{\theta},\boldsymbol{\rho})-1\bigr)^2 \, dt, \]

where \( w(t) \) is a weighting function (for example, \(w(t) = e^{-\alpha t}\) with \(\alpha > 0\)) emphasizing certain time ranges. A parameter sweep evaluates \( J(\boldsymbol{\theta},\boldsymbol{\rho}_k) \) for a grid \( \{\boldsymbol{\rho}_k\}_{k=1}^N \), producing a reproducible table of values.

A practical reporting pipeline for such studies can be sketched as:

flowchart TD
  A["Config grid for plant and controller parameters"] --> B["Batch scripts (Python/C++/Java/MATLAB/MMA)"]
  B --> C["For each sample: simulate, compute metrics, save row to table"]
  C --> D["Aggregate tables and generate plots (heatmaps, histograms)"]
  D --> E["Include generated figures/tables in technical report or thesis"]
        

If every entry in the reported tables is produced by these batch scripts and the scripts, configuration files, and software versions are archived, the robustness study is fully reproducible. In robotics projects, this is essential when comparing alternative linear controllers for robot joints, end-effector motion, or flexible structures.

10. Problems and Solutions

Problem 1 (Controller Matching and Deterministic Design Map): Consider the plant \( G(s) = \dfrac{1}{s(s+1)} \) and a PD controller \( C(s,\boldsymbol{\theta}) = K_p + K_d s \). For unity feedback, derive \( K_p \) and \( K_d \) such that the closed-loop characteristic polynomial equals \( s^2 + 2\zeta\omega_n s + \omega_n^2 \) for prescribed \( \zeta \in (0,1) \) and \( \omega_n > 0 \). Explain why a script implementing this relationship yields a reproducible design.

Solution: The open-loop transfer function is

\[ L(s,\boldsymbol{\theta}) = G(s)C(s,\boldsymbol{\theta}) = \frac{K_p + K_d s}{s(s+1)}. \]

With unity feedback, the closed-loop transfer function is \( T(s,\boldsymbol{\theta}) = \dfrac{L(s,\boldsymbol{\theta})} {1 + L(s,\boldsymbol{\theta})} \), and the denominator of \( T(s,\boldsymbol{\theta}) \) is

\[ D_{\text{cl}}(s) = s(s+1) + (K_p + K_d s) = s^2 + s + K_p + K_d s = s^2 + (1 + K_d)s + K_p. \]

Matching to the desired second-order polynomial \( s^2 + 2\zeta\omega_n s + \omega_n^2 \) yields

\[ \begin{aligned} 1 + K_d &= 2\zeta\omega_n, \\ K_p &= \omega_n^2. \end{aligned} \]

Hence

\[ K_d = 2\zeta\omega_n - 1, \qquad K_p = \omega_n^2. \]

A script that takes \( (\zeta,\omega_n) \) as inputs and computes \( (K_p,K_d) \) via these formulas implements a deterministic mapping \( (\zeta,\omega_n) \mapsto (K_p,K_d) \). If the script, its inputs, and the floating-point precision are fixed and recorded, any future engineer can regenerate exactly the same gains, making the design reproducible.

Problem 2 (Convergence of Discrete Performance Functional): Let \( y(t) \) denote the closed-loop step response of a stable linear system, and define the continuous-time quadratic performance functional

\[ J = \int_0^{\infty} \bigl(y(t) - 1\bigr)^2 \, dt. \]

Suppose we approximate \( J \) numerically on a uniform grid with step \( h > 0 \) by

\[ J_h = h \sum_{k=0}^{\infty} \bigl(y(kh) - 1\bigr)^2. \]

Assume \( y(t) \) is continuous and converges exponentially to 1, so that \( J < \infty \). Show that \( J_h \to J \) as \( h \to 0 \), and explain the implication for reproducible performance evaluation.

Solution: Because the closed-loop system is exponentially stable, there exist constants \( M \ge 0 \) and \( \alpha > 0 \) such that \( |y(t) - 1| \le M e^{-\alpha t} \) for all \( t \ge 0 \). Consequently, \( f(t) = (y(t) - 1)^2 \) is continuous, nonnegative, and satisfies \( f(t) \le M^2 e^{-2\alpha t} \), making the integral \( J \) absolutely convergent.

The Riemann sums \( J_h = h \sum_{k=0}^{\infty} f(kh) \) converge to the integral of \( f \) over \( [0,\infty) \) as \( h \to 0 \), because:

  • on any finite interval \( [0,T] \), the standard Riemann convergence theorem applies to the continuous function \( f \),
  • the tail integral over \( [T,\infty) \) and the tail sum for \( k h \ge T \) can both be made arbitrarily small using the exponential bound, uniformly in \( h \).

Therefore, \( J_h \to J \) as \( h \to 0 \). In a reproducible workflow, specifying the integration step \( h \) and documenting that it is sufficiently small ensures that reported values of \( J_h \) approximate the mathematical quantity \( J \) within a known error bound.

Problem 3 (Well-Defined Settling Time for Exponential Convergence): Let \( y(t) \) be the closed-loop step response of a stable linear system with exponential convergence: \( |y(t) - 1| \le M e^{-\alpha t} \) for some \( M,\alpha > 0 \). For a tolerance \( \varepsilon > 0 \), define the settling time

\[ t_s(\varepsilon) = \inf\bigl\{ t \ge 0 : |y(\tau) - 1| \le \varepsilon \ \forall \tau \ge t \bigr\}. \]

Prove that \( t_s(\varepsilon) \) is finite, and derive an explicit upper bound expressed in terms of \( M \), \( \alpha \), and \( \varepsilon \).

Solution: From the exponential bound \( |y(t) - 1| \le M e^{-\alpha t} \), we want to find \( t \) such that \( M e^{-\alpha t} \le \varepsilon \). Solving for \( t \) gives

\[ M e^{-\alpha t} \le \varepsilon \quad \Rightarrow \quad e^{-\alpha t} \le \frac{\varepsilon}{M} \quad \Rightarrow \quad -\alpha t \le \ln\left(\frac{\varepsilon}{M}\right) \quad \Rightarrow \quad t \ge \frac{1}{\alpha} \ln\left(\frac{M}{\varepsilon}\right). \]

Thus, any \( t \ge \dfrac{1}{\alpha} \ln\left(\dfrac{M}{\varepsilon}\right) \) guarantees \( |y(t) - 1| \le \varepsilon \). Because the bound is monotone in \( t \), it follows that \( t_s(\varepsilon) \) is finite and satisfies

\[ t_s(\varepsilon) \le \frac{1}{\alpha} \ln\left(\frac{M}{\varepsilon}\right). \]

In practice, this inequality can be used to select a simulation horizon long enough to guarantee that the numerical estimate of settling time is meaningful. Recording \( \alpha \) and \( M \) (e.g. via pole locations and initial condition bounds) strengthens the theoretical basis for reproducible reporting of \( t_s(\varepsilon) \).

Problem 4 (Equivalence of Scripted and Manual Root-Locus Gain): Suppose that for a given unity-feedback loop with loop transfer function \( L(s) = K L_0(s) \), a designer manually identifies a gain \( K^\star \) on the root locus such that the dominant closed-loop poles are located at \( s = s_{1,2} \). Show that if the same \( K^\star \) is computed in a script by solving the characteristic equation

\[ 1 + K L_0(s) = 0 \]

at \( s = s_{1,2} \), the resulting closed-loop transfer function is identical. Explain why the scripted approach is preferable from a reproducibility standpoint.

Solution: The closed-loop characteristic equation is \( 1 + K L_0(s) = 0 \). If a manual root-locus procedure picks \( K^\star \) such that \( s_{1,2} \) lies on the root locus, then by definition

\[ 1 + K^\star L_0(s_{1,2}) = 0. \]

A script that solves \( 1 + K L_0(s_{1,2}) = 0 \) for \( K \) yields exactly

\[ K = -\frac{1}{L_0(s_{1,2})} = K^\star, \]

assuming \( L_0(s_{1,2}) \neq 0 \). Thus, the characteristic polynomial, and hence the closed-loop transfer function, are identical between the manual and scripted designs. However, the manual procedure typically leaves no precise record of the graphical construction, whereas the script encodes the exact numerical values and the algebraic steps. Re-running the script with the same \( s_{1,2} \) regenerates \( K^\star \) exactly, satisfying reproducibility requirements.

11. Summary

In this lesson we formalized reproducible computer-aided workflows for linear control design. We modeled the design process as a deterministic mapping from plant, specifications, and configuration to controller parameters, closed-loop transfer functions, and performance metrics. We stressed the importance of scripting all steps (plant definition, controller synthesis, simulation, plotting, and table generation) and of recording configuration metadata such as simulation grids and software versions.

We implemented concrete workflows in Python, C++, Java, MATLAB/Simulink, and Wolfram Mathematica, and sketched how robotics libraries fit naturally into this structure by supplying plant models for robot joints and mechanisms. Finally, we connected these workflows to parameter sweeps and robustness experiments, and we analyzed mathematically the convergence of numerical performance functionals and the well-posedness of settling time. These tools and concepts enable students to build linear control designs that are not only effective but also scientifically reproducible and auditable.

12. References

  1. Claerbout, J.F., & Karrenbach, M. (1992). Electronic documents give reproducible research a new meaning. SEG Technical Program Expanded Abstracts, 27–30.
  2. Donoho, D.L. (2010). An invitation to reproducible computational research. Biostatistics, 11(3), 385–388.
  3. Åström, K.J., & Hägglund, T. (1984). Automatic tuning of simple regulators with specifications on phase and amplitude margins. Automatica, 20(5), 645–651.
  4. Doyle, J.C., Glover, K., Khargonekar, P.P., & Francis, B.A. (1989). State-space solutions to standard “H-infinity” and “H2” control problems. IEEE Transactions on Automatic Control, 34(8), 831–847.
  5. Skogestad, S., & Postlethwaite, I. (1997). Multivariable feedback control: Analysis and design (selected theoretical developments). International Journal of Control, 65(2), 277–314.
  6. Franklin, G.F., Powell, J.D., & Workman, M.L. (1990). Digital control of dynamic systems: Theoretical foundations of computer-aided design. International Journal of Control, 51(4), 821–847.
  7. Ljung, L. (1999). System identification: Theory for the user (theoretical aspects of model-based control design). IEEE Control Systems Magazine, 19(3), 32–41.
  8. Sung, C., & Lee, J.H. (1994). A unified approach to robust PID controller design. Automatica, 30(9), 1585–1589.
  9. Morari, M., & Zafiriou, E. (1989). Robust process control: Fundamental theoretical issues in computer-aided control design. Chemical Engineering Science, 44(8), 1787–1805.
  10. Trefethen, L.N., & Bau, D. (1997). Numerical linear algebra and its applications to control. SIAM Review, 39(3), 423–450.