Indirect Shooting Method

Updated 24 April 2026

Indirect shooting method is a numerical strategy that recasts optimal control problems into boundary value problems using the Pontryagin Maximum Principle.
It employs root-finding techniques and Newton iterations, enhanced by adjoint propagation and continuation strategies, to achieve high-precision solutions.
Extensions to stochastic, delayed, and hybrid systems showcase its practical applications in aerospace, robotics, and trajectory optimization.

An indirect shooting method is a numerical strategy for solving optimal control and two-point boundary value problems in systems governed by differential equations. It reformulates the first-order necessary conditions, typically derived from the Pontryagin Maximum Principle (PMP), as a root-finding problem over a transformed set of unknowns, usually the initial values of costate (adjoint) variables and occasionally additional multipliers or the final time. The method is termed "indirect" because it leverages necessary conditions derived analytically—contrasting with "direct" approaches based solely on discrete parametrizations of controls. Indirect shooting has become a cornerstone of modern computational optimal control in deterministic, stochastic, and delay systems and is foundational for high-precision guidance, trajectory optimization, and robotics applications.

1. Mathematical Foundations and Formulation

The indirect shooting method arises from the recasting of an optimal control problem, such as

$\min_{x(\cdot),\,u(\cdot)} J[x, u] = \phi(x(T)) + \int_0^T L(x,u,t)\,dt, \quad \text{s.t. } \dot{x} = f(x,u,t),\; x(0)=x_0,$

into a boundary-value problem for the Hamiltonian dynamics derived from the PMP: $\begin{cases} \dot{x}(t) = f(x,u,t), \ -\dot{\lambda}(t) = \frac{\partial H}{\partial x}, \ 0 = \frac{\partial H}{\partial u}, \ x(0) = x_0,\quad \lambda(T) = \phi_x(x(T)), \end{cases}$ where $H(x,u,\lambda,t) = L(x,u,t)+\lambda^\top f(x,u,t)$ and $\lambda(t)$ are costates.

The core numerical idea is to select a parameterization of the unknown initial adjoint(s) and any other unspecified conditions (final time, multipliers), propagate the coupled system (states forward, costates backward or jointly), and define a residual function capturing the terminal (or target) constraints. Root-finding (typically via Newton–Kantorovich iterations) targets zeros of this residual, constituting solutions of the original boundary-value problem (Scheiber, 2022, Ross, 2020).

2. Algorithmic Structure and Adjoint System

Given the reformulation above, the method unfolds algorithmically as follows:

Choose an initial guess for all unknowns (e.g., initial adjoints and, when relevant, free final time or multipliers).
Numerically integrate the forward ODEs (state equations), using the current initial adjoint guess, and backward integrate the adjoint ODEs (or jointly propagate them).
Evaluate the terminal residual, i.e., mismatch between achieved and desired terminal/periodic/constraint conditions.
Compute the Jacobian of the residual map. This is often performed using the adjoint system itself, yielding expressions for the Jacobian (e.g., via variational equations or inner-product identities) as documented in (Scheiber, 2022).
Apply a Newton-type update to refine the initial adjoint (and any other unknowns).
Iterate until the residual norm is within tolerance.

The adjunct system's structure allows for substantial analytic simplifications, including costate normalization (to reduce dimensionality), elimination of redundant multipliers (by exploiting homogeneity), and admissible restriction of the shooting domain (e.g., to a positive octant for landing problems, as in (Wang et al., 2024)).

3. Continuation and Homotopy Strategies

A key source of robustness in practical indirect shooting is the use of numerical continuation (homotopy) procedures. Initialization is challenging, especially for stiff or highly nonlinear problems with poor guess robustness. Continuation-based strategies mitigate this by gradually morphing an analytically tractable "order-zero" problem (with known solution) toward the original problem through a family of parametrized intermediate problems (Bonalli et al., 2017, Bonalli et al., 2017, Bonalli et al., 2017, Raff et al., 2024). At each continuation step, the previous solution is used to initialize the next, with corrections applied via Newton iterations.

Applications of this framework include:

Stepwise introduction of physical complexity (e.g., adding gravity, thrust, curvature).
Gradual inclusion of state/control delays, ensuring convergence of the adjoint fields at each stage (Bonalli et al., 2017).
Systematic variation of parameters (e.g., slope or speed in bipedal gait generation (Raff et al., 2024)).
Homotopy traversal in the space of terminal constraints or cost functions.

The existence and convergence of such strategies are underpinned by the continuous dependence of optimal trajectories and multipliers on the homotopy parameter, proven under suitable regularity conditions (Bonalli et al., 2017).

4. Extensions: Stochastic, Delayed, and Hybrid Systems

The indirect shooting method paradigm extends beyond deterministic ODEs:

For stochastic optimal control problems governed by rough differential equations (RDEs), an indirect shooting formulation based on a stochastic version of the PMP (backward RDEs for the adjoint, sample-wise integration over Brownian path realizations) yields a similar root-finding approach, accelerated by handling the whole multipath ensemble in block (Lew, 10 Feb 2025).
For control systems with distributed delays in state or input, the boundary value formulation incorporates advanced/retarded terms in the adjoint ODEs, and the shooting vector includes both adjoint and history terms (Bonalli et al., 2017).
For hybrid/differential-algebraic and periodic problems (e.g., trajectory planning for periodic walking), the shooting vector concatenates all quantities that affect the periodicity and event constraints; Lagrange multipliers and costates associated with reset maps or event surfaces are included explicitly, and root-finding targets all stationarity and transversality conditions simultaneously (Raff et al., 2024).

These extensions preserve the structure of the indirect approach while adapting to the non-classical forward-backward dependencies that arise in stochastic and delay settings.

5. Numerical Characteristics, Conditioning, and Comparison to Direct Strategies

Indirect shooting delivers high-precision solutions, exploiting the full structure of the first-order optimality conditions—particularly the costate feedback relations between control and states. This results in significantly faster or more accurate solutions than direct methods under favorable initialization and regularity.

Key points regarding computational and conditioning properties:

The Newton–Kantorovich iterations for the shooting function show local quadratic convergence under Jacobian nonsingularity and proximity to the solution (Scheiber, 2022).
Condition numbers of the residual Jacobian remain moderate even as the number of collocation points or state dimensions increases, provided the problem is regular and homotopy steps are managed adaptively (Raff et al., 2024).
A central insight is the exact equivalence—at the level of first-order necessary conditions—between direct and indirect shooting methods. In discretized settings, the gradient of the reduced terminal cost with respect to control parameters incorporates the discrete adjoint, and setting this gradient to zero recovers the same equations as the indirect method (Ross, 2020). However, naive application of direct methods using black-box NLP solvers can fail to exploit this structure, leading to ill-conditioning or poor convergence when the underlying Hamiltonian structure is not exposed.
Implementation of Hamiltonian programming (explicitly embedding adjoint propagation and stationarity) improves numerical robustness and aligns with modern automatic differentiation and machine learning frameworks (Ross, 2020).

Empirical studies confirm indirect shooting’s superior speed and accuracy relative to direct methods, especially in challenging parameter regimes or high-precision applications (Bonalli et al., 2017, Lew, 10 Feb 2025).

6. Practical Applications and Empirical Performance

Representative applications include:

Optimal guidance of endo-atmospheric launch vehicles under geometric, mixed control-state constraints, using chart-based formulation and adaptive continuation for high-accuracy interception trajectories (Bonalli et al., 2017, Bonalli et al., 2017).
Physics-informed trajectory optimization for lunar/planetary soft landing, where clever normalization of the adjoint guess dramatically shrinks the shooting domain, ensuring fast and robust convergence (Wang et al., 2024).
Generation of libraries of optimal periodic trajectories for dynamic bipedal robots, where indirect shooting paired with continuation handles parameter sweeps, bifurcation points, and reconstruction of multipliers from passive gaits with high efficiency (Raff et al., 2024).
Nonlinear stochastic optimal control problems, where the indirect strategy achieves order-of-magnitude speedups over direct approaches by solving for a compact set of multipath initial adjoints (Lew, 10 Feb 2025).
Delay-ridden nonlinear control problems, where delay continuation admits robust solutions even for substantial delay magnitudes (Bonalli et al., 2017).

Numerical results report computation times typically in the sub-second to few-second range for complex aerospace or robotics applications, high convergence rates (99%+ in some large-scale Monte Carlo scenarios), and improved regions of convergence versus direct approaches.

7. Limitations, Assumptions, and Ongoing Developments

While indirect shooting is provably effective under regularity and with robust initialization, it remains sensitive to:

The quality of the initial guess, particularly in high-dimensional or stiff systems without a suitable order-zero or analytic starting point.
Nonregular (singular) arcs, though empirical evidence suggests these are infrequent in many practical problems (Bonalli et al., 2017).
The conditioning of the shooting function Jacobian near singularities or folds in the solution manifold, where adaptive continuation and chart switching are essential (Bonalli et al., 2017, Raff et al., 2024).
The extension to large-scale systems or strongly nonsmooth problems, where discretization or hybrid-analytic–direct approaches may be necessary.

Ongoing developments focus on improved analytic initialization, systematic use of homotopy, embedding into structured NLP frameworks, and integration with modern adjoint- and autodiff-based computational platforms, reflecting the strong geometric and algorithmic links between optimal control, machine learning, and dynamical systems (Ross, 2020).

References: