Minimum lap-time trajectory planning; stochastic vehicle dynamics; probabilistic safety certificates; parameter uncertainty; parsimonious formulation; MPC

1 Introduction↩︎

Minimum-lap-time planning (MLTP) is a standard tool in motorsport engineering: by solving a constrained optimal control problem (OCP), it synthesises the time-optimal trajectory and control inputs for a given vehicle racing on a circuit. Engineers use these references offline, to compare vehicle setups and race strategies, and online, as feedforward signals for driver-assistance and autonomous-racing systems. In this context, optimality typically refers to lap-time minimisation. Since MLTP drives the vehicle as fast as the dynamic model and the OCP constraints allow, the resulting trajectory necessarily skirts the tyre-friction limit with essentially no margin to spare.

A trajectory with no margin is a fragile one to follow. Any disturbance—such as a kerb strike or vehicle parameters that deviate from their nominal values—can push the vehicle past the limit the reference assumed. The resulting reference is then hard, or even unsafe, to reproduce, whether followed by an expert driver or a tracking controller.

We address that weakness by planning offline robust minimum-lap-time references. Lap time remains the objective, but we tighten tyre-friction limits probabilistically so that the robust reference carries explicit safety margins precisely where they are most beneficial along the track. Most MLTP formulations do not pursue that goal directly: many works adopt heuristic conservative margins or different assumptions on performance envelopes rather than embedding probabilistic robustness in the optimal control problem itself. Even the frameworks that do embed it consider uncertainty in the vehicle state alone, leaving the vehicle’s own parameters fixed. We close this gap by extending disturbance-aware minimum-lap-time planning to joint state and parameter uncertainty, providing an efficient parsimonious strategy to offset the increased computational burden, and assessing the robust MLTP references through a Monte Carlo validation campaign.

1.1 Related work↩︎

Minimum-lap-time planning spans a broad and well-surveyed landscape [1]. Quasi-steady-state formulations keep the state space small and push vehicle complexity into performance-envelope constraints [2]–[4], whereas high-fidelity multibody models trade computation for accuracy on three-dimensional circuits [5], [6]. Orthogonal efforts speed up these deterministic problems through better OCP formulation [7], [8] or, as the program grows, through parallelization strategies [9], [10]. In contrast to these offline-planning methods, a second branch of research adopts online controllers to recover much of the gap to an offline-optimal reference on a comparable high-fidelity model [11]. Yet across all these studies the MLTP stays nominal.

A smaller group of methods instead accounts for disturbances within planning, trading some lap time for a reference easier to follow under perturbation. Among online-planning studies, [12] plans a nominal lap-time trajectory and then adds a feedback controller that rejects disturbances inferred from a road-roughness model, reducing the steering effort needed to track the reference, thereby enhancing its driveability. The offline-planning framework [13] instead embeds robustness inside a single large-scale MLTP, modelling state uncertainty and propagating the state covariance. That framework enforces a probabilistic back-off on the OCP constraints, tightening the stay-on-track and friction-limit constraints. Those frameworks, however, assume nominal values for the vehicle’s parameters.

Outside motorsport, trajectory optimisation under parametric uncertainty is a far more developed topic, e.g., in robotics and control applications. The dominant approach minimises the trajectory’s sensitivity to the uncertain parameters as a term in the cost [14], an idea refined from open-loop to closed-loop sensitivity [15] and later extended to the control inputs [16]. Treating robustness as a cost penalty suits objectives that tolerate leaving the nominal optimum for a less sensitive one. Minimum-time planning does not: its optimum lies on the friction limit, so penalising sensitivity only shifts the nominal point along that limit without building in any margin. Under a real perturbation, the trajectory then meets the limit with nothing in reserve. The penalty weight, moreover, fixes no probability of respecting the limit. We therefore keep the same modelling stance—key vehicle parameters treated as uncertain—but move robustness from the cost into the constraints, tightening the friction limit by a chance-constrained back-off sized to a chosen probability of staying within it, and costing nothing where that limit is slack.

Extending robustness to parameters is not free. Propagating the joint state-and-parameter covariance and tightening the OCP constraints against it enlarges the nonlinear program markedly. Parsimonious strategies, however, may contain this cost. Concentrating the computational effort where the problem is locally hardest is the principle that drives adaptive refinement in numerical optimal control [17]. Here, we activate the robust machinery only on critical and near-critical track regions, with an original method based on Lagrange multipliers and residuals of the OCP constraints.

Finally, assessing how well a robust reference survives execution inherently requires closed-loop evaluation. Validating disturbance-aware references from the same planning framework with professional drivers on a dynamic driving simulator offers realistic but sparse evidence [18]—a handful of runs per condition. Model predictive control (MPC) closes the loop computationally instead, tracking a planned reference in the presence of disturbances, and Monte Carlo campaigns are a standard way to probe a tracking controller’s robustness under uncertainty [19]–[21]; in those studies, the campaign certifies the controller itself. We take the complementary route: a nonlinear MPC acts as a virtual test driver tracking externally planned references in a Monte Carlo campaign. Survival rates and control effort then become quantitative measures of driveability, and lap-time increase quantifies the cost of robustness.

1.2 Paper’s contributions and organisation↩︎

The paper makes three contributions.

• Parametric uncertainty. We extend the prior disturbance-aware MLTP framework [13] to propagate Gaussian uncertainty in the yaw moment of inertia, centre-of-mass position, and aerodynamic drag coefficient within the same probabilistic back-off as the state disturbances.

• Parsimonious activation. We confine the robust machinery to critical and near-critical nodes, recovering the full-grid-enforcement trajectory at a fraction of the solve time.

• MPC Monte Carlo validation. An MPC virtual driver tracks planned MLTP references over 1000 runs each, showing that the robust references satisfy the friction-limit constraints more reliably, at a moderate lap-time cost and with lower steering effort.

The remainder of the paper is organised as follows. Section 2 introduces the stochastic vehicle-dynamics framework: the single-track model with nonlinear tyres, the continuous- and discrete-time planning problems, and the probabilistic back-off used to tighten the friction-limit constraint. Section 3 formulates the resulting open-loop covariance-propagation planning problem and its parsimonious variant. Section 4 describes the experimental setup—reference planning, MPC design, and the Monte Carlo testing strategy—while Section 5 presents and discusses the results. Section 6 concludes and outlines directions for future work.

2 Handling model uncertainty↩︎

2.1 Stochastic vehicle model↩︎

We adopt a Single Track model [22] with nonlinear axle characteristics, illustrated in Figure 1.

The state vector \(\boldsymbol{x}\), with \(n_x=6\) components, comprises the longitudinal and lateral velocities \(u\) and \(v\) of the centre of mass (CoM) with respect to its body-fixed frame, the yaw rate \(r\), the planar position \((x_G, y_G)\), and the heading \(\psi\) (Figure 1). Four uncertain parameters—the yaw moment of inertia \(J_z\), the CoM height \(h\), the weight balance \(w_b = \frac{a_2}{a_1+a_2}\) (with \(a_1\), \(a_2\) the distances from the CoM to the front and rear axles), and the longitudinal drag coefficient \(C_x\)—form the parameter vector \(\boldsymbol{p}\), with \(n_p=4\). Table 1 reports the nominal values of the main vehicle parameters; hereafter, parameters refers exclusively to these four uncertain quantities (shaded in Table 1), for which the listed values are the distribution means. The full parameter set is available as supplementary material [23].

