Lecture 14 - Introduction to Feedback Linearization

Based on notes created by Sam Coogan and Murat Arcak. Licensed under a “Creative Commons Attribution-NonCommercial-ShareAlike 4.0 International License”

Additional Reading

• Khalil, Chapter 12 (Feedback Control)
• Khalil, Chapter 13 (Feedback Linearization)

Overview

• Control Systems and Control Objectives
• Motivating Feedback Linearization

Motivating Feedback Linearization

Feedback linearization is an extremely powerful tool in nonlinear control because it allows one (under certain conditions) to exactly linearize a nonlinear system via nonlinear control.

Control Systems and Control Objectives

Definition: Control System

Consider an open and connected set D \subseteq \mathbb{R}^n and a set of admissible control inputs U \subset \mathbb{R}^m. A control system is given by a differential equation:

\dot{x} = f(x, u), \quad x \in D, u \in U

where f: D \times U \to \mathbb{R}^n is a C^1 function.

Definition: Affine Control System

A control system is an affine control system (or control-affine) if it can be written in the form:

\dot{x} = f(x) + g(x) u

where f: D \to \mathbb{R}^n and g: D \to \mathbb{R}^{n \times m} are C^1 functions. Here, f is sometimes called the drift and g is sometimes called the actuation matrix or input matrix.

Definition: Outputs

An output is a differentiable function y: D \to \mathbb{R}^k, sometimes written in vector form:

y(x) = \begin{bmatrix} y_1(x) \\ \vdots \\ y_k(x) \end{bmatrix} \in \mathbb{R}^k

Goal of Control: Control objectives can be mathematically encoded through the use of outputs. Explicitly, the goal is to define a feedback control law u: D \to U such that any solution x(t) of the resulting closed loop system:

\dot{x} = f_{cl}(x) = f(x) + g(x) u(x)

drives the outputs to zero (drives y(x(t)) \to 0 as t \to \infty).

Example: One standard control objective is to drive the system to a desired state x_{\text{des}} \in \mathbb{R}^n. In this case:

y(x) = x - x_{\text{des}}

where y: D \to \mathbb{R}^n and thus k = n.

Virtual Constraints: We can achieve “output-based tracking” by defining a set of virtual constraints:

y(x, t) = y_a(x) - y_{\text{des}}(t, \alpha)

where y_a(x) is the actual output and y_{\text{des}}(t, \alpha) is the desired output, which is a function of time and some parameterization \alpha. We can remove the dependence on time using a parameterization of time \tau: D \to \mathbb{R}. This allows us to represent our virtual constraints as:

y(x) = y_a(x) - y_{\text{des}}(\tau(x), \alpha)

Therefore, achieving the objective y \to 0 as t \to \infty implies that y_a \to y_{\text{des}} as t \to \infty.

Bézier Polynomials: In the context of robotic systems, it is often useful to parameterize desired motions using a polynomial parameterization. One of the most popular choices is Bézier polynomials. A Bézier polynomial of degree M is defined as:

\begin{align} y_d(t, \alpha)_i &= \sum_{k=0}^{M} \frac{M!}{k!(M-k)!}\alpha_{k,i} t^k (1-t)^{M-k} \\ &= \sum_{k=0}^M \binom{M}{k} \alpha_{k,i} t^k (1-t)^{M-k} \end{align}

where \alpha_{k,i} are called the control points of the Bézier curve.

Phase-Based Parameterization of Gait

Parameterization of Time: In the case of robotic walking, we can parameterize time as the forward evolution of a walking robot through a step. Explicitly, this is done by defining a function \tau: D \to \mathbb{R}:

\tau(x) = \frac{\theta(q) - \theta^+}{\theta^- - \theta^+}

where \theta: D \to \mathbb{R} is a phase variable quantifying the forward progression of the robot (it must be monotonic), \theta^+ = \theta(x^+) is its value at the “beginning” of the step, and \theta^- = \theta(x^-) is its value at the “end” of the step. One common choice for \theta is the angle of the stance leg with respect to the vertical axis.

Feedback Linearization

We will begin by considering SISO (single input single output) systems. These are systems with k = m = 1 wherein the system takes the form:

\begin{align} \dot{x} &= f(x) + g(x) u \\ y &= h(x) \end{align}

with x \in D \subseteq \mathbb{R}^n, u \in U \subset \mathbb{R}, and h: D \to \mathbb{R}.

Definition: Feedback Linearizable

A control system is feedback linearizable (or input-state linearizable) if there exists a diffeomorphism z = T(x) and a feedback control law u: D \times U \to U (i.e., u(x,v)) such that the closed loop system:

\dot{x} = f(x) + g(x) u(x, v)

with v \in \mathbb{R} being a new control input, renders a linear relationship between the input and the state:

\dot{z} = Az + Bv

Since v is an auxiliary input, it can be chosen to stabilize the system dynamics (x) which are now linear.

However, sometimes (such as the case with tracking control) it is more beneficial to linearize the input-output map rather than the input-state map. This is called input-output linearization.

Definition: Input-Output Linearizable

A control system is input-output linearizable if there exists a feedback control law u: D \times U \to U (i.e., u(x,v)) such that the closed loop system:

\dot{x} = f(x) + g(x) u(x, v)

with v \in \mathbb{R} being a new control input, renders a linear relationship between the input and the output:

y^{(p)} = v

with p denoting the relative degree of the system.

In this case, since v is an auxiliary input, it can be chosen to stabilize the output dynamics (y) which are now linear.

Feedback Linearization Example: Let’s consider an inverted pendulum with torque actuation at the pivot point:

\dot{x} = \begin{bmatrix} x_2 \\ -\frac{g}{l}\sin(x_1) \end{bmatrix} + \begin{bmatrix} 0 \\ 1 \end{bmatrix}u

with x_1 = \theta being the angle (measured from the downward vertical) and x_2 = \dot{\theta} being the angular velocity.

from __future__ import print_function   
from scipy.integrate import odeint 

import time
import math
import numpy as np
import pylab as py

from matplotlib import animation, rc
from IPython.display import HTML
from matplotlib import pyplot as plt

x0 = [-np.pi/2.2, 0]    # initial conditions. 
# x[0] = angle 
# x[1] = angular velocity 

tfinal = 25.0       # Final time. Simulation time = 0 to tfinal.
Nt = 751

def get_params():
  L = 1.4              # length of pendulum 1 (in meter)
  g = 9.8               # gravitatioanl acceleration constant (m/s^2)
  return L, g

# Differential equations describing the system
def pendulum_system(x, t, controller):    
  L, g = get_params()
  dx = np.zeros(2)
  dx[0] = x[1]   # d(theta 1)
  u = controller(x, t)
  dx[1] = -(g/L)*np.sin(x[0]) + u
  return dx

def simulate_pendulum(x0, tVec, controller):
  sol = odeint(pendulum_system, x0, tVec, args=(controller,))
  sol0 = sol[:,0]     # theta_1 
  sol1 = sol[:,1]     # omega 1
  return sol0, sol1, tVec[2]-tVec[1]

