SDC-Based Model Predictive Control: Enhancing Computational Feasibility for Safety-Critical Quadrotor Control

Saber Omidi S. Omidi is a Ph.D. Candidate with the Department of Mechanical Engineering, University of New Hampshire, Durham, NH 03824 USA (e-mail: saber.omidi@unh.edu).Manuscript received [Date]; revised [Date].
Abstract

Nonlinear Model Predictive Control (NMPC) is widely used for controlling high-speed robotic systems such as quadrotors. However, its significant computational demands often hinder real-time feasibility and reliability, particularly in environments requiring robust obstacle avoidance. This paper proposes a novel SDC-Based Model Predictive Control (MPC) framework, which preserves the high-precision performance of NMPC while substantially reducing computational complexity by over 30%. By reformulating the nonlinear quadrotor dynamics through the State-Dependent Coefficient (SDC) method, the original nonlinear program problem is transformed into a sequential quadratic optimization problem. The controller integrates an integral action to eliminate steady-state tracking errors and imposes constraints for safety-critical obstacle avoidance. Additionally, a disturbance estimator is incorporated to enhance robustness against external perturbations. Simulation results demonstrate that the SDC-Based MPC achieves comparable tracking accuracy to NMPC, with greater efficiency in terms of computation times, thereby improving its suitability for real-time applications. Theoretical analysis further establishes the stability and recursive feasibility of the proposed approach.

Index Terms:
Model Predictive Control (MPC); State-Dependent Coefficient (SDC); Nonlinear Control Systems; Quadrotors; Unmanned Aerial Vehicles (UAVs); Computational Efficiency; Safety-Critical Control; Robustness.

I Introduction

Unmanned aerial vehicles (UAVs), particularly quadrotors, have become increasingly important in both academic research and industry [1]. This growth is mainly due to their low cost, simple design, and ability to move with agility [2]. Quadrotors can take off and land vertically and hover in place, which makes them well-suited for tasks that traditional aircraft cannot handle [3]. They are now used in areas like photography, surveillance, search and rescue, precision agriculture, and package delivery [4]. As a result, improving flight control systems is essential to make the most of these vehicles in real-world situations [5]. Effectively controlling quadrotors presents a unique challenge due to their complex, nonlinear, and underactuated dynamics. Over the years, researchers and engineers have explored a wide range of flight control strategies to address these difficulties. Techniques, such as the Proportional-Integral-Derivative (PID) controller, remain popular for basic flight stabilization because they are straightforward to implement and tune [6]. Nevertheless, PID controllers depend on simplified, linearized models and can struggle to maintain performance during aggressive maneuvers or when faced with unexpected disturbances. To achieve better control in such situations, optimal controllers like the Linear Quadratic Regulator (LQR) are used to systematically compute control gains that minimize a predefined cost function. However, LQR and similar methods are typically effective only under nominal, linearized conditions [7]. Recognizing these shortcomings, nonlinear control approaches have gained traction for their ability to deal with the system’s true dynamics. For example, backstepping uses a stepwise design to guarantee stability, while sliding mode control offers strong resilience to model uncertainties and external disturbances by guiding the system along a specific trajectory [8, 9]. The aforementioned control methods are usually organized in a cascade control structure. This approach simplifies the control problem by separating it into a rapid inner loop for attitude stabilization and a slower outer loop for position tracking [10]. In recent years, learning-based controllers—particularly those leveraging reinforcement learning and neural networks—have shown promise by learning optimal control policies directly from real-world experience, eliminating the need for precise mathematical models [11]. Despite their potential, these methods often require large datasets for training and pose challenges in ensuring safe and reliable operation. Ensuring the safety of autonomous systems is a critical consideration in their synthesis and implementation, which has motivated the development of formally verifiable control methodologies. A prominent approach involves the use of Control Barrier Functions (CBFs), which are employed to render a defined safe set forward invariant. A CBF-based controller minimally modifies a nominal performance-oriented law; this intervention occurs only as the system state approaches the boundary of the safe set, thereby guaranteeing the system does not enter the unsafe region [12]. For stronger safety guarantees, albeit often at a higher computational cost, Reachability Analysis is utilized. This methodology entails the computation of the reachable set—the set of all states attainable by the system—and subsequently verifying that this set does not intersect with predefined unsafe regions [13]. Another influential paradigm leverages Lyapunov-based techniques, which are conventionally applied to stability analysis. In this context, a Lyapunov-like function, referred to as a barrier certificate, is constructed such that its value defines a safe sublevel set of the state space. The control law is then synthesized to ensure the time derivative of the barrier certificate is non-positive, thereby formally confining the system’s evolution to this region [14]. Model Predictive Control (MPC) is an advanced control framework that leverages an explicit system model to predict future behavior over a finite time horizon. At each sampling instant, MPC solves an online optimization problem to find an optimal sequence of control inputs, after which only the first input is applied in a strategy known as a receding horizon. Its primary strength lies in its ability to explicitly incorporate system constraints, making it highly suitable for safety-critical applications. The nature of the underlying optimization is dictated by the system model and constraints; linear systems with quadratic costs result in a solvable Quadratic Program (QP), while fully nonlinear systems require solving a computationally intensive Nonlinear Program (NLP). Building upon this core framework, several powerful variants have been developed. Nonlinear MPC (NMPC) directly uses a nonlinear model for high-fidelity prediction at the cost of significant online computation. To address this computational burden, Explicit MPC (EMPC) pre-computes the optimal control law offline as a piecewise affine function of the state, reducing the online effort to a simple function evaluation but is limited to lower-dimensional systems [15]. To handle real-world uncertainties, Robust MPC provides worst-case guarantees. A prominent implementation is Tube MPC, where a nominal trajectory is optimized and a separate feedback controller ensures the actual system state remains within a robust "tube" around this path, containing all possible deviations due to uncertainty. For probabilistic disturbances, Stochastic MPC (SMPC) optimizes expected performance while satisfying chance constraints [16]. Furthermore, modern Learning-based MPC integrates machine learning to learn system dynamics from data, adapting its predictive model over time [17]. Despite its versatility, MPC’s performance is fundamentally reliant on the accuracy of the predictive model, and its significant online computational demand can render it unsuitable for systems with very fast dynamics or limited onboard processing power. This paper addresses the critical challenge of developing safe, robust, and computationally tractable control systems for quadrotor unmanned aerial vehicles (UAVs). A comprehensive simulation framework was established in MATLAB, utilizing CasADi for numerical optimization, to systematically evaluate three distinct Model Predictive Control (MPC) architectures. The controllers were benchmarked against a challenging trajectory tracking mission requiring static obstacle avoidance amidst stochastic external disturbances. The evaluated strategies include a high-fidelity Nonlinear MPC (NMPC), which solves the full nonlinear optimal control problem and serves as a performance benchmark, alongside two novel quasi-linear formulations. The primary contribution of this research is the development of a robust, computationally efficient MPC that leverages a State-Dependent Coefficient (SDC) representation of the quadrotor’s nonlinear dynamics. This technique reformulates the system equations into a pseudo-linear form, which enables the nonlinear optimal control problem to be effectively approximated as a Quadratic Program (QP) at each sampling instant, significantly reducing computational complexity compared to the full NLP. Building upon this efficient foundation, we introduce a Robust SDC-based MPC (R-SDC-MPC). This advanced controller integrates a real-time disturbance observer that estimates lumped model uncertainties and external disturbances from the one-step state prediction error. This disturbance estimate is then used to actively compensate the predictive model, enhancing the controller’s resilience. Through quantitative analysis of tracking accuracy (MSE), computational load (CPU time), and overall cost, our work demonstrates that the proposed R-SDC-MPC offers a compelling trade-off, achieving robustness and performance comparable to NMPC while maintaining the computational tractability of a QP-based approach.

The remainder of this paper is organized as follows. Section II presents the nonlinear dynamical model of the quadrotor and details its reformulation into the State-Dependent Coefficient (SDC) representation. Section III develops the proposed control framework, starting with the NMPC benchmark, followed by the nominal SDC-MPC design, and culminating in the robust controller with disturbance estimation. In Section IV, a comprehensive simulation study is presented to validate the performance of the proposed methods, with a focus on tracking accuracy, computational efficiency, and safety. Finally, Section V concludes the paper with a summary of the key findings and outlines directions for future research.

II Quadrotor Dynamics and Problem Formulation

II-A Nonlinear Quadrotor Model