Tracking the joint state-parameter distribution requires a single covariance matrix; we therefore concatenate \(\boldsymbol{x}\) and \(\boldsymbol{p}\) into the augmented state \(\boldsymbol{\eta}= [\boldsymbol{x};\boldsymbol{p}]\), with \(n_\eta = n_x + n_p\) components. The control inputs are the total longitudinal force \(X = X_1+X_2\) and the front steering angle \(\delta\), collected in the input vector \(\boldsymbol{u}\), with \(n_u=2\).

The rear-wheel-drive configuration of the modelled FSAE vehicle implies that a positive (driving) longitudinal force \(X\) acts only on the rear axle, while a negative (braking) force is distributed between the front and rear axles, according to a fixed brake balance \(\frac{X_1}{X_2}\). Since the input and the fixed brake balance assign the longitudinal forces \(X_1\) and \(X_2\), the computation of vertical load transfers is explicit and allows direct evaluation of the lateral forces \(Y_1\) and \(Y_2\) using Pacejka’s Magic Formula.

The perturbed vehicle dynamics thus takes the explicit form \[\label{eq:dynamics} \dot{\boldsymbol{\eta}}(t) = f(\boldsymbol{\eta}(t),\boldsymbol{u}(t))+\boldsymbol{w}(t),\tag{1}\] where \(\boldsymbol{w}(t)\) is additive Gaussian white noise with zero mean and known covariance \(\boldsymbol{Q}(t)\), i.e., \(\boldsymbol{w}(t) \sim \mathcal{N}(\boldsymbol{0}, \boldsymbol{Q}(t))\). Assuming an initially Gaussian distribution on \(\boldsymbol{\eta}\) and linearising the covariance dynamics to first order, the distribution remains Gaussian throughout, with mean \(\boldsymbol{\mu}(t)\) and covariance \(\boldsymbol{P}(t)\), i.e., \(\boldsymbol{\eta}(t) \sim \mathcal{N}(\boldsymbol{\mu}(t), \boldsymbol{P}(t))\).

Since the vehicle parameters are uncertain but constant, Equation 1 can be expanded into \[\begin{align} \label{eq:expanded95dyn} \begin{bmatrix} \dot{\boldsymbol{x}}(t) \\[0.25em] \dot{\boldsymbol{p}}(t) \end{bmatrix} = \begin{bmatrix} f_x(\boldsymbol{x}(t),\boldsymbol{p}(t),\boldsymbol{u}(t)) \\[0.25em] \mathbf{0} \end{bmatrix} + \begin{bmatrix} \boldsymbol{w}_x(t) \\[0.25em] \mathbf{0} \end{bmatrix}. \end{align}\tag{2}\]

Since the distribution of \(\boldsymbol{\eta}(t)\) remains Gaussian, the mean \(\boldsymbol{\mu}(t)\) and the covariance \(\boldsymbol{P}(t)\) fully characterise its time evolution.

The mean \(\boldsymbol{\mu}(t)\) follows the deterministic dynamics given by \[\label{eq:mean95dyn} \dot{\boldsymbol{\mu}}(t) = f(\boldsymbol{\mu}(t), \boldsymbol{u}(t)).\tag{3}\] The covariance \(\boldsymbol{P}(t)\) evolves according to the Lyapunov matrix differential equation [24] \[\begin{align} \label{eq:dP} \dot{\boldsymbol{P}}(t) = \boldsymbol{A}(t) \boldsymbol{P}(t) + \boldsymbol{P}(t) \boldsymbol{A}^T(t) + \boldsymbol{Q}(t),\quad \boldsymbol{P}(t_0) = \boldsymbol{P}_0 = \boldsymbol{P}_0^T. \end{align}\tag{4}\] In Equation 4 , \(\boldsymbol{P}_0\) is the initial augmented state covariance and \(\boldsymbol{A}(t)\) is the usual shorthand notation for the Jacobian along the mean trajectory \(\boldsymbol{\mu}(t)\), that is \(\boldsymbol{A}(t)=\frac{\partial f(\boldsymbol{\mu}, \boldsymbol{u})}{\partial\boldsymbol{\mu}}\).

To gain insight into 4 , we show how a diagonal \(\boldsymbol{P}_0\)—with state and parameter uncertainties assumed independent—propagates, by partitioning the matrices conformally with \(\boldsymbol{\eta}= [\boldsymbol{x};\boldsymbol{p}]\) and shading the non-zero blocks: \[\begin{align} \label{eq:dP95color1} \dot{\boldsymbol{P}}(0) = \underbrace{ \begin{bNiceMatrix}[margin, columns-width = 0.4em, cell-space-top-limit = 0.4em, cell-space-bottom-limit = 0.4em] \Block[fill=gray!25, tikz={pattern={Lines[angle=45,distance=6pt,line width=0.3pt]}, pattern color=black}]{2-3}{} & & \\ & & \\ \hline \Block[fill=white]{1-3}{} & & \end{bNiceMatrix} }_{\displaystyle \boldsymbol{A}(0)} \underbrace{ \begin{bNiceMatrix}[margin,vlines=3, columns-width = 0.4em, cell-space-top-limit = 0.4em, cell-space-bottom-limit = 0.4em] \Block[fill=red!25]{2-2}{} & & \\ & & \\ \hline & & \Block[fill=blue!25]{1-1}{} \end{bNiceMatrix} }_{\displaystyle \boldsymbol{P}_0} + \underbrace{ \begin{bNiceMatrix}[margin,vlines=3, columns-width = 0.4em, cell-space-top-limit = 0.4em, cell-space-bottom-limit = 0.4em] \Block[fill=red!25]{2-2}{} & & \\ & & \\ \hline & & \Block[fill=blue!25]{1-1}{} \end{bNiceMatrix} }_{\displaystyle \boldsymbol{P}_0} \underbrace{ \begin{bNiceMatrix}[margin,vlines=3, columns-width = 0.4em, cell-space-top-limit = 0.4em, cell-space-bottom-limit = 0.4em] \Block[fill=gray!25, tikz={pattern={Lines[angle=45,distance=6pt,line width=0.3pt]}, pattern color=black}]{3-2}{} & & \Block[fill=white]{3-1}{} \\ & & \\ & & \end{bNiceMatrix} }_{\displaystyle \boldsymbol{A}^T(0)} + \underbrace{ \begin{bNiceMatrix}[margin,vlines=3, columns-width = 0.4em, cell-space-top-limit = 0.4em, cell-space-bottom-limit = 0.4em] \Block[fill=gray!25, tikz={pattern={Lines[angle=135,distance=6pt,line width=0.3pt]}, pattern color=black}]{2-2}{} & & \Block[fill=white]{2-1}{} \\ & & \\ \hline \Block[fill=white]{1-2}{} & & \Block[fill=white]{1-1}{} \end{bNiceMatrix} }_{\displaystyle \boldsymbol{Q}} \end{align}\tag{5}\] The red and blue blocks represent the initial state and parameter uncertainty, respectively. The resulting structure of \(\dot{\boldsymbol{P}}(0)\), \[\begin{align} \label{eq:dP95color2} \underbrace{ \begin{bNiceMatrix}[margin,vlines=3, columns-width = 0.4em, cell-space-top-limit = 0.4em, cell-space-bottom-limit = 0.4em] \Block[fill=red!32, tikz={pattern={Lines[angle=45,distance=6pt,line width=0.3pt]}, pattern color=black}]{2-2}{} & & \Block[fill=blue!32, tikz={pattern={Lines[angle=45,distance=6pt,line width=0.3pt]}, pattern color=black}]{2-1}{} \\ & & \\ \hline & & \end{bNiceMatrix} }_{\displaystyle \boldsymbol{A}(0) \boldsymbol{P}_0} + \underbrace{ \begin{bNiceMatrix}[margin,vlines=3, columns-width = 0.4em, cell-space-top-limit = 0.4em, cell-space-bottom-limit = 0.4em] \Block[fill=red!32, tikz={pattern={Lines[angle=45,distance=6pt,line width=0.3pt]}, pattern color=black}]{2-2}{} & & \\ & & \\ \hline \Block[fill=blue!32, tikz={pattern={Lines[angle=45,distance=6pt,line width=0.3pt]}, pattern color=black}]{1-2}{} & & \end{bNiceMatrix} }_{\displaystyle \boldsymbol{P}_0 \boldsymbol{A}^T(0)} + \underbrace{ \begin{bNiceMatrix}[margin,vlines=3, columns-width = 0.4em, cell-space-top-limit = 0.4em, cell-space-bottom-limit = 0.4em] \Block[fill=gray!25, tikz={pattern={Lines[angle=135,distance=6pt,line width=0.3pt]}, pattern color=black}]{2-2}{} & & \Block[fill=white]{2-1}{} \\ & & \\ \hline \Block[fill=white]{1-2}{} & & \Block[fill=white]{1-1}{} \end{bNiceMatrix} }_{\displaystyle \boldsymbol{Q}} = \underbrace{ \begin{bNiceMatrix}[margin,vlines=3, columns-width = 0.4em, cell-space-top-limit = 0.4em, cell-space-bottom-limit = 0.4em] \Block[fill=red!40, tikz={postaction={pattern={Lines[angle=45,distance=6pt,line width=0.3pt]}, pattern color=black}, postaction={pattern={Lines[angle=-45,distance=6pt,line width=0.3pt]}, pattern color=black}}]{2-2}{} & & \Block[fill=blue!32, tikz={pattern={Lines[angle=45,distance=6pt,line width=0.3pt]}, pattern color=black}]{2-1}{} \\ & & \\ \hline \Block[fill=blue!32, tikz={pattern={Lines[angle=45,distance=6pt,line width=0.3pt]}, pattern color=black}]{1-2}{} & & \end{bNiceMatrix} }_{\displaystyle \dot{\boldsymbol{P}}(0)} \end{align}\tag{6}\] reveals the asymmetric behaviour of 4 : i) state uncertainty drives the state auto-covariance dynamics (top-left, red); ii) parameter uncertainty propagates into the cross-covariance blocks (off-diagonal, blue); iii) the parameter auto-covariance (bottom-right) has zero dynamics. Therefore, parameter uncertainty influences state evolution, but not vice versa.

The analytical solution of 4 , verifiable by differentiation [24], takes the form \[\begin{align} \label{eq:P95from95STM} \boldsymbol{P}(t) = \boldsymbol{\Phi}(t,t_0) \boldsymbol{P}_0 \boldsymbol{\Phi}^T(t,t_0)+\int_{t_0}^{t} \boldsymbol{\Phi}(t,\tau) \boldsymbol{Q}(\tau) \boldsymbol{\Phi}^T(t,\tau) \text{{\sl d}}\tau \quad \boldsymbol{P}(t_0)=\boldsymbol{P}_0, \end{align}\tag{7}\] where \(\boldsymbol{\Phi}(t,t_0)\) is the state transition matrix (STM). This matrix encodes the evolution from \(t_0\) to \(t\) of a perturbation with respect to the mean trajectory \(\boldsymbol{\mu}(t)\) as \(\bar{\boldsymbol{\eta}}(t) = \boldsymbol{\Phi}(t,t_0) \bar{\boldsymbol{\eta}}(t_0)\), with \(\bar{\boldsymbol{\eta}}(\cdot) =\boldsymbol{\eta}(\cdot) - \boldsymbol{\mu}(\cdot)\). The evolution of \(\boldsymbol{\Phi}(t,t_k)\) from a generic \(t_k\) to \(t\) is driven by the following differential equation \[\begin{align} \label{eq:dSTM} \dot{\boldsymbol{\Phi}}(t,t_k) = \boldsymbol{A}(t)\boldsymbol{\Phi}(t,t_k), \quad \boldsymbol{\Phi}(t_k, t_k) = \boldsymbol{I}. \end{align}\tag{8}\] Rather than integrating 4 directly, we propagate the STM via 8 and recover \(\boldsymbol{P}(t)\) a posteriori through 7 , a well-established approach in geometric and orbital mechanics for statistical trajectory determination [25], [26]. The advantage of this choice in our framework will become apparent in Section 3.1.

2.2 Robust friction limit constraint formulation↩︎

The model described in Section 2.1 employs a pure lateral Magic Formula, while longitudinal forces are treated as control inputs. Since this formulation decouples longitudinal and lateral force generation, it does not guarantee that the resulting ground force is physically feasible when both components are required simultaneously. To account for combined-slip behaviour, we enforce that the total ground-reaction forces at each axle remain within the tyre adherence ellipse.

This friction limit constraint is expressed as \[\bar{h}_j(\boldsymbol{\eta},\boldsymbol{u}) = S_j(\boldsymbol{\eta},\boldsymbol{u}) - 1 \leq 0, \qquad{(j=1,2)} \label{eq:friction95limit}\tag{9}\] where \(S_j(\boldsymbol{\eta},\boldsymbol{u})\) is the axle saturation ratio \[S_j(\boldsymbol{\eta},\boldsymbol{u}) = \frac{ \left( \frac{X_j(\boldsymbol{\eta},\boldsymbol{u})}{\mu_{x,j}} \right)^2+ \left( \frac{Y_j(\boldsymbol{\eta},\boldsymbol{u})}{\mu_{y,j}}\right)^2}{Z_j^2(\boldsymbol{\eta},\boldsymbol{u})}. \label{eq:axle95saturation}\tag{10}\] In Equations 9 –10 , \(j=1,2\) refers to the front and rear axles, respectively; \(X_j\) and \(Y_j\) denote the longitudinal and lateral components of the in-plane ground forces (Figure 1); \(Z_j\) is the total vertical load on the axle; \(\mu_{x,j}\) and \(\mu_{y,j}\) are the longitudinal and lateral friction coefficients, both set to \(1.15\) for all axles. The axle saturation ratio is equal to 0 when the in-plane ground force demand on the axle is zero, while it is equal to 1 when the point \(\left(X_j, Y_j\right)\) lies exactly on the friction ellipse, i.e., under full saturation.

