Control Framework Architecture

Classical Control Theory for Analysis and Synthesis

Author

Gil Benezer

Published

January 16, 2026

Overview

The Control Framework provides classical control theory algorithms for analysis and synthesis of dynamical systems. It consists of 3 core modules organized into a clean 2-layer architecture: pure stateless functions and thin composition wrappers.

Framework Components

Pure Function Layer: classical_control_functions.py - Stateless control algorithms
Composition Layer: control_synthesis.py, system_analysis.py - System integration
Type Layer: control_classical.py - TypedDict result definitions

User Interaction Model

Users interact with the control framework through system properties:

# ✓ CORRECT: Access via system properties
system = Pendulum()
A, B = system.linearize(x_eq, u_eq)

# Control design
lqr = system.control.design_lqr(A, B, Q, R, system_type='continuous')
K = lqr['gain']

# System analysis
stability = system.analysis.stability(A, system_type='continuous')
controllability = system.analysis.controllability(A, B)

# ✗ INCORRECT: Direct function imports (not recommended)
from cdesym.control.classical_control_functions import design_lqr
lqr = design_lqr(A, B, Q, R, system_type='continuous')

The system properties (system.control, system.analysis) automatically handle backend management and provide a consistent interface. Direct function calls are possible but bypass the system’s backend configuration.

Architecture Philosophy

Functional Design with Composition - The control framework achieves:

Pure Functions - Stateless algorithms like scipy (design_lqr, analyze_stability)
Thin Wrappers - Minimal composition layer for system integration
Type Safety - TypedDict results for all algorithms
Backend Consistency - Automatic backend handling from parent system
Separation of Concerns - Analysis vs synthesis clearly separated
Mathematical Rigor - Implements classical control theory correctly
Clean Integration - Natural system.control and system.analysis APIs

Design principle: Control algorithms are pure functions that work like scipy—take matrices in, return structured results. Thin wrapper classes provide system integration without adding business logic.

Framework Layers

┌────────────────────────────────────────────────────────────┐
│                   APPLICATION LAYER                        │
│              (ContinuousSystemBase, DiscreteSystemBase)    │
│                                                            │
│  system.control   ─────► ControlSynthesis                 │
│  system.analysis  ─────► SystemAnalysis                   │
└──────────────────┬─────────────────────────────────────────┘
                   │ delegates to
                   ↓
┌────────────────────────────────────────────────────────────┐
│              PURE FUNCTION LAYER                           │
│          classical_control_functions.py                    │
│                                                            │
│  Control Design:          System Analysis:                 │
│  • design_lqr()          • analyze_stability()             │
│  • design_kalman()       • analyze_controllability()       │
│  • design_lqg()          • analyze_observability()         │
│                                                            │
│  All functions are stateless, pure, backend-agnostic       │
└──────────────────┬─────────────────────────────────────────┘
                   │ returns
                   ↓
┌────────────────────────────────────────────────────────────┐
│                   TYPE LAYER                               │
│             control_classical.py                           │
│                                                            │
│  • LQRResult           • StabilityInfo                     │
│  • KalmanFilterResult  • ControllabilityInfo               │
│  • LQGResult           • ObservabilityInfo                 │
└────────────────────────────────────────────────────────────┘

Pure Function Layer

classical_control_functions.py: Stateless Algorithms

File: classical_control_functions.py

The pure function layer provides stateless control algorithms that work like scipy—take matrices in, return structured TypedDict results.

Design philosophy:

Pure functions - No side effects, no state
Backend agnostic - Internal conversion to/from NumPy for scipy
Mathematical correctness - Rigorous implementation of control theory
Comprehensive validation - Dimension checks, positive-definiteness
Clear error handling - Actionable exceptions for infeasible problems

Key categories:

Control Design Functions

LQR Controller Design

Mathematical background:

Continuous-time LQR:

\[J = \int_0^\infty (x'Qx + u'Ru + 2x'Nu) dt\]

Algebraic Riccati Equation (ARE):

\[A'P + PA - PBR^{-1}B'P + Q - N'R^{-1}N = 0\]

Optimal gain:

\[K = R^{-1}(B'P + N')\]

Discrete-time LQR:

\[J = \sum_{k=0}^\infty (x[k]'Qx[k] + u[k]'Ru[k] + 2x[k]'Nu[k])\]

Discrete ARE:

\[P = A'PA - (A'PB + N)(R + B'PB)^{-1}(B'PA + N') + Q\]

Optimal gain:

\[K = (R + B'PB)^{-1}(B'PA + N')\]

Function signature:

def design_lqr(
    A: StateMatrix,
    B: InputMatrix,
    Q: StateMatrix,
    R: InputMatrix,
    N: Optional[InputMatrix] = None,
    system_type: str = "discrete",
    backend: Backend = "numpy"
) -> LQRResult:
    """
    Design Linear Quadratic Regulator (LQR) controller.
    
    Unified interface for continuous and discrete LQR.
    
    Parameters
    ----------
    A : StateMatrix
        State matrix (nx, nx)
    B : InputMatrix
        Input matrix (nx, nu)
    Q : StateMatrix
        State cost matrix (nx, nx), Q ≥ 0, (Q,A) detectable
    R : InputMatrix
        Control cost matrix (nu, nu), R > 0
    N : InputMatrix, optional
        Cross-coupling matrix (nx, nu)
    system_type : str
        'continuous' or 'discrete'
    backend : Backend
        Computational backend
    
    Returns
    -------
    LQRResult
        gain : Optimal feedback gain K (nu, nx)
        cost_to_go : Riccati solution P (nx, nx)
        closed_loop_eigenvalues : eig(A - BK)
        stability_margin : Stability robustness measure
    
    Raises
    ------
    ValueError
        If Q, R dimensions incompatible or conditions violated
    LinAlgError
        If Riccati equation has no stabilizing solution
    """

Usage example:

Code

system = Pendulum()

# Linearize at upright equilibrium
x_eq = np.array([np.pi, 0.0])
u_eq = np.array([0.0])
A, B = system.linearize(x_eq, u_eq)

# Design LQR controller
Q = np.diag([10, 1])   # Penalize angle more than velocity
R = np.array([[0.1]])   # Small control cost

lqr = system.control.design_lqr(A, B, Q, R, system_type='continuous')

# Extract results
K = lqr['gain']
P = lqr['cost_to_go']
eigs = lqr['controller_eigenvalues']

print(f"LQR gain shape: {K.shape}")
print(f"Closed-loop stable: {np.all(np.real(eigs) < 0)}")

Kalman Filter Design

Mathematical background:

System model:

\[x[k+1] = Ax[k] + Bu[k] + w[k], \quad w \sim \mathcal{N}(0, Q)\]

\[y[k] = Cx[k] + v[k], \quad v \sim \mathcal{N}(0, R)\]

Estimator dynamics:

\[\hat{x}[k+1] = A\hat{x}[k] + Bu[k] + L(y[k] - C\hat{x}[k])\]

Optimal gain:

\[L = APC'(CPC' + R)^{-1}\]

Error covariance Riccati equation:

\[P = A(P - PC'(CPC' + R)^{-1}CP)A' + Q\]

Function signature:

def design_kalman_filter(
    A: StateMatrix,
    C: OutputMatrix,
    Q: StateMatrix,
    R: OutputMatrix,
    system_type: str = "discrete",
    backend: Backend = "numpy"
) -> KalmanFilterResult:
    """
    Design Kalman filter for optimal state estimation.
    
    Parameters
    ----------
    A : StateMatrix
        State matrix (nx, nx)
    C : OutputMatrix
        Output matrix (ny, nx)
    Q : StateMatrix
        Process noise covariance (nx, nx), Q ≥ 0
    R : OutputMatrix
        Measurement noise covariance (ny, ny), R > 0
    system_type : str
        'continuous' or 'discrete'
    backend : Backend
        Computational backend
    
    Returns
    -------
    KalmanFilterResult
        gain : Kalman gain L (nx, ny)
        error_covariance : Steady-state P (nx, nx)
        innovation_covariance : S = CPC' + R (ny, ny)
        observer_eigenvalues : eig(A - LC)
    """

Usage example:

Code

# Define measurement model (measure angle only)
C = np.array([[1, 0]])

# Define noise covariances
Q_process = 0.01 * np.eye(2)      # Process noise
R_measurement = np.array([[0.1]])  # Measurement noise

# Design Kalman filter
kalman = system.control.design_kalman(
    A, C, Q_process, R_measurement,
    system_type='continuous'
)

# Extract results
L = kalman['gain']
P = kalman['estimation_error_covariance']
estimator_eigs = kalman['estimator_eigenvalues']

print(f"Kalman gain shape: {L.shape}")
print(f"Observer stable: {np.all(np.real(estimator_eigs) < 0)}")

LQG Controller Design

Mathematical background:

Separation principle:

Design LQR assuming full state feedback: $u = -Kx$
Design Kalman filter for state estimation
Combine via certainty equivalence: $u = -K\hat{x}$

Closed-loop system:

\[\begin{bmatrix} x[k+1] \\ e[k+1] \end{bmatrix} = \begin{bmatrix} A - BK & BK \\ 0 & A - LC \end{bmatrix} \begin{bmatrix} x[k] \\ e[k] \end{bmatrix} + \begin{bmatrix} w[k] \\ w[k] - Lv[k] \end{bmatrix}\]

Eigenvalues:

Controller poles: $\text{eig}(A - BK)$
Observer poles: $\text{eig}(A - LC)$
Combined: $\text{eig}(A - BK) \cup \text{eig}(A - LC)$

Function signature:

def design_lqg(
    A: StateMatrix,
    B: InputMatrix,
    C: OutputMatrix,
    Q_state: StateMatrix,
    R_control: InputMatrix,
    Q_process: StateMatrix,
    R_measurement: OutputMatrix,
    N: Optional[InputMatrix] = None,
    system_type: str = "discrete",
    backend: Backend = "numpy"
) -> LQGResult:
    """
    Design Linear Quadratic Gaussian (LQG) controller.
    
    Combines LQR (optimal control) with Kalman filter (optimal estimation)
    via the separation principle.
    
    Parameters
    ----------
    A : StateMatrix
        State matrix (nx, nx)
    B : InputMatrix
        Input matrix (nx, nu)
    C : OutputMatrix
        Output matrix (ny, nx)
    Q_state : StateMatrix
        LQR state cost (nx, nx)
    R_control : InputMatrix
        LQR control cost (nu, nu)
    Q_process : StateMatrix
        Process noise covariance (nx, nx)
    R_measurement : OutputMatrix
        Measurement noise covariance (ny, ny)
    N : InputMatrix, optional
        Cross-coupling (nx, nu)
    system_type : str
        'continuous' or 'discrete'
    backend : Backend
        Computational backend
    
    Returns
    -------
    LQGResult
        control_gain : LQR gain K (nu, nx)
        estimator_gain : Kalman gain L (nx, ny)
        control_cost_to_go : Controller Riccati P_c
        estimation_error_covariance : Estimator Riccati P_e
        separation_verified : bool
        closed_loop_stable : bool
        controller_eigenvalues : eig(A - BK)
        estimator_eigenvalues : eig(A - LC)
    """

Usage example:

Code

# Design LQG controller
lqg = system.control.design_lqg(
    A, B, C,
    Q_state=np.diag([10, 1]),
    R_control=np.array([[0.1]]),
    Q_process=0.01*np.eye(2),
    R_measurement=np.array([[0.1]]),
    system_type='continuous'
)

# Extract both gains
K_control = lqg['control_gain']
L_estimator = lqg['estimator_gain']

print(f"LQG stable: {lqg['closed_loop_stable']}")
print(f"Separation verified: {lqg['separation_verified']}")
print(f"Controller poles: {lqg['controller_eigenvalues']}")
print(f"Observer poles: {lqg['estimator_eigenvalues']}")

System Analysis Functions

Stability Analysis

Mathematical background:

Continuous-time stability criteria: - Asymptotically stable: All $\text{Re}(\lambda) < 0$ (left half-plane) - Marginally stable: $\max \text{Re}(\lambda) \approx 0$ - Unstable: Any $\text{Re}(\lambda) > 0$

Function signature:

def analyze_stability(
    A: StateMatrix,
    system_type: str = "continuous",
    backend: Backend = "numpy"
) -> StabilityInfo:
    """
    Analyze system stability via eigenvalue placement.
    
    Parameters
    ----------
    A : StateMatrix
        State matrix (nx, nx)
    system_type : str
        'continuous' or 'discrete'
    backend : Backend
        Computational backend
    
    Returns
    -------
    StabilityInfo
        eigenvalues : Complex eigenvalues of A
        magnitudes : |λ| for each eigenvalue
        max_magnitude : max |λ| (spectral radius)
        spectral_radius : Same as max_magnitude
        is_stable : Asymptotic stability
        is_marginally_stable : On boundary
        is_unstable : At least one unstable mode
    """

Usage example:

Code

# Analyze open-loop stability
stability_ol = system.analysis.stability(A, system_type='continuous')

print(f"Open-loop stable: {stability_ol['is_stable']}")
print(f"Spectral radius: {stability_ol['spectral_radius']:.3f}")
print(f"Eigenvalues: {stability_ol['eigenvalues']}")

# Analyze closed-loop stability (with LQR)
A_cl = A - B @ K
stability_cl = system.analysis.stability(A_cl, system_type='continuous')

print(f"\nClosed-loop stable: {stability_cl['is_stable']}")
print(f"Spectral radius: {stability_cl['spectral_radius']:.3f}")

Controllability Analysis

Mathematical background:

Controllability test:

\[\text{rank}(\mathcal{C}) = n_x \quad \text{where} \quad \mathcal{C} = [B, AB, A^2B, \ldots, A^{n_x-1}B]\]

Interpretation: - Controllable: All states can be driven to any value in finite time - Uncontrollable: Some states cannot be influenced by control

Function signature:

def analyze_controllability(
    A: StateMatrix,
    B: InputMatrix,
    backend: Backend = "numpy"
) -> ControllabilityInfo:
    """
    Test system controllability via rank condition.
    
    Parameters
    ----------
    A : StateMatrix
        State matrix (nx, nx)
    B : InputMatrix
        Input matrix (nx, nu)
    backend : Backend
        Computational backend
    
    Returns
    -------
    ControllabilityInfo
        controllability_matrix : C = [B AB ... A^(n-1)B] (nx, nx*nu)
        rank : Rank of controllability matrix
        is_controllable : rank == nx
        uncontrollable_modes : Eigenvalues of uncontrollable subsystem
    """

Usage example:

Code

# Check controllability
ctrl = system.analysis.controllability(A, B)

print(f"Controllable: {ctrl['is_controllable']}")
print(f"Rank: {ctrl['rank']} / {A.shape[0]}")

if ctrl['is_controllable']:
    print("✓ Can design LQR controller")
else:
    print("✗ Cannot design LQR - system not controllable")
    if ctrl.get('uncontrollable_modes') is not None:
        print(f"Uncontrollable modes: {ctrl['uncontrollable_modes']}")

Observability Analysis

Mathematical background:

Observability test:

\[\text{rank}(\mathcal{O}) = n_x \quad \text{where} \quad \mathcal{O} = \begin{bmatrix} C \\ CA \\ CA^2 \\ \vdots \\ CA^{n_x-1} \end{bmatrix}\]

Interpretation: - Observable: Initial state can be determined from output measurements - Unobservable: Some states hidden from measurements

Function signature:

def analyze_observability(
    A: StateMatrix,
    C: OutputMatrix,
    backend: Backend = "numpy"
) -> ObservabilityInfo:
    """
    Test system observability via rank condition.
    
    Parameters
    ----------
    A : StateMatrix
        State matrix (nx, nx)
    C : OutputMatrix
        Output matrix (ny, nx)
    backend : Backend
        Computational backend
    
    Returns
    -------
    ObservabilityInfo
        observability_matrix : O = [C; CA; ...] (nx*ny, nx)
        rank : Rank of observability matrix
        is_observable : rank == nx
        unobservable_modes : Eigenvalues of unobservable subsystem
    """

Usage example:

Code

# Check observability
obs = system.analysis.observability(A, C)

print(f"Observable: {obs['is_observable']}")
print(f"Rank: {obs['rank']} / {A.shape[0]}")

if obs['is_observable']:
    print("✓ Can design Kalman filter")
else:
    print("✗ Cannot design Kalman filter - system not observable")
    if obs.get('unobservable_modes') is not None:
        print(f"Unobservable modes: {obs['unobservable_modes']}")

# Check minimal realization
if ctrl['is_controllable'] and obs['is_observable']:
    print("\n✓ System is minimal (controllable and observable)")
    print("✓ Can design LQG controller")

Composition Wrapper Layer

ControlSynthesis: Control Design Wrapper

File: control_synthesis.py

The ControlSynthesis class provides a thin wrapper for control design algorithms, integrating them with the system’s backend configuration.

Design philosophy:

Composition not inheritance - Utility composed by system, not inherited
No state - Only stores backend setting from parent system
No caching - Delegates immediately to pure functions
Clean API - Methods match control theory terminology

Architecture:

class ControlSynthesis:
    """
    Control synthesis wrapper for system composition.
    
    Thin wrapper that routes to pure control design functions while
    maintaining backend consistency with parent system.
    """
    
    def __init__(self, backend: Backend = "numpy"):
        self.backend = backend
    
    def design_lqr(self, A, B, Q, R, N=None, system_type='discrete'):
        """Route to classical_control_functions.design_lqr()"""
        from cdesym.control.classical_control_functions import design_lqr
        return design_lqr(A, B, Q, R, N, system_type, self.backend)
    
    def design_kalman(self, A, C, Q, R, system_type='discrete'):
        """Route to classical_control_functions.design_kalman_filter()"""
        from cdesym.control.classical_control_functions import design_kalman_filter
        return design_kalman_filter(A, C, Q, R, system_type, self.backend)
    
    def design_lqg(self, A, B, C, Q_state, R_control, Q_process, R_measurement, 
                   N=None, system_type='discrete'):
        """Route to classical_control_functions.design_lqg()"""
        from cdesym.control.classical_control_functions import design_lqg
        return design_lqg(A, B, C, Q_state, R_control, Q_process, R_measurement, 
                         N, system_type, self.backend)

System integration:

# In ContinuousSystemBase / DiscreteSystemBase
@property
def control(self) -> ControlSynthesis:
    """Access control synthesis utilities."""
    if not hasattr(self, '_control_synthesis'):
        from cdesym.control.control_synthesis import ControlSynthesis
        self._control_synthesis = ControlSynthesis(backend=self.backend.default_backend)
    return self._control_synthesis

SystemAnalysis: System Analysis Wrapper

File: system_analysis.py

The SystemAnalysis class provides a thin wrapper for system analysis algorithms, identical in design to ControlSynthesis.

Architecture:

class SystemAnalysis:
    """
    System analysis wrapper for composition.
    
    Thin wrapper that routes to pure system analysis functions while
    maintaining backend consistency with parent system.
    """
    
    def __init__(self, backend: Backend = "numpy"):
        self.backend = backend
    
    def stability(self, A, system_type='continuous'):
        """Route to classical_control_functions.analyze_stability()"""
        from cdesym.control.classical_control_functions import analyze_stability
        return analyze_stability(A, system_type, self.backend)
    
    def controllability(self, A, B):
        """Route to classical_control_functions.analyze_controllability()"""
        from cdesym.control.classical_control_functions import analyze_controllability
        return analyze_controllability(A, B, self.backend)
    
    def observability(self, A, C):
        """Route to classical_control_functions.analyze_observability()"""
        from cdesym.control.classical_control_functions import analyze_observability
        return analyze_observability(A, C, self.backend)

System integration:

# In ContinuousSystemBase / DiscreteSystemBase
@property
def analysis(self) -> SystemAnalysis:
    """Access system analysis utilities."""
    if not hasattr(self, '_system_analysis'):
        from cdesym.control.system_analysis import SystemAnalysis
        self._system_analysis = SystemAnalysis(backend=self.backend.default_backend)
    return self._system_analysis

Design Patterns

Pattern 1: Pure Functions + Thin Wrappers

Principle: Separate stateless algorithms from system integration.

Code

# ✗ ANTI-PATTERN: Methods on system class (violates SRP)
class ContinuousSystemBase:
    def design_lqr(self, Q, R):
        # LQR implementation mixed with system concerns
        pass

# ✓ GOOD PATTERN: Pure function + composition

# Pure function (classical_control_functions.py)
def design_lqr(A, B, Q, R, system_type, backend):
    """Stateless, testable, reusable."""
    # Focus solely on LQR algorithm
    return LQRResult(...)

# Thin wrapper (control_synthesis.py)
class ControlSynthesis:
    def design_lqr(self, A, B, Q, R, system_type):
        return design_lqr(A, B, Q, R, system_type, self.backend)

# System integration (continuous_system_base.py)
@property
def control(self) -> ControlSynthesis:
    return ControlSynthesis(backend=self.backend.default_backend)

Benefits:

✓ Single responsibility per component
✓ Easy unit testing of pure functions
✓ Reusable algorithms outside system context
✓ Maintainability through isolation

Pattern 2: Backend Agnosticism

Principle: Convert to NumPy for scipy, then back to original backend.

def design_lqr(..., backend: Backend):
    """Works with NumPy, PyTorch, JAX transparently."""
    
    # Convert to NumPy for scipy
    A_np = _to_numpy(A, backend)
    B_np = _to_numpy(B, backend)
    Q_np = _to_numpy(Q, backend)
    R_np = _to_numpy(R, backend)
    
    # Solve in NumPy (scipy.linalg)
    P = solve_continuous_are(A_np, B_np, Q_np, R_np)
    K = np.linalg.solve(R_np, B_np.T @ P)
    
    # Convert back to original backend
    K_result = _from_numpy(K, backend)
    P_result = _from_numpy(P, backend)
    
    return LQRResult(gain=K_result, cost_to_go=P_result, ...)

Pattern 3: TypedDict Results

Principle: Structured results with type safety and IDE support.

Code

# All functions return structured TypedDict
result: LQRResult = design_lqr(A, B, Q, R, system_type='continuous')

# IDE autocomplete knows all fields
K = result['gain']                    # ✓ Valid
P = result['cost_to_go']              # ✓ Valid
eigs = result['controller_eigenvalues']  # ✓ Valid
bad = result['nonexistent_key']       # ✗ Type error!

# Type checking prevents errors
def apply_control(result: LQRResult) -> np.ndarray:
    return result['gain']  # ✓ Type checker verifies

Pattern 4: Unified Continuous/Discrete Interface

Principle: Single function handles both continuous and discrete cases.

def design_lqr(A, B, Q, R, N=None, system_type='discrete', backend='numpy'):
    """Unified interface - system_type selects algorithm."""
    
    if system_type == 'continuous':
        # Continuous-time algebraic Riccati equation
        P = solve_continuous_are(A, B, Q, R, s=N)
        K = np.linalg.solve(R, B.T @ P)
    elif system_type == 'discrete':
        # Discrete-time algebraic Riccati equation
        P = solve_discrete_are(A, B, Q, R, s=N)
        K = np.linalg.solve(R + B.T @ P @ B, B.T @ P @ A)
    else:
        raise ValueError(f"Invalid system_type: {system_type}")
    
    # Rest of implementation identical
    closed_loop_eigs = np.linalg.eigvals(A - B @ K)
    
    return LQRResult(gain=K, cost_to_go=P, ...)

Benefits:

✓ Single function for both cases
✓ Less code duplication
✓ Consistent API across system types
✓ Easier maintenance

Usage Workflows

Workflow 1: LQR Controller Design

Complete workflow from linearization to closed-loop simulation:

Code

# 1. Create system
system = Pendulum(m_val=1.0, l_val=0.5, g_val=9.81, beta_val=0.1)

# 2. Define equilibrium (upright position)
x_eq = np.array([np.pi, 0])  # [theta, omega]
u_eq = np.zeros(1)

# 3. Linearize at equilibrium
A, B = system.linearize(x_eq, u_eq)

# 4. Design LQR controller
Q = np.diag([10, 1])   # Penalize angle more than velocity
R = np.array([[0.1]])   # Small control cost

lqr = system.control.design_lqr(A, B, Q, R, system_type='continuous')

# 5. Extract and verify
K = lqr['gain']
eigs = lqr['controller_eigenvalues']

print(f"✓ LQR gain: {K}")
print(f"✓ Closed-loop stable: {np.all(np.real(eigs) < 0)}")

# 6. Implement controller
def lqr_controller(t, x):
    return -K @ (x - x_eq)

# 7. Simulate (closed-loop)
result = system.simulate(
    x0=np.array([np.pi + 0.2, 0]),  # Start near upright
    controller=lqr_controller,
    t_span=(0, 10),
    dt=0.01
)

print(f"✓ Simulation complete: {len(result['time'])} time points")

Workflow 2: Kalman Filter Design

Observer design for partial state measurements:

Code

# 1. Same system and linearization
A, B = system.linearize(x_eq, u_eq)

# 2. Define measurement model (measure angle only)
C = np.array([[1, 0]])  # Observe theta, not omega

# 3. Define noise covariances
Q_process = 0.01 * np.eye(2)      # Process noise
R_measurement = np.array([[0.1]])  # Measurement noise

# 4. Design Kalman filter
kalman = system.control.design_kalman(
    A, C, Q_process, R_measurement,
    system_type='continuous'
)

# 5. Extract gain
L = kalman['gain']
observer_eigs = kalman['estimator_eigenvalues']

print(f"✓ Kalman gain: {L}")
print(f"✓ Observer stable: {np.all(np.real(observer_eigs) < 0)}")

# 6. Implement estimator
x_hat = np.zeros(2)
dt = 0.01
N = 1000

for k in range(N):
    # Prediction
    x_hat_pred = A @ x_hat * dt + x_hat + B @ u[k] * dt
    
    # Measurement update
    y_meas = C @ x_true[k] + np.random.randn() * np.sqrt(0.1)
    innovation = y_meas - C @ x_hat_pred
    x_hat = x_hat_pred + L @ innovation * dt

Workflow 3: LQG Controller Design

Combined optimal control and estimation:

Code

# 1. Define all matrices
A, B = system.linearize(x_eq, u_eq)
C = np.array([[1, 0]])  # Partial state measurement

# 2. Design weights
Q_state = np.diag([10, 1])        # LQR state cost
R_control = np.array([[0.1]])      # LQR control cost
Q_process = 0.01 * np.eye(2)      # Kalman process noise
R_measurement = np.array([[0.1]])  # Kalman measurement noise

# 3. Design LQG
lqg = system.control.design_lqg(
    A, B, C,
    Q_state, R_control,
    Q_process, R_measurement,
    system_type='continuous'
)

# 4. Extract both gains
K = lqg['controller_gain']
L = lqg['estimator_gain']

print(f"✓ LQG stable: {lqg['closed_loop_stable']}")
print(f"✓ Separation verified: {lqg['separation_verified']}")

# 5. Implement LQG controller
x_hat = np.zeros(2)
dt = 0.01

for k in range(N):
    # Control (certainty equivalence)
    u[k] = -K @ (x_hat - x_eq)
    
    # Prediction
    x_hat = A @ x_hat * dt + x_hat + B @ u[k] * dt
    
    # Measurement update
    y_meas = C @ x_true[k] + measurement_noise[k]
    innovation = y_meas - C @ x_hat
    x_hat = x_hat + L @ innovation * dt

Workflow 4: System Analysis

Complete system analysis before controller design:

Code

# 1. Linearize system
A, B = system.linearize(x_eq, u_eq)
C = np.array([[1, 0]])

# 2. Check stability
stability = system.analysis.stability(A, system_type='continuous')
print(f"Open-loop stable: {stability['is_stable']}")
print(f"Eigenvalues: {stability['eigenvalues']}")
print(f"Spectral radius: {stability['spectral_radius']:.3f}")

# 3. Check controllability
ctrl = system.analysis.controllability(A, B)
print(f"\nControllable: {ctrl['is_controllable']}")
print(f"Rank: {ctrl['rank']} / {A.shape[0]}")

if not ctrl['is_controllable']:
    print(f"Uncontrollable modes: {ctrl['uncontrollable_modes']}")

# 4. Check observability
obs = system.analysis.observability(A, C)
print(f"\nObservable: {obs['is_observable']}")
print(f"Rank: {obs['rank']} / {A.shape[0]}")

# 5. Verify conditions for control design
if ctrl['is_controllable']:
    print("\n✓ Can design LQR controller")
else:
    print("\n✗ Cannot design LQR - system not controllable")

if obs['is_observable']:
    print("✓ Can design Kalman filter")
else:
    print("✗ Cannot design Kalman filter - system not observable")

if ctrl['is_controllable'] and obs['is_observable']:
    print("✓ System is minimal - can design LQG controller")

Key Strengths

Control Framework Advantages

Pure Functional Core - Stateless algorithms, easy to test and reuse
Thin Wrappers - Minimal composition layer with no business logic
Type Safety - TypedDict results throughout with IDE support
Backend Agnostic - NumPy/PyTorch/JAX transparency
Separation of Concerns - Analysis vs synthesis clearly separated
Mathematical Rigor - Correct implementation of classical control theory
Clean Integration - Natural system.control and system.analysis APIs
Unified Interface - Single function for continuous/discrete cases
Comprehensive - LQR, Kalman, LQG, stability, controllability, observability
Scipy-like - Familiar API for control engineers

Summary

The Control Framework provides production-quality classical control theory algorithms through a clean two-layer architecture:

Pure Function Layer provides stateless, testable algorithms
Composition Layer integrates algorithms with system backend management
Type Layer ensures type safety with structured TypedDict results

Users interact through intuitive system properties (system.control, system.analysis) that automatically handle backend management and provide consistent interfaces across continuous and discrete systems.

The framework implements classical control theory correctly while maintaining clean software architecture and comprehensive type safety.

--- title: "Control Framework Architecture" subtitle: "Classical Control Theory for Analysis and Synthesis" author: "Gil Benezer" date: today format: html: toc: true toc-depth: 5 code-fold: show code-tools: true theme: cosmo execute: eval: true cache: true warning: false --- ## Overview {#sec-overview} The **Control Framework** provides classical control theory algorithms for analysis and synthesis of dynamical systems. It consists of **3 core modules** organized into a clean **2-layer architecture**: pure stateless functions and thin composition wrappers. ::: {.callout-note} ## Framework Components **Pure Function Layer:** `classical_control_functions.py` - Stateless control algorithms **Composition Layer:** `control_synthesis.py`, `system_analysis.py` - System integration **Type Layer:** `control_classical.py` - TypedDict result definitions ::: ::: {.callout-important} ## User Interaction Model **Users interact with the control framework through system properties:** ```python # ✓ CORRECT: Access via system properties system = Pendulum() A, B = system.linearize(x_eq, u_eq) # Control design lqr = system.control.design_lqr(A, B, Q, R, system_type='continuous') K = lqr['gain'] # System analysis stability = system.analysis.stability(A, system_type='continuous') controllability = system.analysis.controllability(A, B) # ✗ INCORRECT: Direct function imports (not recommended) from cdesym.control.classical_control_functions import design_lqr lqr = design_lqr(A, B, Q, R, system_type='continuous') ``` The system properties (`system.control`, `system.analysis`) automatically handle backend management and provide a consistent interface. Direct function calls are possible but bypass the system's backend configuration. ::: ```{python} #| label: setup-imports #| echo: false #| output: false # Global setup for all code examples import numpy as np from typing import Optional, Tuple, Dict, List, Union import warnings warnings.filterwarnings('ignore', category=UserWarning) # Import for examples from cdesym import Pendulum ``` ## Architecture Philosophy {#sec-architecture-philosophy} **Functional Design with Composition** - The control framework achieves: 1. **Pure Functions** - Stateless algorithms like scipy (design_lqr, analyze_stability) 2. **Thin Wrappers** - Minimal composition layer for system integration 3. **Type Safety** - TypedDict results for all algorithms 4. **Backend Consistency** - Automatic backend handling from parent system 5. **Separation of Concerns** - Analysis vs synthesis clearly separated 6. **Mathematical Rigor** - Implements classical control theory correctly 7. **Clean Integration** - Natural `system.control` and `system.analysis` APIs **Design principle:** Control algorithms are **pure functions** that work like scipy—take matrices in, return structured results. Thin wrapper classes provide system integration without adding business logic. ## Framework Layers {#sec-framework-layers} ``` ┌────────────────────────────────────────────────────────────┐ │ APPLICATION LAYER │ │ (ContinuousSystemBase, DiscreteSystemBase) │ │ │ │ system.control ─────► ControlSynthesis │ │ system.analysis ─────► SystemAnalysis │ └──────────────────┬─────────────────────────────────────────┘ │ delegates to ↓ ┌────────────────────────────────────────────────────────────┐ │ PURE FUNCTION LAYER │ │ classical_control_functions.py │ │ │ │ Control Design: System Analysis: │ │ • design_lqr() • analyze_stability() │ │ • design_kalman() • analyze_controllability() │ │ • design_lqg() • analyze_observability() │ │ │ │ All functions are stateless, pure, backend-agnostic │ └──────────────────┬─────────────────────────────────────────┘ │ returns ↓ ┌────────────────────────────────────────────────────────────┐ │ TYPE LAYER │ │ control_classical.py │ │ │ │ • LQRResult • StabilityInfo │ │ • KalmanFilterResult • ControllabilityInfo │ │ • LQGResult • ObservabilityInfo │ └────────────────────────────────────────────────────────────┘ ``` ## Pure Function Layer {#sec-pure-function-layer} ### classical_control_functions.py: Stateless Algorithms {#sec-classical-control-functions} **File:** `classical_control_functions.py` The pure function layer provides stateless control algorithms that work like scipy—take matrices in, return structured TypedDict results. **Design philosophy:** - **Pure functions** - No side effects, no state - **Backend agnostic** - Internal conversion to/from NumPy for scipy - **Mathematical correctness** - Rigorous implementation of control theory - **Comprehensive validation** - Dimension checks, positive-definiteness - **Clear error handling** - Actionable exceptions for infeasible problems **Key categories:** ### Control Design Functions {#sec-control-design} #### LQR Controller Design **Mathematical background:** **Continuous-time LQR:** $$J = \int_0^\infty (x'Qx + u'Ru + 2x'Nu) dt$$ **Algebraic Riccati Equation (ARE):** $$A'P + PA - PBR^{-1}B'P + Q - N'R^{-1}N = 0$$ **Optimal gain:** $$K = R^{-1}(B'P + N')$$ **Discrete-time LQR:** $$J = \sum_{k=0}^\infty (x[k]'Qx[k] + u[k]'Ru[k] + 2x[k]'Nu[k])$$ **Discrete ARE:** $$P = A'PA - (A'PB + N)(R + B'PB)^{-1}(B'PA + N') + Q$$ **Optimal gain:** $$K = (R + B'PB)^{-1}(B'PA + N')$$ **Function signature:** ```python def design_lqr( A: StateMatrix, B: InputMatrix, Q: StateMatrix, R: InputMatrix, N: Optional[InputMatrix] = None, system_type: str = "discrete", backend: Backend = "numpy" ) -> LQRResult: """ Design Linear Quadratic Regulator (LQR) controller. Unified interface for continuous and discrete LQR. Parameters ---------- A : StateMatrix State matrix (nx, nx) B : InputMatrix Input matrix (nx, nu) Q : StateMatrix State cost matrix (nx, nx), Q ≥ 0, (Q,A) detectable R : InputMatrix Control cost matrix (nu, nu), R > 0 N : InputMatrix, optional Cross-coupling matrix (nx, nu) system_type : str 'continuous' or 'discrete' backend : Backend Computational backend Returns ------- LQRResult gain : Optimal feedback gain K (nu, nx) cost_to_go : Riccati solution P (nx, nx) closed_loop_eigenvalues : eig(A - BK) stability_margin : Stability robustness measure Raises ------ ValueError If Q, R dimensions incompatible or conditions violated LinAlgError If Riccati equation has no stabilizing solution """ ``` **Usage example:** ```{python} #| label: ex-lqr-design system = Pendulum() # Linearize at upright equilibrium x_eq = np.array([np.pi, 0.0]) u_eq = np.array([0.0]) A, B = system.linearize(x_eq, u_eq) # Design LQR controller Q = np.diag([10, 1]) # Penalize angle more than velocity R = np.array([[0.1]]) # Small control cost lqr = system.control.design_lqr(A, B, Q, R, system_type='continuous') # Extract results K = lqr['gain'] P = lqr['cost_to_go'] eigs = lqr['controller_eigenvalues'] print(f"LQR gain shape: {K.shape}") print(f"Closed-loop stable: {np.all(np.real(eigs) < 0)}") ``` #### Kalman Filter Design **Mathematical background:** **System model:** $$x[k+1] = Ax[k] + Bu[k] + w[k], \quad w \sim \mathcal{N}(0, Q)$$ $$y[k] = Cx[k] + v[k], \quad v \sim \mathcal{N}(0, R)$$ **Estimator dynamics:** $$\hat{x}[k+1] = A\hat{x}[k] + Bu[k] + L(y[k] - C\hat{x}[k])$$ **Optimal gain:** $$L = APC'(CPC' + R)^{-1}$$ **Error covariance Riccati equation:** $$P = A(P - PC'(CPC' + R)^{-1}CP)A' + Q$$ **Function signature:** ```python def design_kalman_filter( A: StateMatrix, C: OutputMatrix, Q: StateMatrix, R: OutputMatrix, system_type: str = "discrete", backend: Backend = "numpy" ) -> KalmanFilterResult: """ Design Kalman filter for optimal state estimation. Parameters ---------- A : StateMatrix State matrix (nx, nx) C : OutputMatrix Output matrix (ny, nx) Q : StateMatrix Process noise covariance (nx, nx), Q ≥ 0 R : OutputMatrix Measurement noise covariance (ny, ny), R > 0 system_type : str 'continuous' or 'discrete' backend : Backend Computational backend Returns ------- KalmanFilterResult gain : Kalman gain L (nx, ny) error_covariance : Steady-state P (nx, nx) innovation_covariance : S = CPC' + R (ny, ny) observer_eigenvalues : eig(A - LC) """ ``` **Usage example:** ```{python} #| label: ex-kalman-design # Define measurement model (measure angle only) C = np.array([[1, 0]]) # Define noise covariances Q_process = 0.01 * np.eye(2) # Process noise R_measurement = np.array([[0.1]]) # Measurement noise # Design Kalman filter kalman = system.control.design_kalman( A, C, Q_process, R_measurement, system_type='continuous' ) # Extract results L = kalman['gain'] P = kalman['estimation_error_covariance'] estimator_eigs = kalman['estimator_eigenvalues'] print(f"Kalman gain shape: {L.shape}") print(f"Observer stable: {np.all(np.real(estimator_eigs) < 0)}") ``` #### LQG Controller Design **Mathematical background:** **Separation principle:** 1. Design LQR assuming full state feedback: $u = -Kx$ 2. Design Kalman filter for state estimation 3. Combine via certainty equivalence: $u = -K\hat{x}$ **Closed-loop system:** $$\begin{bmatrix} x[k+1] \\ e[k+1] \end{bmatrix} = \begin{bmatrix} A - BK & BK \\ 0 & A - LC \end{bmatrix} \begin{bmatrix} x[k] \\ e[k] \end{bmatrix} + \begin{bmatrix} w[k] \\ w[k] - Lv[k] \end{bmatrix}$$ **Eigenvalues:** - Controller poles: $\text{eig}(A - BK)$ - Observer poles: $\text{eig}(A - LC)$ - Combined: $\text{eig}(A - BK) \cup \text{eig}(A - LC)$ **Function signature:** ```python def design_lqg( A: StateMatrix, B: InputMatrix, C: OutputMatrix, Q_state: StateMatrix, R_control: InputMatrix, Q_process: StateMatrix, R_measurement: OutputMatrix, N: Optional[InputMatrix] = None, system_type: str = "discrete", backend: Backend = "numpy" ) -> LQGResult: """ Design Linear Quadratic Gaussian (LQG) controller. Combines LQR (optimal control) with Kalman filter (optimal estimation) via the separation principle. Parameters ---------- A : StateMatrix State matrix (nx, nx) B : InputMatrix Input matrix (nx, nu) C : OutputMatrix Output matrix (ny, nx) Q_state : StateMatrix LQR state cost (nx, nx) R_control : InputMatrix LQR control cost (nu, nu) Q_process : StateMatrix Process noise covariance (nx, nx) R_measurement : OutputMatrix Measurement noise covariance (ny, ny) N : InputMatrix, optional Cross-coupling (nx, nu) system_type : str 'continuous' or 'discrete' backend : Backend Computational backend Returns ------- LQGResult control_gain : LQR gain K (nu, nx) estimator_gain : Kalman gain L (nx, ny) control_cost_to_go : Controller Riccati P_c estimation_error_covariance : Estimator Riccati P_e separation_verified : bool closed_loop_stable : bool controller_eigenvalues : eig(A - BK) estimator_eigenvalues : eig(A - LC) """ ``` **Usage example:** ```{python} #| label: ex-lqg-design # Design LQG controller lqg = system.control.design_lqg( A, B, C, Q_state=np.diag([10, 1]), R_control=np.array([[0.1]]), Q_process=0.01*np.eye(2), R_measurement=np.array([[0.1]]), system_type='continuous' ) # Extract both gains K_control = lqg['control_gain'] L_estimator = lqg['estimator_gain'] print(f"LQG stable: {lqg['closed_loop_stable']}") print(f"Separation verified: {lqg['separation_verified']}") print(f"Controller poles: {lqg['controller_eigenvalues']}") print(f"Observer poles: {lqg['estimator_eigenvalues']}") ``` ### System Analysis Functions {#sec-system-analysis} #### Stability Analysis **Mathematical background:** **Continuous-time stability criteria:** - **Asymptotically stable:** All $\text{Re}(\lambda) < 0$ (left half-plane) - **Marginally stable:** $\max \text{Re}(\lambda) \approx 0$ - **Unstable:** Any $\text{Re}(\lambda) > 0$ **Discrete-time stability criteria:** - **Asymptotically stable:** All $|\lambda| < 1$ (inside unit circle) - **Marginally stable:** $\max |\lambda| \approx 1$ - **Unstable:** Any $|\lambda| > 1$ **Function signature:** ```python def analyze_stability( A: StateMatrix, system_type: str = "continuous", backend: Backend = "numpy" ) -> StabilityInfo: """ Analyze system stability via eigenvalue placement. Parameters ---------- A : StateMatrix State matrix (nx, nx) system_type : str 'continuous' or 'discrete' backend : Backend Computational backend Returns ------- StabilityInfo eigenvalues : Complex eigenvalues of A magnitudes : |λ| for each eigenvalue max_magnitude : max |λ| (spectral radius) spectral_radius : Same as max_magnitude is_stable : Asymptotic stability is_marginally_stable : On boundary is_unstable : At least one unstable mode """ ``` **Usage example:** ```{python} #| label: ex-stability-analysis # Analyze open-loop stability stability_ol = system.analysis.stability(A, system_type='continuous') print(f"Open-loop stable: {stability_ol['is_stable']}") print(f"Spectral radius: {stability_ol['spectral_radius']:.3f}") print(f"Eigenvalues: {stability_ol['eigenvalues']}") # Analyze closed-loop stability (with LQR) A_cl = A - B @ K stability_cl = system.analysis.stability(A_cl, system_type='continuous') print(f"\nClosed-loop stable: {stability_cl['is_stable']}") print(f"Spectral radius: {stability_cl['spectral_radius']:.3f}") ``` #### Controllability Analysis **Mathematical background:** **Controllability test:** $$\text{rank}(\mathcal{C}) = n_x \quad \text{where} \quad \mathcal{C} = [B, AB, A^2B, \ldots, A^{n_x-1}B]$$ **Interpretation:** - **Controllable:** All states can be driven to any value in finite time - **Uncontrollable:** Some states cannot be influenced by control **Function signature:** ```python def analyze_controllability( A: StateMatrix, B: InputMatrix, backend: Backend = "numpy" ) -> ControllabilityInfo: """ Test system controllability via rank condition. Parameters ---------- A : StateMatrix State matrix (nx, nx) B : InputMatrix Input matrix (nx, nu) backend : Backend Computational backend Returns ------- ControllabilityInfo controllability_matrix : C = [B AB ... A^(n-1)B] (nx, nx*nu) rank : Rank of controllability matrix is_controllable : rank == nx uncontrollable_modes : Eigenvalues of uncontrollable subsystem """ ``` **Usage example:** ```{python} #| label: ex-controllability-analysis # Check controllability ctrl = system.analysis.controllability(A, B) print(f"Controllable: {ctrl['is_controllable']}") print(f"Rank: {ctrl['rank']} / {A.shape[0]}") if ctrl['is_controllable']: print("✓ Can design LQR controller") else: print("✗ Cannot design LQR - system not controllable") if ctrl.get('uncontrollable_modes') is not None: print(f"Uncontrollable modes: {ctrl['uncontrollable_modes']}") ``` #### Observability Analysis **Mathematical background:** **Observability test:** $$\text{rank}(\mathcal{O}) = n_x \quad \text{where} \quad \mathcal{O} = \begin{bmatrix} C \\ CA \\ CA^2 \\ \vdots \\ CA^{n_x-1} \end{bmatrix}$$ **Interpretation:** - **Observable:** Initial state can be determined from output measurements - **Unobservable:** Some states hidden from measurements **Function signature:** ```python def analyze_observability( A: StateMatrix, C: OutputMatrix, backend: Backend = "numpy" ) -> ObservabilityInfo: """ Test system observability via rank condition. Parameters ---------- A : StateMatrix State matrix (nx, nx) C : OutputMatrix Output matrix (ny, nx) backend : Backend Computational backend Returns ------- ObservabilityInfo observability_matrix : O = [C; CA; ...] (nx*ny, nx) rank : Rank of observability matrix is_observable : rank == nx unobservable_modes : Eigenvalues of unobservable subsystem """ ``` **Usage example:** ```{python} #| label: ex-observability-analysis # Check observability obs = system.analysis.observability(A, C) print(f"Observable: {obs['is_observable']}") print(f"Rank: {obs['rank']} / {A.shape[0]}") if obs['is_observable']: print("✓ Can design Kalman filter") else: print("✗ Cannot design Kalman filter - system not observable") if obs.get('unobservable_modes') is not None: print(f"Unobservable modes: {obs['unobservable_modes']}") # Check minimal realization if ctrl['is_controllable'] and obs['is_observable']: print("\n✓ System is minimal (controllable and observable)") print("✓ Can design LQG controller") ``` ## Composition Wrapper Layer {#sec-composition-wrapper-layer} ### ControlSynthesis: Control Design Wrapper {#sec-control-synthesis} **File:** `control_synthesis.py` The `ControlSynthesis` class provides a thin wrapper for control design algorithms, integrating them with the system's backend configuration. **Design philosophy:** - **Composition not inheritance** - Utility composed by system, not inherited - **No state** - Only stores backend setting from parent system - **No caching** - Delegates immediately to pure functions - **Clean API** - Methods match control theory terminology **Architecture:** ```python class ControlSynthesis: """ Control synthesis wrapper for system composition. Thin wrapper that routes to pure control design functions while maintaining backend consistency with parent system. """ def __init__(self, backend: Backend = "numpy"): self.backend = backend def design_lqr(self, A, B, Q, R, N=None, system_type='discrete'): """Route to classical_control_functions.design_lqr()""" from cdesym.control.classical_control_functions import design_lqr return design_lqr(A, B, Q, R, N, system_type, self.backend) def design_kalman(self, A, C, Q, R, system_type='discrete'): """Route to classical_control_functions.design_kalman_filter()""" from cdesym.control.classical_control_functions import design_kalman_filter return design_kalman_filter(A, C, Q, R, system_type, self.backend) def design_lqg(self, A, B, C, Q_state, R_control, Q_process, R_measurement, N=None, system_type='discrete'): """Route to classical_control_functions.design_lqg()""" from cdesym.control.classical_control_functions import design_lqg return design_lqg(A, B, C, Q_state, R_control, Q_process, R_measurement, N, system_type, self.backend) ``` **System integration:** ```python # In ContinuousSystemBase / DiscreteSystemBase @property def control(self) -> ControlSynthesis: """Access control synthesis utilities.""" if not hasattr(self, '_control_synthesis'): from cdesym.control.control_synthesis import ControlSynthesis self._control_synthesis = ControlSynthesis(backend=self.backend.default_backend) return self._control_synthesis ``` ### SystemAnalysis: System Analysis Wrapper {#sec-system-analysis-wrapper} **File:** `system_analysis.py` The `SystemAnalysis` class provides a thin wrapper for system analysis algorithms, identical in design to `ControlSynthesis`. **Architecture:** ```python class SystemAnalysis: """ System analysis wrapper for composition. Thin wrapper that routes to pure system analysis functions while maintaining backend consistency with parent system. """ def __init__(self, backend: Backend = "numpy"): self.backend = backend def stability(self, A, system_type='continuous'): """Route to classical_control_functions.analyze_stability()""" from cdesym.control.classical_control_functions import analyze_stability return analyze_stability(A, system_type, self.backend) def controllability(self, A, B): """Route to classical_control_functions.analyze_controllability()""" from cdesym.control.classical_control_functions import analyze_controllability return analyze_controllability(A, B, self.backend) def observability(self, A, C): """Route to classical_control_functions.analyze_observability()""" from cdesym.control.classical_control_functions import analyze_observability return analyze_observability(A, C, self.backend) ``` **System integration:** ```python # In ContinuousSystemBase / DiscreteSystemBase @property def analysis(self) -> SystemAnalysis: """Access system analysis utilities.""" if not hasattr(self, '_system_analysis'): from cdesym.control.system_analysis import SystemAnalysis self._system_analysis = SystemAnalysis(backend=self.backend.default_backend) return self._system_analysis ``` ## Design Patterns {#sec-design-patterns} ### Pattern 1: Pure Functions + Thin Wrappers {#sec-pattern-pure-functions} **Principle:** Separate stateless algorithms from system integration. ```{python} #| label: ex-pattern-separation #| eval: false # ✗ ANTI-PATTERN: Methods on system class (violates SRP) class ContinuousSystemBase: def design_lqr(self, Q, R): # LQR implementation mixed with system concerns pass # ✓ GOOD PATTERN: Pure function + composition # Pure function (classical_control_functions.py) def design_lqr(A, B, Q, R, system_type, backend): """Stateless, testable, reusable.""" # Focus solely on LQR algorithm return LQRResult(...) # Thin wrapper (control_synthesis.py) class ControlSynthesis: def design_lqr(self, A, B, Q, R, system_type): return design_lqr(A, B, Q, R, system_type, self.backend) # System integration (continuous_system_base.py) @property def control(self) -> ControlSynthesis: return ControlSynthesis(backend=self.backend.default_backend) ``` **Benefits:** - ✓ Single responsibility per component - ✓ Easy unit testing of pure functions - ✓ Reusable algorithms outside system context - ✓ Maintainability through isolation ### Pattern 2: Backend Agnosticism {#sec-pattern-backend-agnostic} **Principle:** Convert to NumPy for scipy, then back to original backend. ```python def design_lqr(..., backend: Backend): """Works with NumPy, PyTorch, JAX transparently.""" # Convert to NumPy for scipy A_np = _to_numpy(A, backend) B_np = _to_numpy(B, backend) Q_np = _to_numpy(Q, backend) R_np = _to_numpy(R, backend) # Solve in NumPy (scipy.linalg) P = solve_continuous_are(A_np, B_np, Q_np, R_np) K = np.linalg.solve(R_np, B_np.T @ P) # Convert back to original backend K_result = _from_numpy(K, backend) P_result = _from_numpy(P, backend) return LQRResult(gain=K_result, cost_to_go=P_result, ...) ``` ### Pattern 3: TypedDict Results {#sec-pattern-typeddict-results} **Principle:** Structured results with type safety and IDE support. ```{python} #| label: ex-pattern-typeddict #| eval: false # All functions return structured TypedDict result: LQRResult = design_lqr(A, B, Q, R, system_type='continuous') # IDE autocomplete knows all fields K = result['gain'] # ✓ Valid P = result['cost_to_go'] # ✓ Valid eigs = result['controller_eigenvalues'] # ✓ Valid bad = result['nonexistent_key'] # ✗ Type error! # Type checking prevents errors def apply_control(result: LQRResult) -> np.ndarray: return result['gain'] # ✓ Type checker verifies ``` ### Pattern 4: Unified Continuous/Discrete Interface {#sec-pattern-unified-interface} **Principle:** Single function handles both continuous and discrete cases. ```python def design_lqr(A, B, Q, R, N=None, system_type='discrete', backend='numpy'): """Unified interface - system_type selects algorithm.""" if system_type == 'continuous': # Continuous-time algebraic Riccati equation P = solve_continuous_are(A, B, Q, R, s=N) K = np.linalg.solve(R, B.T @ P) elif system_type == 'discrete': # Discrete-time algebraic Riccati equation P = solve_discrete_are(A, B, Q, R, s=N) K = np.linalg.solve(R + B.T @ P @ B, B.T @ P @ A) else: raise ValueError(f"Invalid system_type: {system_type}") # Rest of implementation identical closed_loop_eigs = np.linalg.eigvals(A - B @ K) return LQRResult(gain=K, cost_to_go=P, ...) ``` **Benefits:** - ✓ Single function for both cases - ✓ Less code duplication - ✓ Consistent API across system types - ✓ Easier maintenance ## Usage Workflows {#sec-usage-workflows} ### Workflow 1: LQR Controller Design {#sec-workflow-lqr} **Complete workflow from linearization to closed-loop simulation:** ```{python} #| label: ex-workflow-lqr # 1. Create system system = Pendulum(m_val=1.0, l_val=0.5, g_val=9.81, beta_val=0.1) # 2. Define equilibrium (upright position) x_eq = np.array([np.pi, 0]) # [theta, omega] u_eq = np.zeros(1) # 3. Linearize at equilibrium A, B = system.linearize(x_eq, u_eq) # 4. Design LQR controller Q = np.diag([10, 1]) # Penalize angle more than velocity R = np.array([[0.1]]) # Small control cost lqr = system.control.design_lqr(A, B, Q, R, system_type='continuous') # 5. Extract and verify K = lqr['gain'] eigs = lqr['controller_eigenvalues'] print(f"✓ LQR gain: {K}") print(f"✓ Closed-loop stable: {np.all(np.real(eigs) < 0)}") # 6. Implement controller def lqr_controller(t, x): return -K @ (x - x_eq) # 7. Simulate (closed-loop) result = system.simulate( x0=np.array([np.pi + 0.2, 0]), # Start near upright controller=lqr_controller, t_span=(0, 10), dt=0.01 ) print(f"✓ Simulation complete: {len(result['time'])} time points") ``` ### Workflow 2: Kalman Filter Design {#sec-workflow-kalman} **Observer design for partial state measurements:** ```{python} #| label: ex-workflow-kalman #| eval: false # 1. Same system and linearization A, B = system.linearize(x_eq, u_eq) # 2. Define measurement model (measure angle only) C = np.array([[1, 0]]) # Observe theta, not omega # 3. Define noise covariances Q_process = 0.01 * np.eye(2) # Process noise R_measurement = np.array([[0.1]]) # Measurement noise # 4. Design Kalman filter kalman = system.control.design_kalman( A, C, Q_process, R_measurement, system_type='continuous' ) # 5. Extract gain L = kalman['gain'] observer_eigs = kalman['estimator_eigenvalues'] print(f"✓ Kalman gain: {L}") print(f"✓ Observer stable: {np.all(np.real(observer_eigs) < 0)}") # 6. Implement estimator x_hat = np.zeros(2) dt = 0.01 N = 1000 for k in range(N): # Prediction x_hat_pred = A @ x_hat * dt + x_hat + B @ u[k] * dt # Measurement update y_meas = C @ x_true[k] + np.random.randn() * np.sqrt(0.1) innovation = y_meas - C @ x_hat_pred x_hat = x_hat_pred + L @ innovation * dt ``` ### Workflow 3: LQG Controller Design {#sec-workflow-lqg} **Combined optimal control and estimation:** ```{python} #| label: ex-workflow-lqg #| eval: false # 1. Define all matrices A, B = system.linearize(x_eq, u_eq) C = np.array([[1, 0]]) # Partial state measurement # 2. Design weights Q_state = np.diag([10, 1]) # LQR state cost R_control = np.array([[0.1]]) # LQR control cost Q_process = 0.01 * np.eye(2) # Kalman process noise R_measurement = np.array([[0.1]]) # Kalman measurement noise # 3. Design LQG lqg = system.control.design_lqg( A, B, C, Q_state, R_control, Q_process, R_measurement, system_type='continuous' ) # 4. Extract both gains K = lqg['controller_gain'] L = lqg['estimator_gain'] print(f"✓ LQG stable: {lqg['closed_loop_stable']}") print(f"✓ Separation verified: {lqg['separation_verified']}") # 5. Implement LQG controller x_hat = np.zeros(2) dt = 0.01 for k in range(N): # Control (certainty equivalence) u[k] = -K @ (x_hat - x_eq) # Prediction x_hat = A @ x_hat * dt + x_hat + B @ u[k] * dt # Measurement update y_meas = C @ x_true[k] + measurement_noise[k] innovation = y_meas - C @ x_hat x_hat = x_hat + L @ innovation * dt ``` ### Workflow 4: System Analysis {#sec-workflow-analysis} **Complete system analysis before controller design:** ```{python} #| label: ex-workflow-analysis #| eval: false # 1. Linearize system A, B = system.linearize(x_eq, u_eq) C = np.array([[1, 0]]) # 2. Check stability stability = system.analysis.stability(A, system_type='continuous') print(f"Open-loop stable: {stability['is_stable']}") print(f"Eigenvalues: {stability['eigenvalues']}") print(f"Spectral radius: {stability['spectral_radius']:.3f}") # 3. Check controllability ctrl = system.analysis.controllability(A, B) print(f"\nControllable: {ctrl['is_controllable']}") print(f"Rank: {ctrl['rank']} / {A.shape[0]}") if not ctrl['is_controllable']: print(f"Uncontrollable modes: {ctrl['uncontrollable_modes']}") # 4. Check observability obs = system.analysis.observability(A, C) print(f"\nObservable: {obs['is_observable']}") print(f"Rank: {obs['rank']} / {A.shape[0]}") # 5. Verify conditions for control design if ctrl['is_controllable']: print("\n✓ Can design LQR controller") else: print("\n✗ Cannot design LQR - system not controllable") if obs['is_observable']: print("✓ Can design Kalman filter") else: print("✗ Cannot design Kalman filter - system not observable") if ctrl['is_controllable'] and obs['is_observable']: print("✓ System is minimal - can design LQG controller") ``` ## Key Strengths {#sec-key-strengths} ::: {.callout-tip} ## Control Framework Advantages 1. **Pure Functional Core** - Stateless algorithms, easy to test and reuse 2. **Thin Wrappers** - Minimal composition layer with no business logic 3. **Type Safety** - TypedDict results throughout with IDE support 4. **Backend Agnostic** - NumPy/PyTorch/JAX transparency 5. **Separation of Concerns** - Analysis vs synthesis clearly separated 6. **Mathematical Rigor** - Correct implementation of classical control theory 7. **Clean Integration** - Natural `system.control` and `system.analysis` APIs 8. **Unified Interface** - Single function for continuous/discrete cases 9. **Comprehensive** - LQR, Kalman, LQG, stability, controllability, observability 10. **Scipy-like** - Familiar API for control engineers ::: ## Summary The Control Framework provides **production-quality classical control theory algorithms** through a clean two-layer architecture: - **Pure Function Layer** provides stateless, testable algorithms - **Composition Layer** integrates algorithms with system backend management - **Type Layer** ensures type safety with structured TypedDict results Users interact through intuitive system properties (`system.control`, `system.analysis`) that automatically handle backend management and provide consistent interfaces across continuous and discrete systems. **The framework implements classical control theory correctly while maintaining clean software architecture and comprehensive type safety.**