In this work, we adopt the six-degree-of-freedom (6-DoF) rigid-body model for the quadrotor, as formulated in [6, 18]. Figure 1 shows the 6-DoF quadrotor model, illustrating the inertial and body-fixed frames, the primary forces of thrust and gravity, and the control torques for roll, pitch, and yaw.

Refer to caption
Figure 1: Schematic of the 6-DoF quadrotor model, illustrating the inertial frame {e}\{e\} and the body-fixed frame {b}\{b\}. The quadrotor’s position is described by the vector pp. Control inputs consist of the total thrust force FTF_{T} along the ezbe_{z}^{b} axis and the body torques τϕ\tau_{\phi} (roll), τθ\tau_{\theta} (pitch), and τψ\tau_{\psi} (yaw), which induce the positive rotations ϕ(t)\phi(t), θ(t)\theta(t), and ψ(t)\psi(t), respectively. The gravitational force mgm_{g} acts in the negative ezee_{z}^{e} direction.

The 12-dimensional state vector 𝐱12\mathbf{x}\in\mathbb{R}^{12} is defined as:

𝐱=[𝐩T,𝐯T,𝜼T,𝝎T]T\mathbf{x}=[\mathbf{p}^{T},\mathbf{v}^{T},\bm{\eta}^{T},\bm{\omega}^{T}]^{T} (1)

where 𝐩=[x,y,z]T\mathbf{p}=[x,y,z]^{T} is the position in the inertial frame, 𝐯=[vx,vy,vz]T\mathbf{v}=[v_{x},v_{y},v_{z}]^{T} is the velocity in the inertial frame, 𝜼=[ϕ,θ,ψ]T\bm{\eta}=[\phi,\theta,\psi]^{T} represents the Euler angles (roll, pitch, yaw), and 𝝎=[p,q,r]T\bm{\omega}=[p,q,r]^{T} is the angular velocity vector in the body frame. The control input vector 𝐮4\mathbf{u}\in\mathbb{R}^{4} consists of the total thrust TT and the three body-axis torques 𝝉=[τϕ,τθ,τψ]T\bm{\tau}=[\tau_{\phi},\tau_{\theta},\tau_{\psi}]^{T}:

𝐮=[T,𝝉T]T\mathbf{u}=[T,\bm{\tau}^{T}]^{T} (2)

The translational and rotational dynamics are governed by the Newton-Euler equations:

𝐩˙\displaystyle\dot{\mathbf{p}} =𝐯\displaystyle=\mathbf{v} (3)
𝐯˙\displaystyle\dot{\mathbf{v}} =𝐠+1m𝐑(𝜼)𝐅T\displaystyle=\mathbf{g}+\frac{1}{m}\mathbf{R}(\bm{\eta})\mathbf{F}_{T} (4)
𝜼˙\displaystyle\dot{\bm{\eta}} =𝐖(𝜼)𝝎\displaystyle=\mathbf{W}(\bm{\eta})\bm{\omega} (5)
𝝎˙\displaystyle\dot{\bm{\omega}} =𝐉1(𝝉𝝎×(𝐉𝝎))\displaystyle=\mathbf{J}^{-1}(\bm{\tau}-\bm{\omega}\times(\mathbf{J}\bm{\omega})) (6)

where mm is the total mass, 𝐉3×3\mathbf{J}\in\mathbb{R}^{3\times 3} is the inertia matrix, and 𝐠=[0,0,g]T\mathbf{g}=[0,0,-g]^{T} is the acceleration due to gravity. The thrust vector in the body frame is 𝐅T=[0,0,T]T\mathbf{F}_{T}=[0,0,T]^{T}. The matrix 𝐑(𝜼)SO(3)\mathbf{R}(\bm{\eta})\in SO(3) is the rotation matrix mapping vectors from the body frame to the inertial frame, and 𝐖(𝜼)\mathbf{W}(\bm{\eta}) is the kinematic transformation matrix from body rates to Euler angle rates.

II-B SDC Representation

A general continuous-time nonlinear dynamical system can be expressed as:

𝐱˙(t)=𝐟(𝐱(t),𝐮(t)),𝐱(0)=𝐱0,\dot{\mathbf{x}}(t)=\mathbf{f}(\mathbf{x}(t),\mathbf{u}(t)),\qquad\mathbf{x}(0)=\mathbf{x}_{0}, (5)

where 𝐱(t)n\mathbf{x}(t)\in\mathbb{R}^{n} and 𝐮(t)m\mathbf{u}(t)\in\mathbb{R}^{m} are the system state vector and the control input vector, respectively; 𝐱0n\mathbf{x}_{0}\in\mathbb{R}^{n} is the initial-valued vector of the states. For the purpose of control design, a possible way to express the dynamical system in (5) is as:

𝐱˙(t)\displaystyle\dot{\mathbf{x}}(t) =𝐚(𝐱(t))+𝐁(𝐱(t))𝐮(t),\displaystyle=\mathbf{a}(\mathbf{x}(t))+\mathbf{B}(\mathbf{x}(t))\mathbf{u}(t), (6a)
𝐲(t)\displaystyle\mathbf{y}(t) =𝐜(𝐱(t))\displaystyle=\mathbf{c}(\mathbf{x}(t)) (6b)

where 𝐲(t)q\mathbf{y}(t)\in\mathbb{R}^{q} is the system output vector; 𝐚(𝐱(t)):nn\mathbf{a}(\mathbf{x}(t)):\mathbb{R}^{n}\to\mathbb{R}^{n}, 𝐁(𝐱(t)):nn×m\mathbf{B}(\mathbf{x}(t)):\mathbb{R}^{n}\to\mathbb{R}^{n\times m}, and 𝐜(𝐱(t)):nq\mathbf{c}(\mathbf{x}(t)):\mathbb{R}^{n}\to\mathbb{R}^{q} are assumed to be smooth. The matrix-valued functions in (6) are assumed to be smooth and satisfy 𝐚(0)=𝐜(0)=0\mathbf{a}(0)=\mathbf{c}(0)=0. Therefore, these functions can be written in their pseudo-linearized or SDC forms:

𝐚(𝐱(t))\displaystyle\mathbf{a}(\mathbf{x}(t)) =𝐀(𝐱(t))𝐱(t),\displaystyle=\mathbf{A}(\mathbf{x}(t))\mathbf{x}(t), (8a)
𝐜(𝐱(t))\displaystyle\mathbf{c}(\mathbf{x}(t)) =𝐂(𝐱(t))𝐱(t)\displaystyle=\mathbf{C}(\mathbf{x}(t))\mathbf{x}(t) (8b)

where 𝐀:nn×n\mathbf{A}:\mathbb{R}^{n}\to\mathbb{R}^{n\times n}, and 𝐂:np×n\mathbf{C}:\mathbb{R}^{n}\to\mathbb{R}^{p\times n} are two matrix-valued functions. It is important to note that there are infinite ways to write pseudo-linearized non-scalar systems which allow an additional level of flexibility to select [25]. The specific factorization used in this work is chosen for its straightforward derivation from the Newton-Euler equations and its ability to isolate the control inputs in the 𝐁(𝐱)\mathbf{B}(\mathbf{x}) matrix, which is advantageous for control design. Furthermore, this structure has been shown to maintain the necessary stabilizability and detectability properties for a wide range of operating conditions. Now that the matrices are expressed in SDC form, the system can be written as follows:

𝐱˙(t)\displaystyle\dot{\mathbf{x}}(t) =𝐀(𝐱(t))𝐱(t)+𝐁(𝐱(t))𝐮(t),\displaystyle=\mathbf{A}(\mathbf{x}(t))\mathbf{x}(t)+\mathbf{B}(\mathbf{x}(t))\mathbf{u}(t), (9a)
𝐲(t)\displaystyle\mathbf{y}(t) =𝐂(𝐱(t))𝐱(t).\displaystyle=\mathbf{C}(\mathbf{x}(t))\mathbf{x}(t). (9b)

we consider a full state-feedback control problem, where the entire state vector 𝐱(t)\mathbf{x}(t) is assumed to be available to the controller at each time instant. Consequently, the system output is the state itself, rendering a separate output equation unnecessary for the control formulation. The dynamics of the quadrotor are therefore represented by the following SDC-based pseudo-linear model:

𝐁(𝐱)=[𝟎3×1𝟎3×31m𝐑(𝜼)𝐞3𝟎3×3𝟎3×1𝟎3×3𝟎3×1𝐉1]\mathbf{B}(\mathbf{x})=\begin{bmatrix}\mathbf{0}_{3\times 1}&\mathbf{0}_{3\times 3}\\ \frac{1}{m}\mathbf{R}(\bm{\eta})\mathbf{e}_{3}&\mathbf{0}_{3\times 3}\\ \mathbf{0}_{3\times 1}&\mathbf{0}_{3\times 3}\\ \mathbf{0}_{3\times 1}&\mathbf{J}^{-1}\end{bmatrix} (7)

where 𝐞3=[0,0,1]T\mathbf{e}_{3}=[0,0,1]^{T}. The system matrix 𝐀(𝐱)12×12\mathbf{A}(\mathbf{x})\in\mathbb{R}^{12\times 12} is:

𝐀(𝐱)=[𝟎3×3𝐈3×3𝟎3×3𝟎3×3𝟎3×3𝟎3×3𝐀23(𝜼,T)𝟎3×3𝟎3×3𝟎3×3𝐀33(𝜼,𝝎)𝐖(𝜼)𝟎3×3𝟎3×3𝟎3×3𝐀44(𝝎)]\mathbf{A}(\mathbf{x})=\begin{bmatrix}\mathbf{0}_{3\times 3}&\mathbf{I}_{3\times 3}&\mathbf{0}_{3\times 3}&\mathbf{0}_{3\times 3}\\ \mathbf{0}_{3\times 3}&\mathbf{0}_{3\times 3}&\mathbf{A}_{23}(\bm{\eta},T)&\mathbf{0}_{3\times 3}\\ \mathbf{0}_{3\times 3}&\mathbf{0}_{3\times 3}&\mathbf{A}_{33}(\bm{\eta},\bm{\omega})&\mathbf{W}(\bm{\eta})\\ \mathbf{0}_{3\times 3}&\mathbf{0}_{3\times 3}&\mathbf{0}_{3\times 3}&\mathbf{A}_{44}(\bm{\omega})\end{bmatrix} (8)

where the non-zero, state-dependent blocks are defined as:

𝐀23(𝜼,T)\displaystyle\mathbf{A}_{23}(\bm{\eta},T) =Tm𝐑(𝜼)𝜼,\displaystyle=\frac{T}{m}\frac{\partial\mathbf{R}(\bm{\eta})}{\partial\bm{\eta}}, (9)
𝐀33(𝜼,𝝎)\displaystyle\mathbf{A}_{33}(\bm{\eta},\bm{\omega}) =(𝐖(𝜼)𝝎)𝜼,\displaystyle=\frac{\partial(\mathbf{W}(\bm{\eta})\bm{\omega})}{\partial\bm{\eta}}, (10)
𝐀44(𝝎)\displaystyle\mathbf{A}_{44}(\bm{\omega}) =𝐉1(𝐒(𝐉𝝎)𝐒(𝝎)𝐉).\displaystyle=\mathbf{J}^{-1}(\mathbf{S}(\mathbf{J}\bm{\omega})-\mathbf{S}(\bm{\omega})\mathbf{J}). (11)

Here 𝐈3×3\mathbf{I}_{3\times 3} is the identity matrix and 𝟎3×3\mathbf{0}_{3\times 3} is zero matrix. 𝐑(𝜼)𝜼\frac{\partial\mathbf{R}(\bm{\eta})}{\partial\bm{\eta}} is the partial derivative of the rotation matrix with respect to the Euler angles. (𝐖(𝜼)𝝎)𝜼\frac{\partial(\mathbf{W}(\bm{\eta})\bm{\omega})}{\partial\bm{\eta}} is the partial derivative of the kinematic transformation. 𝐒()\mathbf{S}(\cdot) is the skew-symmetric matrix operator. For a comprehensive review of related quasi-linear control methodologies, including the closely allied State-Dependent Riccati Equation (SDRE) framework, the reader is referred to the works in [19, 20, 21, 22].

II-C SDC-Based Optimal Control Problem Formulation

The main idea of the proposed SDC-based MPC framework is an online optimization problem solved at each sampling instant tt. This problem aims to compute an optimal control input sequence over a finite prediction horizon, TpT_{p}, that steers the quadrotor along a time-varying reference trajectory, 𝐱ref(τ)\mathbf{x}_{\text{ref}}(\tau). Given the current state 𝐱(t)\mathbf{x}(t), the finite-horizon optimal control problem is formulated as follows:

min𝐮()J(.)\displaystyle\min_{\mathbf{u}(\cdot)}J(.) =tt+Tp(𝐱(τ)𝐱ref(τ)𝐐2+𝐮(τ)𝐑2)𝑑τ\displaystyle\!=\!\int_{t}^{t+T_{p}}\!\left(\!\|\mathbf{x}(\tau)\!-\!\mathbf{x}_{\text{ref}}(\tau)\|^{2}_{\mathbf{Q}}\!+\!\|\mathbf{u}(\tau)\|^{2}_{\mathbf{R}}\right)\!d\tau
+𝐱(t+Tp)𝐱ref(t+Tp)𝐏(𝐱)2\displaystyle+\|\mathbf{x}(t+T_{p})-\mathbf{x}_{\text{ref}}(t+T_{p})\|^{2}_{\mathbf{P}(\mathbf{x})} (12a)
s.t.𝐱˙(τ)\displaystyle\text{s.t.}\quad\dot{\mathbf{x}}(\tau) =𝐀(𝐱(τ))𝐱(τ)+𝐁(𝐱(τ))𝐮(τ),\displaystyle=\mathbf{A}(\mathbf{x}(\tau))\mathbf{x}(\tau)+\mathbf{B}(\mathbf{x}(\tau))\mathbf{u}(\tau), (12b)
𝐱(τ)\displaystyle\mathbf{x}(\tau) 𝕏,\displaystyle\in\mathbb{X}, (12c)
𝐮(τ)\displaystyle\mathbf{u}(\tau) 𝕌,\displaystyle\in\mathbb{U}, (12d)
𝐱(t+Tp)\displaystyle\mathbf{x}(t+T_{p}) 𝕏T(t),τ[t,t+Tp]\displaystyle\in\mathbb{X}_{T}(t),\quad\forall\tau\in[t,t+T_{p}] (12e)

where the optimization is performed over the control input sequence 𝐮()\mathbf{u}(\cdot) for all τ[t,t+Tp]\tau\in[t,t+T_{p}]. The objective function in (12a) is composed of a running cost and a terminal cost. The integral term represents the running cost, which penalizes the deviation of the predicted state trajectory 𝐱(τ)\mathbf{x}(\tau) from the reference 𝐱ref(τ)\mathbf{x}_{\text{ref}}(\tau) as well as the magnitude of the control input 𝐮(τ)\mathbf{u}(\tau). The semi-positive-definite and positive-definite weighting matrices 𝐐n×n\mathbf{Q}\in\mathbb{R}^{n\times n} and 𝐑m×m\mathbf{R}\in\mathbb{R}^{m\times m}, respecetivly, allow for tuning the trade-off between tracking accuracy and control effort. The terminal cost, weighted by the matrix 𝐏(𝐱)\mathbf{P}(\mathbf{x}), penalizes the final state error and is crucial for proving closed-loop stability. The optimization is subject to several constraints that define the system’s behavior and limitations.

  • The differential equation (12b) enforces the system’s SDC dynamics throughout the prediction horizon.

  • Constraints (12c) and (12d) confine the state and control inputs to their respective admissible sets. The input constraint set, 𝕌\mathbb{U}, directly represents the physical limitations of the actuators (e.g., maximum thrust and torques). The state constraint set, 𝕏\mathbb{X}, is paramount for defining the safe operating envelope. For safety-critical tasks such as obstacle avoidance, this set is explicitly structured to exclude unsafe regions of the state-space. By imposing such constraints, the MPC controller is forced to find a control sequence that guarantees the predicted trajectory remains entirely within the safe region.

  • Finally, the terminal inequality constraint(12e) compels the predicted state to lie within a terminal set 𝕏T(t)\mathbb{X}_{T}(t) at the end of the horizon. This set is typically defined relative to the reference state 𝐱ref(t+Tp)\mathbf{x}_{\text{ref}}(t+T_{p}) and is essential for guaranteeing the recursive feasibility and stability of the closed-loop system.

Finally, a primary objective of this work is to improve the computational efficiency of the online optimization, ensuring the controller is feasible for real-time implementation without sacrificing performance or violating system constraints.

III Design of the SDC-Based MPC Framework