Since \(\boldsymbol{\eta}\) is a random variable, \(\bar{h}_j(\boldsymbol{\eta},\boldsymbol{u})\) is itself random. We therefore impose the chance constraint \(\mathop{\mathrm{Pr}}\big\{\bar{h}_j(\boldsymbol{\eta},\boldsymbol{u}) \leq 0\big\} \geq p\), where \(p\) is a prescribed confidence level. We reformulate this chance constraint as a deterministic constraint on the mean \(\boldsymbol{\mu}\), at the cost of introducing a back-off term that tightens the nominal constraint. The back-off is defined as \(\beta_j = \gamma\sigma_j\), where \(\gamma=\Phi^{-1}(p)\) is the standard-normal quantile corresponding to the required probability \(p\) of satisfying the constraint, and \(\sigma_j(\boldsymbol{\mu}, \boldsymbol{u}, \boldsymbol{P}) = \big[\nabla_{\boldsymbol{\mu}} \bar{h}_j(\boldsymbol{\mu})^T \, \boldsymbol{P}\, \nabla_{\boldsymbol{\mu}} \bar{h}_j(\boldsymbol{\mu})\big]^{\frac{1}{2}}\) is the standard deviation of \(\bar{h}_j\) linearised around the mean \(\boldsymbol{\mu}\). The robust constraint is \[h_j(\boldsymbol{\mu},\boldsymbol{u}) = \bar{h}_j(\boldsymbol{\mu},\boldsymbol{u}) + \beta_j(\boldsymbol{\mu}, \boldsymbol{P}, \boldsymbol{u}) \leq 0, \qquad j=1,2. \label{eq:friction95limit95backoff}\tag{11}\]

3 Robust minimum-lap-time planning↩︎

We formulate the Optimal Control Problem (OCP) encoding robust minimum-lap-time planning on a motorsport track. We transcribe it into a discrete nonlinear programme (NLP) by applying the direct collocation approach [27]. The track is parameterised along its centreline by a curvilinear coordinate \(\alpha \in [0,1]\) and sampled at \(N+1\) grid nodes \(\alpha_0, \ldots, \alpha_N\). Within each interval \([\alpha_{k-1}, \alpha_k]\), the state trajectories are approximated by polynomials defined at \(d\) collocation points on the unit interval.

3.1 Multiple short-horizon propagation↩︎

Planning under state and parameter uncertainty requires propagating the covariance \(\boldsymbol{P}\). Since we tackle an open-loop planning problem, i.e., without feedback policies, position and orientation uncertainty accumulates along the track. Over a medium-to-long spatial horizon, e.g., a circuit sector, the covariance grows without bound and produces an overly conservative back-off.

Therefore, we propagate the covariance matrix for a short prediction horizon, starting from each grid point \(k\) of the discretised trajectory and accumulating uncertainty only over the following \(H\) steps, with \(H\) chosen small enough to capture the disturbance evolution before any corrective driver action, which the planning framework does not model.

Under this strategy, each grid point carries multiple versions of the covariance matrix. To this end, we introduce the notation \(\boldsymbol{P}_k^j\), where \(j = 0, \ldots, H\), and \(\boldsymbol{P}_k^j\) represents the version of the covariance matrix at step \(k\) that was initialised \(j\) steps earlier. Figure 2 illustrates the multiple instances of the covariance matrix. The dashed rectangle highlights the \(k\)-th discretisation step. At each grid point, two particular versions of the covariance matrix are emphasised: \(\boldsymbol{P}^0_k\) (red node), representing the instance to be initialised at that step, and \(\boldsymbol{P}^H_k\) (green node), corresponding to the version that has been propagated over the full horizon of \(H\) steps. This end-of-horizon instance \(\boldsymbol{P}^H_k\) is precisely the covariance matrix used to tighten at step \(k\) the friction limit constraint introduced in Section 2.2.

3.2 OCP formulation↩︎

We introduce the following notation for the decision variables in the OCP. At each grid node \(k\), \(\boldsymbol{\mu}_k\) denotes the mean augmented state and \(\boldsymbol{\Phi}_k=\boldsymbol{\Phi}(t_k,t_{k-1})\) the state transition matrix; the corresponding values at the \(d\) collocation points within interval \(k\) are collected in \(\boldsymbol{\xi}_k = (\boldsymbol{\xi}_{k,1}, \ldots, \boldsymbol{\xi}_{k,d})\) and \(\boldsymbol{\Sigma}_k = (\boldsymbol{\Sigma}_{k,1}, \ldots, \boldsymbol{\Sigma}_{k,d})\), respectively. Here, \(\boldsymbol{\Sigma}_{k,i}=\boldsymbol{\Phi}(t_{k,i},t_{k-1})\) is the STM at the \(i\)-th collocation point within the \(k\)-th interval \([t_{k-1},t_k]\). The inputs \(\boldsymbol{u}_k\) and algebraic variables \(\boldsymbol{z}_k\)—including ground reaction forces and the covariance matrices—are defined at the grid nodes only and held piecewise constant over each interval.

The resulting OCP has the following formulation: \[\tag{12} \begin{alignat}{3} \llap{\displaystyle\underset{\substack{\boldsymbol{\mu}_k,\boldsymbol{\xi}_k, \boldsymbol{u}_k, \\ \boldsymbol{\Phi}_k, \boldsymbol{\Sigma}_k, \boldsymbol{z}_k}}{\text{minimise}}}\, & & & J_k(\boldsymbol{\mu}_k,\boldsymbol{\xi}_k, \boldsymbol{u}_k) & & \tag{13} \\ \text{s.t.} \quad & \boldsymbol{0}& = & \; \boldsymbol{\Psi}^{\mu}_k(\boldsymbol{\mu}_{k-1},\boldsymbol{\mu}_k,\boldsymbol{\xi}_k, \boldsymbol{u}_k,\boldsymbol{z}_k), & \quad & k \in [1,N] \tag{14} \\ & \boldsymbol{\mu}_0 & = & \; \bar{\boldsymbol{\mu}}_0 & & \tag{15} \\ & \boldsymbol{0}& = & \; \boldsymbol{\Psi}^{\text{\Phi}}_k(\boldsymbol{\xi}_k,\boldsymbol{\Phi}_k,\boldsymbol{\Sigma}_k,\boldsymbol{u}_k,\boldsymbol{z}_k), & \quad & k \in [1,N] \tag{16} \\ & \boldsymbol{0}& = & \; \boldsymbol{\Omega}_k(\boldsymbol{\mu}_k,\boldsymbol{\xi}_k,\boldsymbol{\Phi}_k,\boldsymbol{\Sigma}_k, \boldsymbol{u}_k,\boldsymbol{z}_k), & \quad & k \in [0,N] \tag{17} \\ & 0 & \geq & \; h_i(\boldsymbol{\mu}_k, \boldsymbol{u}_k, \boldsymbol{z}_k) + \beta_i(\boldsymbol{\mu}_k, \boldsymbol{u}_k, \boldsymbol{z}_k), & \quad & k \in [1,N],\; i \in \mathcal{I}\quad (\lambda_{k,i} \geq 0) \tag{18} \end{alignat}\] The cost 13 encodes the minimum-lap-time objective. Equation 14 includes the collocation and continuity conditions for the mean dynamics 3 , initialised at \(\bar{\boldsymbol{\mu}}_0\) in 15 . Equation 16 encodes the collocation and continuity conditions for the dynamics of the STM 8 . Equation 17 collects the path equality constraints, which include the algebraic recovery of the covariance instances \(\boldsymbol{P}_k^j\) from \(\boldsymbol{\Phi}_k\) via 7 . Finally, 18 enforces the robust inequality constraints, with \(\mathcal{I}\) the set of indices of the inequality constraints, and \(\lambda_{k,i}\) the Lagrange multiplier associated with the \(i\)-th constraint at the \(k\)-th grid node. In our case, 18 includes the robust friction limit constraints defined in 11 . The back-off \(\beta_i\) is evaluated using the end-of-horizon covariance \(\boldsymbol{P}_k^H\), whose elements are included in \(\boldsymbol{z}_k\).