def make_animation(dt, theta_vec):
  
  def get_mass_pos(theta):
    L, g = get_params()
    x = L*np.sin(theta)
    y = -L*np.cos(theta)
    return x, y
  
  def init():
    line1.set_data([], [])
    line4.set_data([], [])
    line5.set_data([], [])
    time_string.set_text('')
    return  line4, line5, line1, time_string
  
  def animate(i):
    trail1 = 6              
    line1.set_data(posx[i:max(1,i-trail1):-1], posy[i:max(1,i-trail1):-1])   
    line4.set_data([posx[i], 0], [posy[i],0])                
    line5.set_data([0, 0], [0, 0])
    time_string.set_text(time_template % (i*dt))
    return  line4, line5, line1, time_string
  
  posx, posy = get_mass_pos(theta_vec)
  
  fig = plt.figure()
  L, g = get_params()
  ax = plt.axes(xlim=(-L-0.5, L+0.5), ylim=(-2.5, 1.5))
  line1, = ax.plot([], [], 'o-', color='#d2eeff', markersize=12, markerfacecolor='#0077BE', lw=2, markevery=10000, markeredgecolor='k')   
  line4, = ax.plot([], [], color='k', linestyle='-', linewidth=2)
  line5, = ax.plot([], [], 'o', color='k', markersize=10)
  time_template = 'Time = %.1f s'
  time_string = ax.text(0.05, 0.9, '', transform=ax.transAxes)
  
  ax.get_xaxis().set_ticks([])    
  ax.get_yaxis().set_ticks([])    
  ax.set_aspect('equal', adjustable='box')
  
  anim = animation.FuncAnimation(fig, animate, init_func=init,
                 frames=len(theta_vec)-1, interval=1000*dt, blit=True)
  plt.close(fig)
  return anim

# Example controller function
def controller(x, t):
  return 0  # No control input

tVec = np.linspace(0, tfinal, Nt)
thetavec, thetadotvec, dt = simulate_pendulum(x0, tVec, controller)
anim = make_animation(dt, thetavec)
HTML(anim.to_jshtml())

Once Loop Reflect

Our control objective is to drive the pendulum to the upright position \theta = \pi, encoded by:

y = h(x) = x_1 - \pi

This system can be feedback linearized by choosing the control law u = \frac{g}{l}\sin(x_1) + v. Plugging this into the system dynamics yields the closed-loop system:

\dot{x} = \begin{bmatrix} x_2 \\ v \end{bmatrix} = \begin{bmatrix} 0 & 1 \\ 0 & 0 \end{bmatrix}x + \begin{bmatrix} 0 \\ 1 \end{bmatrix}v

To go further, we can select v = -k_p y - k_d \dot{y} to both stabilize the closed-loop system (for the shifted system \tilde{x} = x - (\pi, 0) such that the equilibrium point is at the origin) and to drive the output to zero. This results in the closed-loop system:

\dot{\tilde{x}} = \begin{bmatrix} 0 & 1 \\ -k_p & -k_d \end{bmatrix}\tilde{x}

Thus, we can stabilize our system by choosing k_p and k_d such that the eigenvalues have negative real parts.

We can also analyze the behavior of the output dynamics:

\begin{align} y &= x_1 - \pi \\ \dot{y} &= \dot{x}_1 = x_2 \\ \ddot{y} &= \dot{x}_2 = v \end{align}

Thus, the same conditions can be enforced on k_p and k_d to stabilize the second-order output dynamics:

\ddot{y} + k_d \dot{y} + k_p y = 0

Simulation Results for Input-Output Linearized System

def controller(x,t):
    L,g = get_params()
    Kp = 1.0
    Kd = 2.0
    u = g/L*np.sin(x[0]) - Kp*(x[0] - np.pi) - Kd*x[1]
    return u

# Animation
theta_vec, thetadot_vec, dt = simulate_pendulum(x0, tVec, controller)
anim = make_animation(dt, theta_vec)
plt.close()
HTML(anim.to_jshtml())

Once Loop Reflect

We can see that the system now resembles a linear system:

fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(8, 4))
L, g = get_params()
Kp = 1.0
Kd = 2.0

# Create a grid for the streamplot
x1 = np.linspace(-np.pi, 2*np.pi, 20)
x2 = np.linspace(-5, 5, 20)
X1, X2 = np.meshgrid(x1, x2)

# Left plot: without a controller
U_uncontrolled = X2
V_uncontrolled = -(g/L)*np.sin(X1)
ax1.streamplot(X1, X2, U_uncontrolled, V_uncontrolled, color='coral', density=1.5)
ax1.set_xlabel('$x_1$ (angle, rad)')
ax1.set_ylabel('$x_2$ (angular velocity, rad/s)')
ax1.set_title('Uncontrolled Pendulum')
ax1.grid(True, alpha=0.3)