III-A NMPC Benchmark

The NMPC benchmark is based on the quasi-infinite horizon scheme proposed in [23]. At each sampling instant tt, the controller solves the following finite-horizon optimal control problem for the current state x(t)x(t):

min𝐮()J(.)\displaystyle\min_{\mathbf{u}(\cdot)}J(.) =tt+Tp(𝐱(τ)𝐱ref(τ)𝐐2+𝐮(τ)𝐑2)𝑑τ\displaystyle\!=\!\int_{t}^{t+T_{p}}\!\left(\!\|\mathbf{x}(\tau)\!-\!\mathbf{x}_{\text{ref}}(\tau)\|^{2}_{\mathbf{Q}}\!+\!\|\mathbf{u}(\tau)\|^{2}_{\mathbf{R}}\right)\!d\tau
+𝐱(t+Tp)𝐱ref(t+Tp)𝐐t2\displaystyle+\|\mathbf{x}(t+T_{p})-\mathbf{x}_{\text{ref}}(t+T_{p})\|^{2}_{\mathbf{Q}_{t}} (13a)
s.t.𝐱˙(τ)\displaystyle\text{s.t.}\quad\dot{\mathbf{x}}(\tau) =𝐟(𝐱(τ),𝐮(τ)),\displaystyle=\mathbf{f}(\mathbf{x}(\tau),\mathbf{u}(\tau)), (13b)
𝐱(τ)\displaystyle\mathbf{x}(\tau) 𝕏,\displaystyle\in\mathbb{X}, (13c)
𝐮(τ)\displaystyle\mathbf{u}(\tau) 𝕌,\displaystyle\in\mathbb{U}, (13d)
𝐱(t+Tp)\displaystyle\mathbf{x}(t+T_{p}) Ω,τ[t,t+Tp]\displaystyle\in\Omega,\quad\forall\tau\in[t,t+T_{p}] (13e)

where the terminal penalty matrix 𝐐tn×n\mathbf{Q}_{t}\in\mathbb{R}^{n\times n} is a positive-definite, symmetric matrix chosen to guarantee stability [23]. It is computed off-line as the unique solution to the Lyapunov equation:

(𝐀K+κ𝐈)T𝐒+𝐒(𝐀K+κ𝐈)=(𝐐+𝐊T𝐑𝐊)\displaystyle(\mathbf{A}_{K}+\kappa\mathbf{I})^{T}\mathbf{S}\!+\!\mathbf{S}(\mathbf{A}_{K}+\kappa\mathbf{I})\!=\!-(\mathbf{Q}+\mathbf{K}^{T}\mathbf{R}\mathbf{K}) (14)

where 𝐀K=𝐀+𝐁𝐊\mathbf{A}_{K}=\mathbf{A}+\mathbf{BK} is the closed-loop matrix of the system’s Jacobian linearization, and κ0\kappa\geq 0 is a tuning parameter [23]. 𝐊\mathbf{K} is a locally stabilizing state feedback gain determined linearized system, typically by solving the LQR problem [24]. The terminal region Ω\Omega is an invariant set for the nonlinear system under a local linear state feedback law 𝐮=𝐊𝐱(t)\mathbf{u}=\mathbf{K}\mathbf{x}(t) [23]. It is defined as a level set of the Lyapunov function, Ω:={𝐱n𝐱T𝐒𝐱α}\Omega:=\{\mathbf{x}\in\mathbb{R}^{n}\mid\mathbf{x}^{T}\mathbf{S}\mathbf{x}\leq\alpha\}, where α>0\alpha>0 is chosen to ensure that for any state within Ω\Omega, the corresponding linear control input satisfies the constraint 𝐊𝐱(t)𝕌\mathbf{K}\mathbf{x}(t)\in\mathbb{U} [23]. The explicit use of a general nonlinear dynamics model, fundamentally casts the optimal control problem as a non-convex NLP. Solving such an NLP is computationally demanding and challenging when the predicted trajectory interacts with the boundaries of the constraint sets. This is particularly for safety-critical applications where complex state constraints, 𝐱(τ)Ω\mathbf{x}(\tau)\in\Omega, are activated to enforce collision avoidance by proscribing certain regions of the state-space. The computational burden can render real-time implementation infeasible. Moreover, the reliance on a local linearization to prove stability introduces significant conservatism for systems with prominent nonlinearities. The terminal region where the linear approximation is valid can be impractically small [23]. This necessitates a longer prediction horizon, TpT_{p}, to ensure the optimization problem remains feasible, which in turn increases the computational burden and can severely restrict the controller’s effective region of attraction.

III-B Nominal SDC-Based MPC Formulation

To address the computational challenges inherent to the NMPC benchmark, the nominal SDC-MPC scheme is proposed. The primary advantage of this approach is its ability to reformulate the computationally intensive NLP into a sequence of QP with utilization of the SDC representation into the dynamic constraint. By parameterizing the system matrices with the most recent state measurement, the dynamic constraint for the optimal control problem is rendered as a affine system. This model is unique to each time step, as the matrices 𝐀(𝐱)\mathbf{A}(\mathbf{x}) and 𝐁(𝐱)\mathbf{B}(\mathbf{x}) are re-calculated based on the new state measurement. This transformation of the dynamics into a set of linear equality constraints is critical. When combined with a quadratic objective function and convex state and input constraints, it ensures that the overall optimal control problem becomes a convex QP. This QP can be solved efficiently and reliably, avoiding the complexities and computational burden of the original NLP. Another key consequence of SDC formulation is that the terminal penalty matrix, 𝐏(𝐱)\mathbf{P}(\mathbf{x}), also becomes state-dependent. Consequently, a unique terminal penalty matrix, 𝐏(𝐱)\mathbf{P}(\mathbf{x}), must be computed at each time step based on the current state measurement to ensure stability of the system. With considering this, at each state 𝐱(t)\mathbf{x}(t), the State-Dependent Lyapunov Equation (SDLE) is formulated to find the positive-definite matrix 𝐏(𝐱)\mathbf{P}(\mathbf{x}) that satisfies:

(𝐀\displaystyle\big(\mathbf{A} (𝐱)𝐁(𝐱)𝐊(𝐱))T𝐏(𝐱)\displaystyle(\mathbf{x})-\mathbf{B}(\mathbf{x})\mathbf{K}(\mathbf{x})\big)^{T}\mathbf{P}(\mathbf{x})
+𝐏(𝐱)(𝐀\displaystyle+\mathbf{P}(\mathbf{x})\big(\mathbf{A} (𝐱)𝐁(𝐱)𝐊(𝐱))\displaystyle(\mathbf{x})-\mathbf{B}(\mathbf{x})\mathbf{K}(\mathbf{x})\big)
+(𝐐\displaystyle+\big(\mathbf{Q} (𝐱)+𝐊(𝐱)T𝐑𝐊(𝐱))=𝟎\displaystyle(\mathbf{x})+\mathbf{K}(\mathbf{x})^{T}\mathbf{R}\mathbf{K}(\mathbf{x})\big)=\mathbf{0} (15)

where the existence of a positive-definite solution 𝐏(𝐱)\mathbf{P}(\mathbf{x}) for all 𝐱\mathbf{x} in a given domain is a sufficient condition for local stability. The existence of a unique, positive-definite solution 𝐏(𝐱)\mathbf{P}(\mathbf{x}) to SDLE is contingent upon the stability of the closed-loop system, which is directly determined by the feedback gain matrix 𝐊(𝐱)\mathbf{K}(\mathbf{x}). This gain is the sub-optimal state feedback controller, computed as

𝐊(𝐱)=𝐑1𝐁(𝐱)T𝐒(𝐱),\displaystyle\mathbf{K}(\mathbf{x})=\mathbf{R}^{-1}\mathbf{B}(\mathbf{x})^{T}\mathbf{S}(\mathbf{x}), (16)

and is designed to stabilize the system at the state 𝐱\mathbf{x}. The gain itself is derived from the solution 𝐒(𝐱)\mathbf{S}(\mathbf{x}) of the State-Dependent Riccati Equation (SDRE), given by:

𝐀(𝐱)T𝐒(𝐱)+𝐒(𝐱)𝐀(𝐱)\displaystyle\mathbf{A}(\mathbf{x})^{T}\mathbf{S}(\mathbf{x})+\mathbf{S}(\mathbf{x})\mathbf{A}(\mathbf{x})
\displaystyle- 𝐒(𝐱)𝐁(𝐱)𝐑1𝐁(𝐱)T𝐒(𝐱)+𝐐(𝐱)=𝟎\displaystyle\mathbf{S}(\mathbf{x})\mathbf{B}(\mathbf{x})\mathbf{R}^{-1}\mathbf{B}(\mathbf{x})^{T}\mathbf{S}(\mathbf{x})+\mathbf{Q}(\mathbf{x})=\mathbf{0} (17)

Therefore, the existence of a positive-definite solution to the Lyapunov equation depends on finding a stabilizing solution to the Riccati equation. While solving the SDRE and SDLE at each time step introduces a computational overhead, this process is typically very fast for systems of this dimension and its cost is included in the total per-iteration CPU times reported in Section IV. For the quadrotor model, this overhead was found to be a small fraction of the total optimization time. Under the following theorem, we can guarantee that the SDRE solution is feasible [19].

Theorem 1 (Existence of a Stabilizing SDRE Solution).

For a nonlinear system in SDC form, a unique, positive-definite stabilizing solution 𝐒(𝐱)\mathbf{S}(\mathbf{x}) to the State-Dependent Riccati Equation (17) exists at a state 𝐱\mathbf{x} if the following state-dependent pairs are pointwise stabilizable and detectable, respectively [19, 25]:

  1. 1.

    The pair (𝐀(𝐱),𝐁(𝐱))(\mathbf{A}(\mathbf{x}),\mathbf{B}(\mathbf{x})) is stabilizable.

  2. 2.

    The pair (𝐀(𝐱),𝐐(𝐱)1/2)(\mathbf{A}(\mathbf{x}),\mathbf{Q}(\mathbf{x})^{1/2}) is detectable.

The fundamental argument for establishing recursive feasibility is conceptually consistent for SDC based MPC and linear MPC schemes. This argument relies on constructing a feasible candidate solution for the current time step using the optimal solution from the preceding step. Therefore, the proof of recursive feasibility is contingent upon a key assumption regarding the system’s local dynamics and a constructive lemma that guarantees the existence of the necessary terminal ingredients. Nevertheless, the theoretical underpinnings of the proof and its practical implications diverge significantly due to the state-dependent, nonlinear nature of the SDC framework.

Assumption 1 (Stabilizability of the Linearization).

The Jacobian linearization of the system at the origin, represented by the pair of matrices (𝐀,𝐁)(\mathbf{A},\mathbf{B}), is stabilizable [23].

This assumption is fundamental, as it ensures that the system can be locally stabilized, which is a prerequisite for constructing a valid terminal set for the MPC. This leads to the following critical lemma.

Lemma 2 (Existence of Terminal Ingredients).

If Assumption 1 holds, then a stabilizing linear feedback gain 𝐊\mathbf{K}, a positive-definite terminal penalty matrix 𝐏\mathbf{P}, and a terminal region 𝕏T\mathbb{X}_{T} can be constructed. This terminal region 𝕏T\mathbb{X}_{T} is a positively invariant set for the full nonlinear system under the local linear feedback law 𝐮=𝐊𝐱\mathbf{u}=\mathbf{Kx} [23, 26].

With this foundation, we can formally state and prove the recursive feasibility of the controller.

Theorem 3 (Recursive Feasibility).

If Assumption 1 holds and the optimization problem is feasible at the initial time t=0t=0, then a feasible solution is guaranteed to exist for all subsequent time steps t>0t>0.

Proof.

The proof is by induction. Assume a feasible optimal solution exists at time tt, yielding the control sequence 𝐮(;t)\mathbf{u}^{*}(\cdot;t) and state trajectory 𝐱(;t)\mathbf{x}^{*}(\cdot;t), which satisfies the terminal constraint 𝐱(t+Tp;t)𝕏T\mathbf{x}^{*}(t+T_{p};t)\in\mathbb{X}_{T}. At the next sampling instant, t+δt+\delta, we construct a candidate solution 𝐮^()\hat{\mathbf{u}}(\cdot) for the new horizon by using the tail of the previous solution and appending the terminal control law:

𝐮^(τ)={𝐮(τ;t)for τ[t+δ,t+Tp]𝐊𝐱^(τ)for τ(t+Tp,t+δ+Tp]\hat{\mathbf{u}}(\tau)=\begin{cases}\mathbf{u}^{*}(\tau;t)&\text{for }\tau\in[t+\delta,t+T_{p}]\\ \mathbf{K}\hat{\mathbf{x}}(\tau)&\text{for }\tau\in(t+T_{p},t+\delta+T_{p}]\end{cases} (18)

where the existence of the stabilizing gain 𝐊\mathbf{K} and the terminal set 𝕏T\mathbb{X}_{T} is guaranteed by Lemma 2. This candidate is feasible because the terminal state 𝐱(t+Tp;t)\mathbf{x}^{*}(t+T_{p};t) lies in 𝕏T\mathbb{X}_{T}, which is a positively invariant set under the control law 𝐮=𝐊𝐱\mathbf{u}=\mathbf{Kx}. Therefore, the new terminal state at t+δ+Tpt+\delta+T_{p} is also guaranteed to be in 𝕏T\mathbb{X}_{T}, satisfying the terminal constraint. Since a feasible solution can be constructed, the optimization problem remains feasible for all subsequent steps [23, 27]. ∎

The determination of a locally stabilizing feedback gain, derived from the Jacobian linearization of the system at its equilibrium, constitutes the foundational prerequisite for the recursive feasibility proof. This concept is illustrated conceptually in Figure 2. This off-line synthesis of a stabilizing terminal controller and an associated invariant terminal set enables the construction of a feasible candidate solution at each subsequent sampling instant, thereby satisfying the conditions required to guarantee recursive feasibility.

Refer to caption
Figure 2: Conceptual schematic of recursive feasibility with a time-varying terminal set. The solid blue line shows the optimal trajectory 𝐱(;t)\mathbf{x}^{*}(\cdot;t) computed at time tt, which terminates in the invariant set 𝕏T(t)\mathbb{X}_{T}(t). At time t+δt+\delta, the dashed red line shows a feasible candidate trajectory. The solid green line represents the new optimal trajectory, 𝐱(;t+δ)\mathbf{x}^{*}(\cdot;t+\delta), which terminates in the updated terminal set 𝕏T(t+δ)\mathbb{X}_{T}(t+\delta).

III-C Robust SDC-Based MPC with Disturbance Estimation

To enhance the controller’s resilience against inevitable model-plant mismatch and external disturbances, we extend the nominal SDC-based MPC framework with a robust formulation. This is achieved by designing a disturbance observer to actively estimate and compensate for uncertainties. The true system dynamics are modeled as being affected by a lumped disturbance term, 𝐝(t)\mathbf{d}(t), which represents the combined effects of parametric uncertainties, unmodeled dynamics, and exogenous forces:

𝐱˙(t)=𝐀(𝐱)𝐱(t)+𝐁(𝐱)𝐮(t)+𝐄d𝐝(t)\dot{\mathbf{x}}(t)=\mathbf{A}(\mathbf{x})\mathbf{x}(t)+\mathbf{B}(\mathbf{x})\mathbf{u}(t)+\mathbf{E}_{d}\mathbf{d}(t) (19)

where 𝐝(t)n\mathbf{d}(t)\in\mathbb{R}^{n} is the vector of unknown disturbance forces and torques, and 𝐄d\mathbf{E}_{d} is the disturbance input matrix. To ensure the tractability of the estimation problem, a standard assumption is made on the nature of this disturbance.

Assumption 2 (Disturbance Characteristics).

The lumped disturbance 𝐝(t)\mathbf{d}(t) is assumed to be bounded and slowly time-varying, such that its rate of change is negligible over a single sampling interval, i.e., 𝐝˙(t)𝟎\dot{\mathbf{d}}(t)\approx\mathbf{0} for t[tk,tk+1)t\in[t_{k},t_{k+1}) [28].

Under this assumption, a discrete-time state-observer is designed to estimate the disturbance in real-time. The estimation process at each sampling instant tkt_{k} involves computing the one-step state prediction error, 𝐞k\mathbf{e}_{k}, which is the discrepancy between the measured state, 𝐱meas(tk)\mathbf{x}_{\text{meas}}(t_{k}), and the state predicted from the previous control action using the nominal model:

𝐞k=𝐱meas(tk)𝐱pred(tk)\mathbf{e}_{k}=\mathbf{x}_{\text{meas}}(t_{k})-\mathbf{x}_{\text{pred}}(t_{k}) (20)

This prediction error is used to infer the disturbance that acted on the system, yielding a raw measurement, 𝐝^meas(tk)\hat{\mathbf{d}}_{\text{meas}}(t_{k}). To reduce sensitivity to measurement noise, this estimate is passed through a first-order low-pass filter:

𝐝^(tk)=α𝐝^(tk1)+(1α)𝐝^meas(tk)\hat{\mathbf{d}}(t_{k})=\alpha\hat{\mathbf{d}}(t_{k-1})+(1-\alpha)\hat{\mathbf{d}}_{\text{meas}}(t_{k}) (21)

The observer gain α[0,1)\alpha\in[0,1) determines the trade-off between the observer’s convergence rate and its noise sensitivity. The gain was selected empirically as 0.9 to provide a balance between responsiveness to changing disturbances and rejection of high-frequency measurement noise.

Lemma 4 (Observer Error Boundedness).

If Assumption 2 holds and the system states remain bounded, the estimation error of the observer (21), defined as 𝐝~(tk)=𝐝(tk)𝐝^(tk)\tilde{\mathbf{d}}(t_{k})=\mathbf{d}(t_{k})-\hat{\mathbf{d}}(t_{k}), is bounded for a suitable choice of the observer gain α\alpha [29].

This disturbance estimate, 𝐝^(tk)\hat{\mathbf{d}}(t_{k}), is subsequently injected into the MPC’s internal prediction model as a corrective offset term. For the optimization performed at time tkt_{k}, the pseudo-linear prediction dynamics are given by:

𝐱i+1|k=𝐀k𝐱i|k+𝐁k𝐮i|k+𝐂k+𝐄d𝐝^(tk)\mathbf{x}_{i+1|k}=\mathbf{A}_{k}\mathbf{x}_{i|k}+\mathbf{B}_{k}\mathbf{u}_{i|k}+\mathbf{C}_{k}+\mathbf{E}_{d}\hat{\mathbf{d}}(t_{k}) (22)
Remark 1 (On the SDC Prediction Model).

The prediction model (22) is an affine system for the duration of a single optimization instance. The matrices 𝐀k\mathbf{A}_{k}, 𝐁k\mathbf{B}_{k}, and 𝐂k\mathbf{C}_{k} are constant within this context. They are obtained by evaluating the state-dependent matrices 𝐀(𝐱)\mathbf{A}(\mathbf{x}) and 𝐁(𝐱)\mathbf{B}(\mathbf{x}) at a specific operating point, such as the current state 𝐱(tk)\mathbf{x}(t_{k}). This approach differs fundamentally from Jacobian linearization; it is an exact algebraic restructuring that fully captures the nonlinearities at that specific point, rather than an approximation that discards higher-order terms.

By incorporating the disturbance estimate, the controller actively accommodates for the disturbance over its prediction horizon. The stability of the composite observer-controller system can be analyzed by considering the separation principle, where the stability of the nominal MPC and the boundedness of the observer error together ensure the ultimate boundedness of the closed-loop system states [30]. The formal verification of recursive feasibility and stability is essential for this robust control architecture. In the presence of persistent disturbances, the objective shifts from proving nominal asymptotic stability to ensuring that the optimization remains feasible despite the model-plant mismatch and that the system state is ultimately bounded within a small neighborhood of the reference trajectory. These properties are typically established by constructing a robustly invariant terminal set and employing an Input-to-State Stability (ISS) framework in the Lyapunov analysis [30, 26]. The implementation of the proposed robust SDC-based MPC is formally outlined in Algorithm 1. The algorithm proceeds iteratively at each sampling instant tkt_{k}; it first updates the disturbance estimate using the measured state and the one-step prediction error from the nominal model.

Algorithm 1 Robust SDC-Based MPC Algorithm
1:Initialization: Set initial state 𝐱0\mathbf{x}_{0}, disturbance estimate 𝐝^0=𝟎\hat{\mathbf{d}}_{0}=\mathbf{0}, time step k=0k=0.
2:repeat
3:  Measure the current state of the system, 𝐱meas(tk)\mathbf{x}_{\text{meas}}(t_{k}).
4:  if k>0k>0 then
5:   Predict the state using the nominal model from the previous step:
6:  𝐱pred(tk)=Fnominal(𝐱meas(tk1),𝐮(tk1))\mathbf{x}_{\text{pred}}(t_{k})=F_{\text{nominal}}(\mathbf{x}_{\text{meas}}(t_{k-1}),\mathbf{u}(t_{k-1})).
7:   Calculate the one-step prediction error: 𝐞k=𝐱meas(tk)𝐱pred(tk)\mathbf{e}_{k}=\mathbf{x}_{\text{meas}}(t_{k})-\mathbf{x}_{\text{pred}}(t_{k}).
8:   Infer the disturbance measurement, 𝐝^meas(tk)\hat{\mathbf{d}}_{\text{meas}}(t_{k}), from the error 𝐞k\mathbf{e}_{k}.
9:   Update the disturbance estimate using the low-pass filter:
10:  𝐝^(tk)=α𝐝^(tk1)+(1α)𝐝^meas(tk)\hat{\mathbf{d}}(t_{k})=\alpha\hat{\mathbf{d}}(t_{k-1})+(1-\alpha)\hat{\mathbf{d}}_{\text{meas}}(t_{k}).
11:  end if
12:  Define the current linearization point 𝐱lin,k,𝐮lin,k\mathbf{x}_{\text{lin},k},\mathbf{u}_{\text{lin},k} (e.g., from the reference trajectory).
13:  Evaluate the state-dependent matrices 𝐀k=𝐀(𝐱lin,k)\mathbf{A}_{k}=\mathbf{A}(\mathbf{x}_{\text{lin},k}) and 𝐁k=𝐁(𝐱lin,k)\mathbf{B}_{k}=\mathbf{B}(\mathbf{x}_{\text{lin},k}).
14:  Compute the nominal offset term 𝐂k=f(𝐱lin,k,𝐮lin,k)𝐀k𝐱lin,k𝐁k𝐮lin,k\mathbf{C}_{k}=f(\mathbf{x}_{\text{lin},k},\mathbf{u}_{\text{lin},k})-\mathbf{A}_{k}\mathbf{x}_{\text{lin},k}-\mathbf{B}_{k}\mathbf{u}_{\text{lin},k}.
15:  Form the robust prediction model for the QP solver:
16:  𝐱i+1|k=𝐀k𝐱i|k+𝐁k𝐮i|k+𝐂k+𝐄d𝐝^(tk)\mathbf{x}_{i+1|k}=\mathbf{A}_{k}\mathbf{x}_{i|k}+\mathbf{B}_{k}\mathbf{u}_{i|k}+\mathbf{C}_{k}+\mathbf{E}_{d}\hat{\mathbf{d}}(t_{k}).
17:  Solve the MPC optimization problem to find the optimal control sequence 𝐔k={𝐮0|k,,𝐮N1|k}\mathbf{U}^{*}_{k}=\{\mathbf{u}^{*}_{0|k},...,\mathbf{u}^{*}_{N-1|k}\}.
18:  Apply the first element of the sequence to the plant: 𝐮(tk)=𝐮0|k\mathbf{u}(t_{k})=\mathbf{u}^{*}_{0|k}.
19:  Set kk+1k\leftarrow k+1.
20:until end of simulation

IV Simulation Results and Performance Analysis

To validate the theoretical claims and quantitatively assess the performance of the proposed control frameworks, a comprehensive simulation study was conducted. The study benchmarks three distinct controllers: the high-fidelity NMPC, the nominal SDC-MPC, and the proposed Robust SDC-based MPC. The controllers were evaluated on a challenging, safety-critical mission requiring a quadrotor to track a dynamic reference trajectory while actively avoiding a static cylindrical obstacle. To simulate realistic operating conditions, the system was subjected to stochastic external disturbances affecting both translational and rotational dynamics. All simulations were implemented in MATLAB R2025b, utilizing the CasADi optimization framework with the IPOPT solver for both the NLP and QP problems. The IPOPT solver was used for all optimization problems to ensure a fair comparison focused on the complexity of the problem formulations, independent of solver-specific optimizations. It is noted, however, that employing a dedicated QP solver for the SDC-based methods could potentially yield even greater computational speed-ups. The computations were executed on a personal laptop equipped with an Apple M4 processor and 24 GB of RAM, providing a consistent hardware baseline for comparing the computational performance of each controller. A consistent set of parameters was used across the simulation to ensure a fair comparison. The quadrotor model was configured with a mass m=1.0m=1.0 kg and moments of inertia Jx=Jy=0.029J_{x}=J_{y}=0.029 kg\cdotm2 and Jz=0.055J_{z}=0.055 kg\cdotm2. For the MPC design, a sampling time of T=0.1T=0.1 s and a prediction horizon of N=20N=20 steps were chosen. The objective function weighting matrices were defined as 𝐐=diag([50,50,80,20,20,20,10,10,10,2,2,2])\mathbf{Q}=\text{diag}([50,50,80,20,20,20,10,10,10,2,2,2]) and 𝐑=diag([0.1,0.5,0.5,0.2])\mathbf{R}=\text{diag}([0.1,0.5,0.5,0.2]). The control inputs were constrained to 𝐮[[0,20]T,[1,1]T,[1,1]T,[0.5,0.5]T]\mathbf{u}\in[[0,20]^{T},[-1,1]^{T},[-1,1]^{T},[-0.5,0.5]^{T}]. The static obstacle was modeled as a vertical cylinder centered at [1,1][-1,-1] with a radius of 0.50.5 m. Finally, the external disturbances were modeled as zero-mean Gaussian noise with standard deviations of 0.30.3 for forces and 0.10.1 for torques. Figure 3 presents the trajectory of the quadrotor, providing a qualitative assessment of tracking performance and safety constraint satisfaction for the three controllers. A primary observation is that all three methods successfully navigate the static obstacle, deviating from the reference path to maintain a safe distance as enforced by the MPC constraints; this confirms the efficacy of the obstacle avoidance formulation. A comparative analysis of the trajectories reveals the superior performance of the proposed Robust SDC-based MPC, which maintains significantly tighter adherence to the reference path despite the presence of stochastic disturbances. In contrast, both the NMPC and the nominal SDC-MPC exhibit more pronounced oscillations and tracking errors. The inset provides a magnified view of the trajectories near the barrier, confirming that all controllers respect the required safety margin.

Refer to caption
Figure 3: The quadrotor trajectory for the safety-critical tracking mission. A primary observation is that all three controllers successfully deviate from the reference path (black dotted line) to satisfy the safety constraint imposed by the static obstacle (yellow circle), confirming the efficacy of the obstacle avoidance formulation. A comparative analysis of tracking performance reveals that the proposed Robust SDC-based MPC (green solid line) exhibits superior adherence to the reference trajectory, particularly in the presence of stochastic disturbances. The NMPC (blue solid line) and the nominal SDC-MPC (red dashed line) show more pronounced deviations from the desired path. The inset provides a magnified view of the trajectories in the vicinity of the obstacle, confirming that all controllers maintain the required safety margin.

A quantitative comparison of the controllers’ performance is summarized in Table I. The table presents three key metrics: the average CPU time required per control step, the Mean Squared Error (MSE) for position tracking, and the total accumulated cost over the entire mission. These metrics provide a clear and concise assessment of the trade-offs between computational efficiency, tracking accuracy, and overall control effort for each of the benchmarked methods.

TABLE I: Quantitative Performance Comparison of MPC Controllers
Metric NMPC SDC-MPC Robust SDC-MPC
Avg. CPU (ms) 25.85 16.28 16.40
MSE (m2m^{2}) 0.0304 0.0302 0.0321
Total Cost 5.47e+04 5.49e+04 5.44e+04
Refer to caption
Figure 4: Per-iteration CPU performance for the three benchmarked controllers. The plot illustrates the computational load required to solve the optimization problem at each time step. Both SDC-based controllers (red and green) consistently exhibit lower average computation times than the NMPC (blue), underscoring the efficiency gained by the proposed method. A significant spike in CPU time is observable for all controllers around the 11-13 second mark, which corresponds to the challenging maneuver where the quadrotor is actively avoiding the obstacle. Even during this computationally intensive phase, the SDC-based methods remain significantly faster than the NMPC, highlighting their advantage for real-time, safety-critical applications.

The per-iteration computational performance, depicted in Figure 4, further reinforces the advantages of the proposed SDC-based methods. The plot illustrates the CPU time required to solve the optimization problem at each control step. A consistent trend is observed where both the nominal and robust SDC-based controllers exhibit a significantly lower computational burden compared to the full NMPC benchmark. This efficiency is particularly critical during the most demanding phase of the mission, from approximately 11 to 13 seconds, where the safety constraints for obstacle avoidance become highly active. During this period, the solution time for all controllers increases, yet the SDC-based methods remain substantially faster, confirming their suitability for real-time implementation in computationally constrained, safety-critical scenarios. Figure 5 illustrates the evolution of the optimal cost function value, JJ^{*}, at each time step for the three controllers. The cost value serves as a proxy for the controller’s perceived difficulty in satisfying its objectives; higher costs indicate a greater combined penalty from tracking errors, control effort, and constraint violations. The plot reveals several distinct phases: an initial transient phase where the cost decreases as the quadrotor converges to the path, a significant increase in cost during the obstacle avoidance maneuver (around 11-13 seconds) where the controller must balance tracking with safety, and subsequent periods of stabilization. Notably, the cost profiles for all three controllers are qualitatively similar, indicating that the SDC-based methods achieve comparable overall performance to the NMPC benchmark, as corroborated by the total cost values in Table I.

Refer to caption
Figure 5: Evolution of the optimal cost function value at each sampling instant. The logarithmic scale on the y-axis highlights the relative changes in the cost. The cost profiles are comparable across all controllers, with significant increases corresponding to challenging maneuvers, particularly the obstacle avoidance phase.

The position tracking performance in the X, Y, and Z axes is detailed in Figure 6. All controllers effectively track the time-varying reference trajectory across all three dimensions. The top subplot shows the X-position, where the step change in the reference after 4 seconds is handled well by all methods, though the Robust SDC-MPC demonstrates the least overshoot. The middle and bottom subplots for the Y and Z positions further highlight the superior tracking accuracy of the Robust SDC-MPC, which exhibits smaller deviations from the reference, especially during the dynamic circular portion of the trajectory. This improved performance is directly attributable to the disturbance observer, which actively compensates for unmodeled dynamics and external perturbations, leading to more precise tracking.

Refer to caption
Figure 6: Position tracking performance for the X, Y, and Z axes. The black dotted line is the reference trajectory. All controllers demonstrate effective tracking, with the Robust SDC-MPC (green) showing the highest fidelity to the reference path throughout the maneuver.

To further quantify the tracking accuracy, the squared position error for each controller is plotted over time in Figure 7. This plot corroborates the qualitative observations from previous figures and provides a clearer view of the performance differences, particularly due to the logarithmic scale. The Robust SDC-MPC consistently achieves the lowest tracking error throughout the simulation, with error levels often an order of magnitude smaller than the other controllers, especially during less dynamic phases of the trajectory (e.g., around 10 seconds and 15 seconds). This demonstrates the effectiveness of the disturbance observer in actively rejecting perturbations and reducing steady-state tracking error, a conclusion strongly supported by the lower overall MSE for this controller reported in Table I.

Refer to caption
Figure 7: Squared position tracking error over time for the three controllers. The logarithmic y-axis highlights the superior performance of the Robust SDC-MPC, which consistently maintains a lower error profile.

IV-A Discussion of Limitations

While the proposed SDC-based framework demonstrates significant advantages, it is important to acknowledge its limitations. The performance of the robust controller relies on the assumption that disturbances are slowly time-varying. Its effectiveness may be reduced in the presence of high-frequency or abrupt disturbances that violate this assumption. Additionally, the SDC method is an exact algebraic transformation at a specific point in time. For highly aggressive maneuvers over a long prediction horizon, the system state may diverge significantly from the point of linearization, potentially degrading the accuracy of the prediction model. The stability and recursive feasibility guarantees also rely on a terminal set derived from a local linearization, which, while standard, may be conservative. Future work will seek to address these limitations, potentially through the use of more advanced disturbance observers or adaptive prediction models.

V Conclusion

This paper addressed the challenge of implementing high-performance, safety-critical control for quadrotors under significant computational constraints. We proposed a novel SDC-based MPC framework and a robust extension that incorporates a disturbance observer. The SDC reformulation allows the nonlinear optimal control problem to be solved as a sequence of computationally efficient QPs, directly addressing the limitations of traditional NMPC. The simulation results conclusively demonstrate the efficacy of the proposed approach. The Robust SDC-based MPC not only matched the safety and tracking performance of the NMPC benchmark but, through active disturbance compensation, achieved superior tracking accuracy. Critically, this enhanced performance was realized with a significant reduction in computational time, with the SDC-based methods proving to be over 30% faster on average than the NMPC. The key contribution of this work is the demonstration that the SDC-based methodology successfully achieves a compelling balance between safety-critical performance and computational efficiency, establishing it as a highly feasible alternative to NMPC for real-time quadrotor applications. Future research will prioritize the experimental validation of the proposed controller on a physical quadrotor platform, which will be crucial for assessing its real-world performance under conditions of sensor noise and communication latency. Further avenues of investigation include the exploration of more advanced disturbance estimation techniques, the development of adaptive SDC formulations to handle highly aggressive maneuvers, and the extension of this computationally efficient framework to cooperative control of multi-agent aerial systems.

References

  • [1] J. Peksa and D. Mamchur, “A review on the state of the art in copter drones and flight control systems,” Sensors, vol. 24, no. 11, p. 3349, 2024.
  • [2] S. Allan and M. Barczyk, “A low-cost experimental quadcopter drone design for autonomous search-and-rescue missions in gnss-denied environments,” Drones, vol. 9, no. 8, p. 523, 2025.
  • [3] Y. Wang, H. Ji, Q. Kang, H. Qi, and J. Wen, “Autonomous trajectory control for quadrotor evtol in hover and low-speed flight via the integration of model predictive and following control,” Drones, vol. 9, no. 8, p. 537, 2025.
  • [4] S. S. Unde, V. Kurkute, S. S. Chavan, D. D. Mohite, A. A. Harale, and A. Chougle, “The expanding role of multirotor uavs in precision agriculture with applications ai integration and future prospects,” Discover Mechanical Engineering, vol. 4, no. 1, pp. 1–26, 2025.
  • [5] M. Rinaldi, S. Primatesta, and G. Guglieri, “A comparative study for control of quadrotor uavs,” Applied Sciences, vol. 13, no. 6, p. 3464, 2023.
  • [6] C. Massé, O. Gougeon, D.-T. N’Guyen, and D. Saussié, “Modeling and control of a quadcopter flying in a wind field: A comparison between lqr and structured h\infty control techniques,” in 2018 International Conference on Unmanned Aircraft Systems (ICUAS). IEEE, 2018, pp. 1408–1417.
  • [7] C.-C. Chen and Y.-T. Chen, “Feedback linearized optimal control design for quadrotor with multi-performances,” IEEE Access, vol. 9, pp. 26 674–26 695, 2021.
  • [8] J. Wang, K. A. Alattas, Y. Bouteraa, O. Mofid, and S. Mobayen, “Adaptive finite-time backstepping control tracker for quadrotor uav with model uncertainty and external disturbance,” Aerospace Science and Technology, vol. 133, p. 108088, 2023.
  • [9] B. Mu, K. Zhang, and Y. Shi, “Integral sliding mode flight controller design for a quadrotor and the application in a heterogeneous multi-agent system,” IEEE Transactions on Industrial Electronics, vol. 64, no. 12, pp. 9389–9398, 2017.
  • [10] S. Bouabdallah, A. Noth, and R. Siegwart, “PID vs LQ control techniques applied to an indoor micro quadrotor,” in 2004 IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS)(IEEE Cat. No. 04CH37566), vol. 1. IEEE, 2004, pp. 345–350.
  • [11] J. Hwangbo, I. Sa, R. Siegwart, and M. Hutter, “Control of a quadrotor with reinforcement learning,” IEEE Robotics and Automation Letters, vol. 2, no. 4, pp. 2096–2103, 2017.
  • [12] A. D. Ames, S. Coogan, M. Egerstedt, G. Notomista, K. Sreenath, and P. Tabuada, “Control Barrier Functions: Asymptotic Safety and Stability,” IEEE Transactions on Automatic Control, vol. 62, no. 7, pp. 3161–3176, 2017.
  • [13] I. M. Mitchell, “A toolbox of level set methods,” in Proceedings of the 2005 international conference on Svetlogorsk. Citeseer, 2005, pp. 1–12.
  • [14] S. Prajna, A. Jadbabaie, and G. J. Pappas, “A framework for analysis of nonlinear systems using barrier certificates,” IEEE Transactions on Automatic Control, vol. 52, no. 1, pp. 23–35, 2007.
  • [15] A. Bemporad, M. Morari, V. Dua, and E. N. Pistikopoulos, “The explicit linear quadratic regulator for constrained systems,” in Automatica, vol. 38, no. 1. Elsevier, 2002, pp. 3–20.
  • [16] M. Farina, L. Giulioni, and R. Scattolini, “Stochastic model predictive control: An overview and perspectives for future research,” Annual Reviews in Control, vol. 42, pp. 115–128, 2016.
  • [17] U. Rosolia and F. Borrelli, “Learning model predictive control for iterative tasks,” in 2017 IEEE 56th Annual Conference on Decision and Control (CDC). IEEE, 2017, pp. 3272–3277.
  • [18] B. Xu, A. Suleman, and Y. Shi, “A multi-rate hierarchical fault-tolerant adaptive model predictive control framework: Theory and design for quadrotors,” Automatica, vol. 153, p. 111015, 2023.
  • [19] T. Cimen, “State-dependent riccati equation (sdre) control: A survey,” in Proceedings of the 17th world congress, the international federation of automatic control, vol. 17, no. 1, 2008, pp. 3761–3775.
  • [20] J. R. Cloutier, “State-dependent riccati equation techniques: an overview,” in Proceedings of the 1997 American Control Conference, vol. 2. IEEE, 1997, pp. 932–936.
  • [21] C. Mracek and J. Cloutier, “Control designs for the nonlinear benchmark problem via the state-dependent riccati equation method,” International Journal of Robust and Nonlinear Control: IFAC-Affiliated Journal, vol. 8, no. 4-5, pp. 401–416, 1998.
  • [22] K. D. Hammett, N. Harl, H. An, and E. Zuazua, “State-dependent coefficient factorization for nonlinear optimal feedback control,” IEEE Transactions on Automatic Control, vol. 62, no. 5, pp. 2497–2502, 2016.
  • [23] H. Chen and F. Allgöwer, “A Quasi-Infinite Horizon Nonlinear Model Predictive Control Scheme with Guaranteed Stability,” Automatica, vol. 34, no. 10, pp. 1205–1217, 1998.
  • [24] J. B. Mare and J. A. De Doná, “Solution of the input-constrained lqr problem using dynamic programming,” Systems & control letters, vol. 56, no. 5, pp. 342–348, 2007.
  • [25] Y. Batmani, M. Davoodi, and N. Meskin, “Nonlinear suboptimal tracking controller design using state-dependent riccati equation technique,” IEEE Transactions on Control Systems Technology, vol. 25, no. 5, pp. 1833–1839, 2016.
  • [26] J. B. Rawlings, D. Q. Mayne, and M. M. Diehl, Model Predictive Control: Theory, Computation, and Design. Nob Hill Publishing, 2017, vol. 2.
  • [27] J. B. Rawlings and K. R. Muske, “The stability of constrained receding horizon control,” IEEE Transactions on Automatic Control, vol. 38, no. 10, pp. 1512–1516, 1993.
  • [28] W.-H. Chen, J. Yang, L. Guo, and S. Li, “Disturbance-observer-based control and related methods—an overview,” IEEE Transactions on Industrial Electronics, vol. 63, no. 2, pp. 1083–1095, 2016.
  • [29] H. Wang, D. Zhou, and P. Li, “Disturbance observer-based robust control for nonlinear systems,” Acta Automatica Sinica, vol. 35, no. 6, pp. 721–727, 2009.
  • [30] D. Q. Mayne, J. B. Rawlings, C. V. Rao, and P. O. Scokaert, Constrained model predictive control: Stability and optimality. Nob Hill Publishing, 2000, vol. 2.
[Uncaptioned image] Saber Omidi (Student Member, IEEE) received his B.Sc. degree in Mechanical Engineering from Qazvin Islamic Azad University, Qazvin, Iran, in 2017, and his M.Sc. degree in Mechanical Engineering from University of Kurdistan, Sanandaj, Iran, in 2021. He is currently pursuing his Ph.D. degree in Mechanical Engineering at University of New Hampshire. His research interests include safety-critical control, learning-based control, with applications in autonomous systems and mechanical systems.