To recover the nominal—non-stochastic—MLTP formulation, it suffices to remove the STM dynamics from 12 and the back-off terms from 18 . With respect to the nominal MLTP problem, the robust formulation in 12 introduces the STMs as additional decision variables at all grid and collocation nodes—each STM instance adding \(n_x(n_x+n_p)\) scalar variables, since the lower \(n_p\) rows of each STM are structurally constant and equal to \([\boldsymbol{0}\mid \boldsymbol{I}_{n_p}]\), by an argument analogous to 5 –6 . Thus, the robust OCP incurs a considerably higher computational cost than the nominal MLTP formulation.

3.3 Parsimonious formulation↩︎

Starting from the nominal MLTP solution, we build a parsimonious variant of 12 that activates the robust machinery (state-transition matrices, covariance propagation, and constraint back-off) only at nodes where the nominal trajectory operates near the friction limit. Elsewhere, the back-off would tighten constraints already satisfied with ample margin, yielding no safety benefit. Confining robustness to the nodes at or near saturation therefore reduces the number of decision variables in the OCP, avoiding unnecessary STM and covariance instances and thereby lowering the computational cost relative to the full robust formulation.

Let \(\mathcal{C}=\{1,\dots,N\}\) be the set of grid nodes. Among the inequality constraints \(i\in\mathcal{I}\), we focus on the friction-limit constraints 9 , indexed by axle \(j\in\{1,2\}\) (front and rear). At each node \(k\), we read from the nominal solution two indicators of these constraints.

The first is the nominal Lagrange multiplier [28] \(\lambda_{k,j}^{(\mathrm{nom})}\ge 0\), which is positive only on the active set, i.e.,where the axle lies exactly on the friction boundary. At the NLP optimum, its magnitude is the shadow price of that constraint: the larger \(\lambda_{k,j}^{(\mathrm{nom})}\), the larger the lap-time reduction attainable by relaxing the corresponding friction limit; it therefore also ranks the saturated nodes according to how strongly they constrain performance.

The second is the constraint residual \({\bar{h}}_{k,j}=\bar{h}_j(\boldsymbol{\mu}_k,\boldsymbol{u}_k,\boldsymbol{z}_k)\), the function 9 evaluated at node \(k\) (hence the added subscript \(k\)): zero at saturation and negative inside the friction ellipse—the more negative, the larger the margin from the friction limit. We summarise each node by its worst axle, \[\label{eq:node95residual} {\bar{h}}_k = \max_{j\in\{1,2\}} {\bar{h}}_{k,j} ,\tag{19}\] so that \({\bar{h}}_k\) retains the larger residual of the front (\(j=1\)) and rear (\(j=2\)) axles, yielding a single residual per node.

We build the parsimonious set \(\mathcal{C}_{\mathrm{par}}\) in two steps. First, the critical (saturated) nodes are those carrying an active friction constraint, \[\label{eq:Csat} \mathcal{C}_{\mathrm{sat}}= \bigl\{\, k \in \mathcal{C}: \exists\, j\in\{1,2\},\; \lambda_{k,j}^{(\mathrm{nom})} > 0 \,\bigr\} ,\tag{20}\] where, in practice, “positive” means above a small numerical tolerance, since interior-point solvers return slightly positive multipliers. Having no residual margin, these nodes are robustified unconditionally.

Critical nodes alone, however, leave unprotected the nodes that operate near the friction limit without reaching it, where a disturbance could push the axle to saturation. We therefore introduce a small additional budget of near-critical nodes, sized as a fraction \(\rho\in(0,1)\) of the grid, \(n_\rho=\mathrm{ceil}(\rho N)\), where \(\mathrm{ceil}(\cdot)\) rounds up to the next integer so that the budget comprises at least one node; for instance, \(\rho=0.05\) (a \(5\%\) budget) on a grid of \(N=100\) nodes yields \(n_\rho=5\) nodes. This budget is allocated to the non-critical nodes closest to the limit, that is, the \(n_\rho\) nodes of \(\mathcal{C}\setminus\mathcal{C}_{\mathrm{sat}}\) with the largest (less negative) residual, \[\label{eq:Cthr} \mathcal{C}_{\mathrm{thr}}= \bigl\{\, k \in \mathcal{C}\setminus\mathcal{C}_{\mathrm{sat}}: {\bar{h}}_k \ge {\bar{h}}_{\mathrm{thr}} \,\bigr\} ,\tag{21}\] where \({\bar{h}}_{\mathrm{thr}}\) is the \(n_\rho\)-th largest residual among those nodes. The fraction \(\rho\) thus acts as a robustness budget, setting how many nodes, beyond the strictly saturated ones, are treated as robust.

The parsimonious OCP keeps the full structure of 12 but carries the STM and covariance variables, and the back-off, only on \[\label{eq:Cpar} \mathcal{C}_{\mathrm{par}}= \mathcal{C}_{\mathrm{sat}}\cup \mathcal{C}_{\mathrm{thr}}\subseteq \mathcal{C},\tag{22}\] reducing the number of decision variables relative to the full robust formulation, where they appear at every node. The result is a computationally lighter OCP that retains the solution quality of the full robust formulation, as Section 4 will detail.

4 Experiments setup↩︎

We validate the proposed framework through a simulation campaign. A model predictive controller (MPC) serves as a virtual driver, tracking MLTP references with different robustness settings over multiple closed-loop runs on a simulated FSAE vehicle; each run includes impulsive disturbances and parameter variations to probe reference driveability. Section 4.1 details the reference planning, Section 4.2 the MPC design, and Section 4.3 the testing strategy.

4.1 Reference trajectory planning↩︎

We apply the robust MLTP framework introduced in Section 3 to a FSAE vehicle driving on a representative sector of the Catalunya circuit (Figure 3). The entire circuit is parameterised by the curvilinear parameter \(\alpha \in [0,1]\), and the selected sector corresponds to the interval \(\left[0.70, 0.77\right]\). This sector includes two distinct corners: a low-speed turn for \(\alpha\in[0.72,0.73]\), and a high-speed turn for \(\alpha\in[0.75,0.76]\). We discretise the track section into \(N=140\) spatial intervals—corresponding to a spatial resolution \(\Delta s \approx 0.43\) m—and build a grid of \(N+1\) nodes for the OCP introduced in Problem 12 . The direct-collocation approach employs Gauss-Legendre collocation points and \(d=2\) as collocation degree.