# Right plot: with IO linearization
U_io = X2
V_io = - Kp*(X1 - np.pi) - Kd*X2
ax2.streamplot(X1, X2, U_io, V_io, color='steelblue', density=1.5)
ax2.set_xlabel('$x_1$ (angle, rad)')
ax2.set_ylabel('$x_2$ (angular velocity, rad/s)')
ax2.set_title('With Input-Output Linearization')
ax2.grid(True, alpha=0.3)

plt.tight_layout()
plt.show()

Comparison with PD Control (without input-output linearization)

For comparison, the following would be the system using a simple PD controller without feedback linearization:

def controller(x,t):
    L,g = get_params()
    Kp = 1.0
    Kd = 2.0
    u = - Kp*(x[0] - np.pi) - Kd*x[1]
    return u

# Animation
theta_vec, thetadot_vec, dt = simulate_pendulum(x0, tVec, controller)
anim = make_animation(dt, theta_vec)
plt.close()
HTML(anim.to_jshtml())

Once Loop Reflect

However, if we increase the gains K_p and K_d high enough the system can be stabilized without feedback linearization, albeit with a much smaller region of attraction. This is because the nonlinear term -(g/L)\sin(x_1) dominates the dynamics when x_1 is large, and thus the system behaves more like the uncontrolled system in those regions of the state space.

def controller(x,t):
    L,g = get_params()
    Kp = 10.0
    Kd = 2.0
    u = - Kp*(x[0] - np.pi) - Kd*x[1]
    return u

# Animation
theta_vec, thetadot_vec, dt = simulate_pendulum(x0, tVec, controller)
anim = make_animation(dt, theta_vec)
plt.close()
HTML(anim.to_jshtml())

Once Loop Reflect

fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(8, 4))
L, g = get_params()
Kp1 = 1.0
Kd1 = 2.0

Kp2 = 10.0
Kd2 = 2.0

# Create a grid for the streamplot
x1 = np.linspace(-np.pi, 2*np.pi, 20)
x2 = np.linspace(-5, 5, 20)
X1, X2 = np.meshgrid(x1, x2)

# Left plot: low-gain PD Control
U_uncontrolled = X2
V_uncontrolled = -(g/L)*np.sin(X1) - Kp1*(X1 - np.pi) - Kd1*X2
ax1.streamplot(X1, X2, U_uncontrolled, V_uncontrolled, color='coral', density=1.5)
ax1.set_xlabel('$x_1$ (angle, rad)')
ax1.set_ylabel('$x_2$ (angular velocity, rad/s)')
ax1.set_title('Low-Gain PD Control')
ax1.grid(True, alpha=0.3)

# Right plot: high-gain PD Control
U_io = X2
V_io = -(g/L)*np.sin(X1) - Kp2*(X1 - np.pi) - Kd2*X2
ax2.streamplot(X1, X2, U_io, V_io, color='steelblue', density=1.5)
ax2.set_xlabel('$x_1$ (angle, rad)')
ax2.set_ylabel('$x_2$ (angular velocity, rad/s)')
ax2.set_title('High-Gain PD Control')
ax2.grid(True, alpha=0.3)

plt.tight_layout()
plt.show()

Actuation Limits

Notice that for all systems, as soon as we impose actuation limits (i.e., u \in U with U being a bounded set), the system may no longer be globally feedback linearizable. In other words, there may not exist a globally valid feedback linearizing control law due to actuation limits. For the pendulum, this means that it cannot be stabilized without a swing-up controller that can drive the pendulum to a region of the state space where the controller can stabilize it.

def controller(x,t):
    L,g = get_params()
    Kp = 1.0
    Kd = 2.0
    u = g/L*np.sin(x[0]) - Kp*(x[0] - np.pi) - Kd*x[1]
    u = np.clip(u, -5, 5)  # Actuation limits
    return u

# Animation
theta_vec, thetadot_vec, dt = simulate_pendulum(x0, tVec, controller)
anim = make_animation(dt, theta_vec)
plt.close()
HTML(anim.to_jshtml())

Once Loop Reflect