We consider four MLTP references. The nominal reference (NOM), which serves as the baseline, solves the deterministic MLTP with no uncertainty on states or parameters; it is recoverable from 12 by removing the STM dynamics and back-off terms. ROB-S is robust with respect to the vehicle states only: the multiple short-horizon propagation scheme (Section 3.1) runs over all the grid nodes \(k\in[1,N]\), and the robust inequality constraints, i.e., with back-off terms, hold on the whole track sector. PAR-S retains state-only uncertainty but exploits the parsimonious formulation of Section 3.3, enforcing the robust constraints only on the parsimonious set \(\mathcal{C}_{\mathrm{par}}\) introduced in 22 ; hence, for each node \(k\in\mathcal{C}_{\mathrm{par}}\), only the STM instances over the \(H\) preceding steps and the resulting end-of-horizon covariance \(\boldsymbol{P}_k^H\) enter the OCP as decision variables. PAR-SP extends PAR-S to joint state-and-parameter uncertainty.

In ROB-S, PAR-S and PAR-SP we adopt \(H=5\) steps for the multiple short-horizon propagation—approximately \(2\) m of track, short enough to capture the disturbance evolution before a corrective driver response—and a confidence level \(p=0.90\) (i.e., \(\gamma=1.28\)), a moderately conservative design choice for the back-off terms. Each one of the short-horizons of Figure 2 is initialised at \(\boldsymbol{P}_k^0=\boldsymbol{P}_0\), which controls the uncertainty on states and parameters together with the covariance matrix \(\boldsymbol{Q}\) of the additive noise \(\boldsymbol{w}\) defined in 1 . Assuming independent initial uncertainties on states and parameters, we impose a diagonal \(\boldsymbol{P}_0 = \mathop{\mathrm{diag}}(\bar{\boldsymbol{\sigma}}_{\eta})^2\), with \(\bar{\boldsymbol{\sigma}}_{\eta}\) the standard deviations of the augmented state \(\boldsymbol{\eta}\) reported in Table 2; for comparison, the nominal (mean) parameter values are listed in Table 1 (Section 2.1). The additive noise \(\boldsymbol{w}\) models small, random disturbances to the vehicle dynamics—grip micro-variations, unmodelled aerodynamic effects—acting directly on the time derivatives of the velocity states \((u,v,r)\). We consequently set \(\boldsymbol{Q}= \mathop{\mathrm{diag}}([\bar{\boldsymbol{\sigma}}_w,\,\boldsymbol{0}])^2\) with \(\bar{\boldsymbol{\sigma}}_w = [0.055,\;0.032,\;0.05]\) (m/s\(^2\), m/s\(^2\), rad/s\(^2\)) and zero entries for all remaining states and parameters.

4pt

We transcribe the four OCPs in the CasADi–MATLAB environment [29], and solve the resulting nonlinear programs with IPOPT’s interior-point algorithm [30]. The nominal solution serves as warm start for the other three MLTPs and allows us to derive the critical and near-critical node sets \(\mathcal{C}_{\mathrm{sat}}\) and \(\mathcal{C}_{\mathrm{thr}}\) introduced in Section 3.3 for the parsimonious formulations used in PAR-S and PAR-SP. Figure 3 shows the critical nodes \(\mathcal{C}_{\mathrm{sat}}\) (orange) and the near-critical nodes \(\mathcal{C}_{\mathrm{thr}}\) (blue) computed from NOM, with the optimal trajectory (left panel), maximum Lagrange multiplier of friction limit constraints (top right panel), and maximum constraint residual (bottom right panel) plotted as a function of the curvilinear abscissa \(\alpha\). Before solving each robust OCP, we run an intermediate warm-start solve: starting from NOM, we augment the decision-variable set with the STM instances of 12 but leave the robust back-off terms inactive, so its solution both initialises the STM-related variables and supplies the converged starting point for the final, fully-constrained solve of ROB-S, PAR-S and PAR-SP. Because this solve only needs to seed the final problem, we run it at looser IPOPT tolerances, which is why it converges faster than NOM in Table 3 despite its larger problem size.

Figure 4 shows the front and rear axle saturation indices—\(S_1\) and \(S_2\), respectively—along the track curvilinear abscissa \(\alpha\), corresponding to the solution of NOM (solid black lines), ROB-S (solid red lines), PAR-S (dashed yellow lines), and PAR-SP (solid blue lines). We report in Table 3 the problem size, sparsity, IPOPT iterations and wall time for each MLTP, broken down into the warm-start solve and the final solve described above. The references are generated on a laptop with an Intel Core i7-12700H CPU at 2.30 GHz.

4pt

Returning to Figure 4, we can now observe the effect of the robust friction-limit constraints in ROB-S, PAR-S and PAR-SP: when the vehicle negotiates the low-speed corner at \(\alpha \in [0.71, 0.73]\), NOM reaches full saturation on the front and rear axles, while the three robust references preserve a margin from the saturation bound. To achieve this, all the robust MLTPs begin braking slightly earlier than NOM ahead of the corner entry, at \(\alpha \in [0.70, 0.71]\). The four references coincide for \(\alpha \in [0.735, 0.77]\), where the reduced tyre usage leaves the friction-limit constraint inactive. Interestingly, ROB-S and PAR-S show identical saturation profiles, yet PAR-S requires only \(109.8\) s against the \(251.3\) s of ROB-S (Table 3): the parsimonious formulation preserves solution accuracy at less than half the computational cost. Finally, PAR-SP introduces a slightly larger margin from full saturation than ROB-S and PAR-S, consistent with the additional parameter uncertainty it accounts for.

4.2 MPC design↩︎

Using an MPC as the virtual driver, rather than a human in a driving simulator, makes the closed-loop trials reproducible and allows a statistically meaningful number of runs. The MPC internal model is the same Single-Track model used in the MLTP (Section 2.1), with state \(\boldsymbol{x}\in\mathbb{R}^{n_x}\) and control \(\boldsymbol{u}\in\mathbb{R}^{n_u}\) of dimensions \(n_x=6\) and \(n_u=2\); the uncertain vehicle parameters are held at their nominal values (Table 1). At each sampling instant \(k\), the MPC predicts the vehicle trajectory and plans an optimal control sequence over a horizon of \(N_{\mathrm{pred}}=10\) stages spaced by \(h=0.01\) s.

The pairs \((\boldsymbol{x}_{k,i},\boldsymbol{u}_{k,i})\) of predicted states and planned controls—the first subscript \(k\) referring to the current sampling instant, the second \(i=0,\ldots,N_{\mathrm{pred}}\) to the prediction stage—constitute the set of decision variables of the MPC, while their target counterparts are the MLTP reference pairs \((\boldsymbol{x}_{k,i}^{\mathrm{ref}},\boldsymbol{u}_{k,i}^{\mathrm{ref}})\). To track the reference, the MPC solves an OCP with the quadratic cost \[\sum_{i=0}^{N_{\mathrm{pred}}-1} \Big( \|\boldsymbol{x}_{k,i}-\boldsymbol{x}_{k,i}^{\mathrm{ref}}\|_Q^2 + \|\boldsymbol{u}_{k,i}-\boldsymbol{u}_{k,i}^{\mathrm{ref}}\|_R^2 \Big) + \|\boldsymbol{x}_{k,N_{\mathrm{pred}}}-\boldsymbol{x}_{k,N_{\mathrm{pred}}}^{\mathrm{ref}}\|_{Q_N}^2, \label{eq:mpc95cost}\tag{23}\] where \(Q\), \(R\) and \(Q_N\) weight matrices. The OCP is subject to the internal-model dynamics and to the nominal friction-limit constraints 9 enforced along the horizon. Following the receding-horizon principle, only the first optimal control \(\boldsymbol{u}_{k,0}\) is applied to the simulated FSAE vehicle before the problem is solved again at \(k+1\) with updated state feedback and reference window. The problem is solved online with the Advanced-Step Real-Time Iteration (AS-RTI) scheme [31], [32]—an extension of the real-time iteration approach [33]—in the acados framework [34], performing a single warm-started iteration per instant rather than minimising 23 to convergence. This real-time-capable scheme would allow the same controller to be deployed on an actual vehicle in future work.

The core of the tracking scheme resides in how we assemble the targets \((\boldsymbol{x}_{k,i}^{\mathrm{ref}},\boldsymbol{u}_{k,i}^{\mathrm{ref}})\) online from the MLTP reference (Section 4.1), as shown in Figure 5. At sampling instant \(k\), the vehicle’s abscissa \(\alpha_k\) along the track synchronises the MPC horizon (dashed black lines) with the reference (dashed red lines): it locates the leading target \((\boldsymbol{x}_{k,0}^{\mathrm{ref}},\boldsymbol{u}_{k,0}^{\mathrm{ref}})\) and its reference time \(\tilde{t}_k\). From \(\tilde{t}_k\) we sample forwards the subsequent targets at the stage temporal spacing \(h\), \[\boldsymbol{x}^{\mathrm{ref}}_{k,i} = \boldsymbol{x}^{\mathrm{ref}}(\tilde{t}_k + i h), \qquad i = 0, \ldots, N_{\mathrm{pred}}, \label{eq:mpc95ref95sampling}\tag{24}\] and apply the same scheme for the controls. In summary, in Figure 5 the black car is the simulated FSAE vehicle, while the red vehicle is the ghost car being followed, corresponding to the MLTP reference. In our experiments the MPC does not search for a new optimal trajectory; it faithfully tracks the MLTP reference planned offline, which is the object of validation.

4.3 Testing strategy↩︎

We assess each reference with a Monte Carlo campaign of \(1000\) runs, the MPC of Section 4.2 driving the simulated FSAE vehicle along the planned reference (Section 4.1). The same controller drives all four references, so outcome differences reflect reference driveability alone. Each run applies three kinds of perturbations mirroring the uncertainty of the robust framework (Section 3): (i) a force–moment pulse realising the state uncertainty, (ii) a random scatter of the vehicle parameters, and (iii) the additive process noise \(\boldsymbol{w}\). The paired-sample design with a common seed makes the \(i\)-th run of every batch face identical disturbances.

The impulse (i) mirrors the model’s state uncertainty: we take the current state as the centre of a distribution with covariance \(\boldsymbol{P}_0\) (the MLTP design values), sample a random state jump, and inject the longitudinal/lateral force and yaw-moment pulses that produce it. We apply the pulses at \(\alpha_P=0.725\), midway through the low-speed corner, where both axles run close to the friction limit in the NOM reference (Figure 4). The pulses have duration \(T_p=0.1\) s, are centred at instant \(t_p\), and share a raised-cosine profile \[1+\cos\left[ \frac{2\pi(t-t_p)}{T_p}\right], \qquad t\in\left[t_p-\frac{T_p}{2},\;t_p+\frac{T_p}{2}\right], \label{eq:impulse95shape}\tag{25}\] scaled on each of its longitudinal-force, lateral-force and yaw-moment components so that the delivered impulse matches the target state change [35]. The parameter scatter (ii) draws the yaw inertia, centre-of-mass position and aerodynamic drag coefficient (Section 2.1) of the simulated vehicle from a Gaussian about their nominal values, with the standard deviations in the parameter block of \(\boldsymbol{P}_0\) (Table 2). The sampled parameters are held fixed along the run while the MPC internal model keeps the nominal values; only PAR-SP plans against this uncertainty, yet all four references are tested against it. The process noise (iii) draws \(\boldsymbol{w}(t_k)\) at every step from the planning covariance \(\boldsymbol{Q}\) of the stochastic dynamics 1 .

A run is classified as physically failed if any plausibility bound is violated: body sideslip greater than \(0.5\,\mathrm{rad}\), yaw rate greater than \(3.0\,\mathrm{rad/s}\), or either-axle slip angle greater than \(0.3\,\mathrm{rad}\). A run that does not trigger this physical-failure criterion reaches the end of the track sector and is therefore classified as completed. A completed run is further classified as survived if the saturation ratio \(S_j\) never exceeds \(0.999\) on either axle for a continuous time window longer than the dwell threshold \(T_{\mathrm{thr}}\). Hence, for each reference, the batch of \(1000\) runs gives rise to a completed cohort and a survived subcohort. In Section 5, we analyse how \(T_{\mathrm{thr}}\) influences the resulting survival counts, i.e., the number of runs in the survived subcohort.

5 Results and discussion↩︎

5.1 Driveability analysis↩︎

We assess the driveability of each MLTP reference (Section 4.1) by means of the fraction of survived runs (Section 4.3) out of the \(1000\) perturbed runs the MPC virtual driver (Section 4.2) drives tracking that reference. The paired campaign applies identical disturbances to all four references, so differences in failure rate are attributable to the references alone. A run is declared failed based on the criteria introduced in Section 4.3; we record the track position where each failed run first meets one of these criteria. For completed but failed runs, the failure location corresponds to the end of the \(T_{\mathrm{thr}}\) adherence-dwell window.

Figure 6 reports the fraction of survived runs against the track abscissa \(\alpha\), with the vertical line at \(\alpha=\alpha_P\) indicating the location of the impulsive disturbance actions. For \(T_{\mathrm{thr}}=0.10\) s the sector-end survival counts are NOM \(111/1000\) (\(11.1\%\), black line), PAR-S \(345/1000\) (\(34.5\%\), yellow line), and PAR-SP \(585/1000\) (\(58.5\%\), blue line); the shaded regions show how a different choice of \(T_{\mathrm{thr}}\) would affect the survival counts, with the upper bound of each band corresponding to the more permissive \(T_{\mathrm{thr}}=0.15\) s and the lower bound to the stricter threshold \(T_{\mathrm{thr}}=0.05\) s. ROB-S is omitted since its counts coincide with PAR-S. Completion, by contrast, is nearly universal: ROB-S, PAR-S and PAR-SP complete every one of the \(1000\) realisations, while NOM completes \(984/1000\), the \(16\) shortfall runs leaving their plausibility bound before reaching the sector end. The common completed cohort (Section 4.3) therefore coincides with NOM’s own completed set, \(984/1000\), and underpins the analyses of Sections 5.1 and 5.2.

To visualise the different behaviour of the survived and failed runs, we plot the trajectories of a representative sample of the Monte Carlo campaign for the NOM and PAR-SP (Figure 7 (a)–(b)). Each panel draws \(100\) runs independently from that reference’s own campaign, in proportion to its end-section survival rate, so the blue (survived) to red (failed) ratio mirrors the counts above.

Figure 8 illustrates the adherence behaviour of the whole batch relative to each reference. Here, we analyse the dispersion of the axle saturation indices over the common completed cohort of \(984/1000\) runs, so that every band spans the full track sector for the same paired-seed runs. At each abscissa \(\alpha\) we take the 10th, 25th, 50th, 75th and 90th percentiles of \(S_1\) and \(S_2\) (10 ) and plot them for NOM (black) and PAR-SP (blue): solid lines are the medians, the darker bands the 25th–75th percentile range, the lighter bands the 10th–90th. For clarity we omit ROB-S and PAR-S, whose distributions nearly overlap that of PAR-SP.

The adherence profiles in Figure 8 reflect the MLTP reference curves of Figure 4, yet the bands reveal how some runs tracking a robust reference may still hit tyre saturation under certain disturbance realisations. For instance, although the PAR-SP reference keeps a margin to the friction limit, the 25th–75th percentile band of \(S_2\) for PAR-SP reaches saturation over \(\alpha\in[0.70,0.71]\) during braking into the slow-speed corner; this explains the survival-rate drop at the same \(\alpha\) range in Figure 6. Each remaining survival-rate drop in Figure 6 corresponds to one of the percentile bands reaching the friction limit. More trajectory samples and further telemetry signals are available in the supplementary material [23].

5.2 Sector-time cost and steering effort↩︎

Robustness carries a cost: the back-off terms that hold the planned trajectory off the friction boundary entail a higher sector time. Figure 9 (a) shows the sector-time distribution over the common completed cohort (Section 5.1). For each reference, the boxplot includes the median (thick line), the 25th–75th percentile box, the 10th–90th whiskers, while the dot indicates the planned MLTP sector time. As expected, median sector time rises monotonically with robustness: NOM \(10.50\) s, ROB-S and PAR-S \(10.62\) s (\(+120\) ms over NOM), and PAR-SP \(10.67\) s (\(+170\) ms over NOM, \(+50\) ms over ROB-S/PAR-S). ROB-S and PAR-S share identical distributions because they plan the same reference (Figure 4).

On the other hand, the steering effort \(\int \dot{\delta}^2\,\mathrm{d}t\)—a measure of the overall steering activity required to follow the reference—decreases with robustness, as shown in Figure 9 (b). The boxplots describe the steering effort distribution for each reference, while the dot indicates the planned MLTP effort. The plot is cropped to \([1,3]\times10^{-5}\) for readability; NOM’s 90th-percentile whisker extends to \(6.6\times10^{-5}\), beyond the plotted range.

5.3 Discussion↩︎

The campaign confirms that the robust MLTP framework (Section 3) produces trajectory and control references that keep the vehicle away from tyre saturation, as Figure 6 shows. Although the exact number of survivors depends on the dwell window \(T_{\mathrm{thr}}\) over which saturation is tolerated, the shaded band shows a consistent trend: accounting for state uncertainty (ROB-S/PAR-S) raises the survived fraction over NOM, and adding parameter uncertainty (PAR-SP) raises it further still.

The representative samples of Figure 7 illustrate this difference. The NOM runs deviate markedly off-line over \(\alpha\in[0.73,0.75]\) and spread over a wider corridor than the robust PAR-SP runs, signalling a partial loss of vehicle control and a higher sensitivity to the applied disturbances.

Since the underlying problem is one of minimum-time planning, this robustness has a natural counterpart in the sector time (Figure 9 (a)): moving from NOM to ROB-S/PAR-S costs a median \(120\) ms, while accounting also for parameter uncertainty (PAR-SP) adds a further \(50\) ms. In all four cases the MPC driver does not, on average, reproduce the planned MLTP sector time (dots), settling at higher median values (thick lines); this gap, however, narrows as robustness increases—for PAR-SP, more than \(10\%\) of runs beat the planned time, with the MLTP dot lying above the 10th-percentile whisker.

The steering-effort panel (Figure 9 (b)) further measures driveability—understood here as the ease with which a driver reproduces the MLTP reference. The NOM runs require repeated steering corrections that skew their distribution asymmetrically towards high values, with a 90th percentile markedly larger than those of the three robust references, which remain close to one another at a far lower effort.

Although the campaign relies on an MPC virtual driver, the higher driveability of the robust references—greater survival and lower steering effort—plausibly carries over to a human driver, a hypothesis that further testing on a driving simulator or on track could confirm.

6 Conclusion↩︎

The minimum-time racing line is, by construction, fragile: it runs along the friction-limit boundary with no margin to spare, so a small disturbance, be it a gust, a kerb strike, a momentary loss of grip or a drift in the vehicle parameters, is enough to make the nominal MLTP solution hard or unsafe to follow. To plan a reference that keeps a conservative margin exactly where it is needed, we presented and validated a parsimonious disturbance-aware framework for minimum-lap-time planning under parametric uncertainty. It extends a prior disturbance-aware MLTP framework along three lines: parameter uncertainty in the covariance propagation, a parsimonious activation strategy for the robust constraints, and validation of the planned references by an MPC virtual driver in a Monte Carlo campaign.

Targeting tyre saturation as the critical condition, we enforced robust friction-limit constraints in the OCP and confined them, through a spatially selective activation strategy, to the circuit segments where they matter most (only \(34.3\%\) of the track grid nodes in the tested case). This parsimony costs no accuracy: under state uncertainty alone, the parsimonious reference (PAR-S) reproduces the same trajectory and controls as its full counterpart (ROB-S), in which the robust constraints hold at every node, yet PAR-S solves in \(109.8\) s against \(251.3\) s. This efficiency keeps the formulation tractable once parameter uncertainty enters (PAR-SP).

Validation employed MPC as a virtual test driver: for each of the four planned references (NOM, ROB-S, PAR-S, PAR-SP), the same controller drove a simulated FSAE vehicle over 1000 runs on the Catalunya sector, under randomly drawn impulsive disturbances, scattered vehicle parameters and process noise. With a run counted as failed once it remains at the friction limit longer than the dwell threshold \(T_{\mathrm{thr}}=0.10\) s, the survived fraction rose from \(11.1\%\) (NOM) to \(34.5\%\) (ROB-S/PAR-S) and \(58.5\%\) (PAR-SP), a roughly fivefold safety gain for the full state-and-parameter formulation. The cost of this safety was small: a median of \(120\) ms (ROB-S/PAR-S) and \(170\) ms (PAR-SP) over a \(10.50\) s nominal sector, alongside tighter trajectory dispersion and lower steering effort.

These results rest on a single-track vehicle model, one track sector and a single disturbance location, with the MPC standing in for the driver; within these bounds, the higher driveability of the robust references (greater survival, tighter dispersion, lower steering effort) plausibly carries over to a human driver. The natural next steps are to confirm this on a driving simulator or on track, and to extend the framework to more detailed vehicle models while containing its computational cost.

Disclosure statement↩︎

No potential conflict of interest was reported by the author(s).

Data availability statement↩︎

The full parameter set of the vehicle model used in this paper, together with an interactive dashboard illustrating the closed-loop Monte Carlo results, are provided as supplementary material and openly available at [23]. The dashboard can also be accessed directly at https://martinogulisano.github.io/data-parsimonious-disturbance-aware/.

References↩︎