LGL Costate Estimation and General Pseudospectral Transcription
This section completes the treatment of costate estimation for the three principal Legendre pseudospectral schemes:
Legendre–Gauss (LG),
Legendre–Gauss–Radau (LGR), and
Legendre–Gauss–Lobatto (LGL).
The first part derives the LGL costate equations from the KKT conditions of the LGL nonlinear program. The second part compares the LG, LGR, and LGL costate mappings. The final part extends the pseudospectral transcription from the simplified Mayer problem on [ − 1 , 1 ] [-1,1] [ − 1 , 1 ] [ t 0 , t f ] [t_0,t_f] [ t 0 , t f ]
Legendre–Gauss–Lobatto Grid ¶ Let
τ 1 , … , τ N \tau_1,\ldots,\tau_N τ 1 , … , τ N be the LGL points. Unlike LG and LGR grids, the LGL grid includes both endpoints:
τ 1 = − 1 , τ N = 1. \boxed{
\tau_1=-1,
\qquad
\tau_N=1.
} τ 1 = − 1 , τ N = 1. The interior points are roots of the derivative of an appropriate Legendre polynomial.
LGL State Approximation ¶ Approximate the state by
X ( τ ) = ∑ i = 1 N X i L i ( τ ) , \boxed{
\boldsymbol{X}(\tau)
=
\sum_{i=1}^{N}\boldsymbol{X}_iL_i(\tau),
} X ( τ ) = i = 1 ∑ N X i L i ( τ ) , where
L i ( τ ) = ∏ j = 1 j ≠ i N τ − τ j τ i − τ j . L_i(\tau)
=
\prod_{\substack{j=1\\j\neq i}}^{N}
\frac{\tau-\tau_j}{\tau_i-\tau_j}. L i ( τ ) = j = 1 j = i ∏ N τ i − τ j τ − τ j . This interpolation polynomial has degree at most N − 1 N-1 N − 1
Its derivative is
X ˙ ( τ ) = ∑ i = 1 N X i L i ′ ( τ ) , \dot{\boldsymbol{X}}(\tau)
=
\sum_{i=1}^{N}
\boldsymbol{X}_iL_i'(\tau), X ˙ ( τ ) = i = 1 ∑ N X i L i ′ ( τ ) , which has degree at most N − 2 N-2 N − 2
LGL Differentiation Matrix ¶ Define
D j i = L i ′ ( τ j ) , i , j = 1 , … , N . D_{ji}
=
L_i'(\tau_j),
\qquad
i,j=1,\ldots,N. D ji = L i ′ ( τ j ) , i , j = 1 , … , N . Thus,
D ∈ R N × N . \boxed{
\boldsymbol{D}\in\mathbb{R}^{N\times N}.
} D ∈ R N × N . Collocating the dynamics at all LGL points gives
∑ i = 1 N D j i X i = f ( X j , U j ) , j = 1 , … , N . \sum_{i=1}^{N}
D_{ji}\boldsymbol{X}_i
=
\boldsymbol{f}(\boldsymbol{X}_j,\boldsymbol{U}_j),
\qquad
j=1,\ldots,N. i = 1 ∑ N D ji X i = f ( X j , U j ) , j = 1 , … , N . In matrix form,
D X 1 : N = F 1 : N . \boxed{
\boldsymbol{D}\boldsymbol{X}_{1:N}
=
\boldsymbol{F}_{1:N}.
} D X 1 : N = F 1 : N . Over-Collocation and Singularity ¶ The derivative polynomial has degree at most N − 2 N-2 N − 2 N − 1 N-1 N − 1 N N N
As a result,
D is singular . \boxed{
\boldsymbol{D}\ \text{is singular}.
} D is singular . This fact has important consequences for costate estimation and for the relationship between differential and integral pseudospectral formulations.
Simplified Mayer Problem ¶ Consider
min x ( ⋅ ) , u ( ⋅ ) Φ ( x ( 1 ) ) , \min_{\boldsymbol{x}(\cdot),\boldsymbol{u}(\cdot)}
\Phi\bigl(\boldsymbol{x}(1)\bigr), x ( ⋅ ) , u ( ⋅ ) min Φ ( x ( 1 ) ) , subject to
x ˙ ( τ ) = f ( x ( τ ) , u ( τ ) ) , τ ∈ [ − 1 , 1 ] , \dot{\boldsymbol{x}}(\tau)
=
\boldsymbol{f}\bigl(\boldsymbol{x}(\tau),\boldsymbol{u}(\tau)\bigr),
\qquad
\tau\in[-1,1], x ˙ ( τ ) = f ( x ( τ ) , u ( τ ) ) , τ ∈ [ − 1 , 1 ] , and
x ( − 1 ) = x 0 . \boldsymbol{x}(-1)=\boldsymbol{x}_0. x ( − 1 ) = x 0 . The Hamiltonian is
H ( x , u , λ ) = λ T f ( x , u ) . H(\boldsymbol{x},\boldsymbol{u},\boldsymbol{\lambda})
=
\boldsymbol{\lambda}^{\mathsf{T}}\boldsymbol{f}(\boldsymbol{x},\boldsymbol{u}). H ( x , u , λ ) = λ T f ( x , u ) . The continuous optimality conditions are
λ ( 1 ) = ∇ x Φ ( x ( 1 ) ) , \boldsymbol{\lambda}(1)
=
\nabla_{\boldsymbol{x}}\Phi\bigl(\boldsymbol{x}(1)\bigr), λ ( 1 ) = ∇ x Φ ( x ( 1 ) ) , λ ˙ = − ∇ x ( λ T f ) , \dot{\boldsymbol{\lambda}}
=
-
\nabla_{\boldsymbol{x}}
\left(
\boldsymbol{\lambda}^{\mathsf{T}}\boldsymbol{f}
\right), λ ˙ = − ∇ x ( λ T f ) , and
∇ u ( λ T f ) = 0 . \nabla_{\boldsymbol{u}}
\left(
\boldsymbol{\lambda}^{\mathsf{T}}\boldsymbol{f}
\right)
=
\boldsymbol{0}. ∇ u ( λ T f ) = 0 . LGL Nonlinear Program ¶ The LGL transcription is
min X , U Φ ( X N ) subject to F 1 : N − D X 1 : N = 0 , x 0 − X 1 = 0 . \boxed{
\begin{aligned}
\min_{\boldsymbol{X},\boldsymbol{U}}\quad
&
\Phi(\boldsymbol{X}_N)
\\
\text{subject to}\quad
&
\boldsymbol{F}_{1:N}
-
\boldsymbol{D}\boldsymbol{X}_{1:N}
=
\boldsymbol{0},
\\
&
\boldsymbol{x}_0-\boldsymbol{X}_1
=
\boldsymbol{0}.
\end{aligned}
} X , U min subject to Φ ( X N ) F 1 : N − D X 1 : N = 0 , x 0 − X 1 = 0 . The NLP variables are
X 1 , … , X N , U 1 , … , U N . \boldsymbol{X}_1,\ldots,\boldsymbol{X}_N,
\qquad
\boldsymbol{U}_1,\ldots,\boldsymbol{U}_N. X 1 , … , X N , U 1 , … , U N . LGL NLP Lagrangian ¶ Let
Λ 1 : N = [ Λ 1 ⋮ Λ N ] \boldsymbol{\Lambda}_{1:N}
=
\begin{bmatrix}
\boldsymbol{\Lambda}_1\\
\vdots\\
\boldsymbol{\Lambda}_N
\end{bmatrix} Λ 1 : N = ⎣ ⎡ Λ 1 ⋮ Λ N ⎦ ⎤ denote the multipliers associated with the collocation equations.
Let
denote the multiplier associated with the initial condition.
The NLP Lagrangian is
L = Φ ( X N ) + ⟨ Λ 1 : N , F 1 : N − D X 1 : N ⟩ + ⟨ μ , x 0 − X 1 ⟩ . \begin{aligned}
\mathcal{L}
={}&
\Phi(\boldsymbol{X}_N)
+
\left\langle
\boldsymbol{\Lambda}_{1:N},
\boldsymbol{F}_{1:N}
-
\boldsymbol{D}\boldsymbol{X}_{1:N}
\right\rangle
\nonumber\\
&+
\left\langle
\boldsymbol{\mu},
\boldsymbol{x}_0-\boldsymbol{X}_1
\right\rangle.
\end{aligned} L = Φ ( X N ) + ⟨ Λ 1 : N , F 1 : N − D X 1 : N ⟩ + ⟨ μ , x 0 − X 1 ⟩ . KKT Conditions ¶ Because both endpoints are collocation points, the initial and terminal state equations must be treated separately from the interior equations.
Initial Point ¶ Stationarity with respect to X 1 \boldsymbol{X}_1 X 1
∇ X [ Λ 1 T f ( X 1 , U 1 ) ] − d 1 T Λ 1 : N = μ , \boxed{
\nabla_{\boldsymbol{X}}
\left[
\boldsymbol{\Lambda}_1^{\mathsf{T}}
\boldsymbol{f}(\boldsymbol{X}_1,\boldsymbol{U}_1)
\right]
-
\boldsymbol{d}_1^{\mathsf{T}}\boldsymbol{\Lambda}_{1:N}
=
\boldsymbol{\mu},
} ∇ X [ Λ 1 T f ( X 1 , U 1 ) ] − d 1 T Λ 1 : N = μ , where d 1 \boldsymbol{d}_1 d 1 D \boldsymbol{D} D
Terminal Point ¶ Stationarity with respect to X N \boldsymbol{X}_N X N
d N T Λ 1 : N − ∇ X [ Λ N T f ( X N , U N ) ] = ∇ X Φ ( X N ) . \boxed{
\boldsymbol{d}_N^{\mathsf{T}}\boldsymbol{\Lambda}_{1:N}
-
\nabla_{\boldsymbol{X}}
\left[
\boldsymbol{\Lambda}_N^{\mathsf{T}}
\boldsymbol{f}(\boldsymbol{X}_N,\boldsymbol{U}_N)
\right]
=
\nabla_{\boldsymbol{X}}\Phi(\boldsymbol{X}_N).
} d N T Λ 1 : N − ∇ X [ Λ N T f ( X N , U N ) ] = ∇ X Φ ( X N ) . Interior Points ¶ For
j = 2 , … , N − 1 , j=2,\ldots,N-1, j = 2 , … , N − 1 , stationarity gives
∇ X [ Λ j T f ( X j , U j ) ] = d j T Λ 1 : N , \boxed{
\nabla_{\boldsymbol{X}}
\left[
\boldsymbol{\Lambda}_j^{\mathsf{T}}
\boldsymbol{f}(\boldsymbol{X}_j,\boldsymbol{U}_j)
\right]
=
\boldsymbol{d}_j^{\mathsf{T}}\boldsymbol{\Lambda}_{1:N},
} ∇ X [ Λ j T f ( X j , U j ) ] = d j T Λ 1 : N , where d j \boldsymbol{d}_j d j j j j D \boldsymbol{D} D
Control Stationarity ¶ For all LGL points,
∇ U ⟨ Λ 1 : N , F 1 : N ⟩ = 0 . \boxed{
\nabla_{\boldsymbol{U}}
\left\langle
\boldsymbol{\Lambda}_{1:N},
\boldsymbol{F}_{1:N}
\right\rangle
=
\boldsymbol{0}.
} ∇ U ⟨ Λ 1 : N , F 1 : N ⟩ = 0 . LGL Costate Mapping ¶ Let
W = diag ( w 1 , … , w N ) , \boldsymbol{W}
=
\operatorname{diag}(w_1,\ldots,w_N), W = diag ( w 1 , … , w N ) , where w i w_i w i
Define the costate estimate by
λ 1 : N = W − 1 Λ 1 : N . \boxed{
\boldsymbol{\lambda}_{1:N}
=
\boldsymbol{W}^{-1}\boldsymbol{\Lambda}_{1:N}.
} λ 1 : N = W − 1 Λ 1 : N . Componentwise,
λ i = Λ i w i . \boldsymbol{\lambda}_i
=
\frac{\boldsymbol{\Lambda}_i}{w_i}. λ i = w i Λ i . This is algebraically the simplest of the three mappings.
LGL Adjoint Differentiation Matrix ¶ Define
D † ∈ R N × N \boldsymbol{D}^\dagger\in\mathbb{R}^{N\times N} D † ∈ R N × N such that
D i i † = D i i D_{ii}^\dagger=D_{ii} D ii † = D ii for interior diagonal entries, while the endpoint diagonal entries are modified:
D 11 † = − D 11 − 1 w 1 , \boxed{
D_{11}^\dagger
=
-D_{11}-\frac{1}{w_1},
} D 11 † = − D 11 − w 1 1 , and
D N N † = − D N N + 1 w N . \boxed{
D_{NN}^\dagger
=
-D_{NN}+\frac{1}{w_N}.
} D NN † = − D NN + w N 1 . For
define
D i j † = − w j w i D j i . \boxed{
D_{ij}^\dagger
=
-
\frac{w_j}{w_i}D_{ji}.
} D ij † = − w i w j D ji . Compact Discrete Adjoint Equation ¶ Let
e 1 = [ 1 0 ⋯ 0 ] T , e N = [ 0 ⋯ 0 1 ] T . \boldsymbol{e}_1
=
\begin{bmatrix}
1&0&\cdots&0
\end{bmatrix}^{\mathsf{T}},
\qquad
\boldsymbol{e}_N
=
\begin{bmatrix}
0&\cdots&0&1
\end{bmatrix}^{\mathsf{T}}. e 1 = [ 1 0 ⋯ 0 ] T , e N = [ 0 ⋯ 0 1 ] T . Then the transformed KKT system can be written as
D † λ 1 : N = − ∇ X ⟨ λ 1 : N , F 1 : N ⟩ + e 1 w 1 ( μ − λ 1 ) + e N w N [ λ N − ∇ X Φ ( X N ) ] . \boxed{
\begin{aligned}
\boldsymbol{D}^\dagger\boldsymbol{\lambda}_{1:N}
={}&
-
\nabla_{\boldsymbol{X}}
\left\langle
\boldsymbol{\lambda}_{1:N},
\boldsymbol{F}_{1:N}
\right\rangle
\\
&+
\frac{\boldsymbol{e}_1}{w_1}
\left(
\boldsymbol{\mu}-\boldsymbol{\lambda}_1
\right)
+
\frac{\boldsymbol{e}_N}{w_N}
\left[
\boldsymbol{\lambda}_N-
\nabla_{\boldsymbol{X}}\Phi(\boldsymbol{X}_N)
\right].
\end{aligned}
} D † λ 1 : N = − ∇ X ⟨ λ 1 : N , F 1 : N ⟩ + w 1 e 1 ( μ − λ 1 ) + w N e N [ λ N − ∇ X Φ ( X N ) ] . The control stationarity equation becomes
∇ U ⟨ λ 1 : N , F 1 : N ⟩ = 0 . \boxed{
\nabla_{\boldsymbol{U}}
\left\langle
\boldsymbol{\lambda}_{1:N},
\boldsymbol{F}_{1:N}
\right\rangle
=
\boldsymbol{0}.
} ∇ U ⟨ λ 1 : N , F 1 : N ⟩ = 0 . Coupling of Boundary Conditions and Costate Dynamics ¶ For a pure discrete adjoint equation, one would prefer
D † λ 1 : N = − ∇ X ⟨ λ 1 : N , F 1 : N ⟩ . \boldsymbol{D}^\dagger\boldsymbol{\lambda}_{1:N}
=
-
\nabla_{\boldsymbol{X}}
\left\langle
\boldsymbol{\lambda}_{1:N},
\boldsymbol{F}_{1:N}
\right\rangle. D † λ 1 : N = − ∇ X ⟨ λ 1 : N , F 1 : N ⟩ . Instead, the LGL equation contains two additional endpoint terms:
e 1 w 1 ( μ − λ 1 ) \frac{\boldsymbol{e}_1}{w_1}
\left(
\boldsymbol{\mu}-\boldsymbol{\lambda}_1
\right) w 1 e 1 ( μ − λ 1 ) and
e N w N [ λ N − ∇ X Φ ( X N ) ] . \frac{\boldsymbol{e}_N}{w_N}
\left[
\boldsymbol{\lambda}_N-
\nabla_{\boldsymbol{X}}\Phi(\boldsymbol{X}_N)
\right]. w N e N [ λ N − ∇ X Φ ( X N ) ] . Therefore,
This differs from LG, where both endpoint conditions are decoupled, and LGR, where only one endpoint requires special treatment.
Why LGL Costate Estimates Can Be Less Accurate ¶ The LGL differentiation matrix is square and singular. Consequently:
there is no ordinary inverse of D \boldsymbol{D} D
the differential form has no uniquely equivalent integration matrix;
a null-space component can enter the discrete adjoint representation;
endpoint residuals remain embedded in the costate dynamics;
the resulting costate estimate may be less accurate than LG or LGR estimates.
The primal state and control solution may still be accurate. The limitation is specifically associated with recovering the continuous adjoint from the NLP multipliers.
Summary of Costate Mappings ¶ LG Collocation ¶ Let Λ 1 : N \boldsymbol{\Lambda}_{1:N} Λ 1 : N Λ N + 1 \boldsymbol{\Lambda}_{N+1} Λ N + 1
The costate at the Gauss points is
λ 1 : N = W − 1 Λ 1 : N + 1 Λ N + 1 . \boxed{
\boldsymbol{\lambda}_{1:N}
=
\boldsymbol{W}^{-1}\boldsymbol{\Lambda}_{1:N}
+
\boldsymbol{1}\boldsymbol{\Lambda}_{N+1}.
} λ 1 : N = W − 1 Λ 1 : N + 1 Λ N + 1 . The final costate is
λ N + 1 = Λ N + 1 . \boxed{
\boldsymbol{\lambda}_{N+1}
=
\boldsymbol{\Lambda}_{N+1}.
} λ N + 1 = Λ N + 1 . The initial costate is recovered from
λ 0 = λ N + 1 − d 0 T Λ 1 : N , \boxed{
\boldsymbol{\lambda}_0
=
\boldsymbol{\lambda}_{N+1}
-
\boldsymbol{d}_0^{\mathsf{T}}\boldsymbol{\Lambda}_{1:N},
} λ 0 = λ N + 1 − d 0 T Λ 1 : N , with the precise column notation chosen consistently with the LG differentiation matrix.
LGR Collocation ¶ The costate at the LGR points is
λ 1 : N = W − 1 Λ 1 : N . \boxed{
\boldsymbol{\lambda}_{1:N}
=
\boldsymbol{W}^{-1}\boldsymbol{\Lambda}_{1:N}.
} λ 1 : N = W − 1 Λ 1 : N . The terminal costate is
λ N + 1 = d N + 1 T Λ 1 : N . \boxed{
\boldsymbol{\lambda}_{N+1}
=
\boldsymbol{d}_{N+1}^{\mathsf{T}}\boldsymbol{\Lambda}_{1:N}.
} λ N + 1 = d N + 1 T Λ 1 : N . For a left-Radau grid,
so λ 1 \boldsymbol{\lambda}_1 λ 1
LGL Collocation ¶ The costate at all LGL points, including both endpoints, is
λ 1 : N = W − 1 Λ 1 : N . \boxed{
\boldsymbol{\lambda}_{1:N}
=
\boldsymbol{W}^{-1}\boldsymbol{\Lambda}_{1:N}.
} λ 1 : N = W − 1 Λ 1 : N . However, the endpoint transversality residuals remain coupled to the discrete adjoint dynamics.
Comparison of LG, LGR, and LGL ¶ Feature LG LGR LGL Included endpoints Neither One Both Differentiation matrix Rectangular Rectangular Square Matrix singularity No relevant square singularity No relevant square singularity Singular Equivalent integral form Yes Yes Not uniquely invertible Endpoint coupling Decoupled One endpoint partly coupled Both endpoints coupled Interior mapping W − 1 Λ + 1 Λ N + 1 \boldsymbol{W}^{-1}\boldsymbol{\Lambda}+\boldsymbol{1}\boldsymbol{\Lambda}_{N+1} W − 1 Λ + 1 Λ N + 1 W − 1 Λ \boldsymbol{W}^{-1}\boldsymbol{\Lambda} W − 1 Λ W − 1 Λ \boldsymbol{W}^{-1}\boldsymbol{\Lambda} W − 1 Λ Typical costate quality High High Can be less accurate
Comparison of Legendre pseudospectral costate estimation methods.
General Bolza Optimal-Control Problem ¶ Consider
min x ( ⋅ ) , u ( ⋅ ) , t 0 , t f J = Φ ( x ( t 0 ) , t 0 , x ( t f ) , t f ) + ∫ t 0 t f L ( x ( t ) , u ( t ) , t ) d t . \begin{aligned}
\min_{\boldsymbol{x}(\cdot),\boldsymbol{u}(\cdot),t_0,t_f}
J
={}&
\Phi
\left(
\boldsymbol{x}(t_0),t_0,
\boldsymbol{x}(t_f),t_f
\right)
\nonumber\\
&+
\int_{t_0}^{t_f}
L\left(
\boldsymbol{x}(t),\boldsymbol{u}(t),t
\right)\,\mathrm{d} t.
\end{aligned} x ( ⋅ ) , u ( ⋅ ) , t 0 , t f min J = Φ ( x ( t 0 ) , t 0 , x ( t f ) , t f ) + ∫ t 0 t f L ( x ( t ) , u ( t ) , t ) d t . The dynamics are
x ˙ ( t ) = f ( x ( t ) , u ( t ) , t ) . \dot{\boldsymbol{x}}(t)
=
\boldsymbol{f}\left(
\boldsymbol{x}(t),\boldsymbol{u}(t),t
\right). x ˙ ( t ) = f ( x ( t ) , u ( t ) , t ) . Endpoint constraints may be written as
ϕ ( x ( t 0 ) , t 0 , x ( t f ) , t f ) = 0 . \boldsymbol{\phi}
\left(
\boldsymbol{x}(t_0),t_0,
\boldsymbol{x}(t_f),t_f
\right)
=
\boldsymbol{0}. ϕ ( x ( t 0 ) , t 0 , x ( t f ) , t f ) = 0 . Path constraints may be written as
c min ≤ c ( x ( t ) , u ( t ) , t ) ≤ c max . \boldsymbol{c}_{\min}
\leq
\boldsymbol{c}\left(
\boldsymbol{x}(t),\boldsymbol{u}(t),t
\right)
\leq
\boldsymbol{c}_{\max}. c m i n ≤ c ( x ( t ) , u ( t ) , t ) ≤ c m a x . Mapping the Time Interval to [ − 1 , 1 ] [-1,1] [ − 1 , 1 ] ¶ Define
τ = 2 ( t − t 0 ) t f − t 0 − 1. \boxed{
\tau
=
\frac{2(t-t_0)}{t_f-t_0}-1.
} τ = t f − t 0 2 ( t − t 0 ) − 1. The inverse transformation is
t = t f − t 0 2 τ + t f + t 0 2 . \boxed{
t
=
\frac{t_f-t_0}{2}\tau
+
\frac{t_f+t_0}{2}.
} t = 2 t f − t 0 τ + 2 t f + t 0 . Differentiation gives
d τ d t = 2 t f − t 0 , \frac{\,\mathrm{d}\tau}{\,\mathrm{d} t}
=
\frac{2}{t_f-t_0}, d t d τ = t f − t 0 2 , and
d t d τ = t f − t 0 2 . \boxed{
\frac{\,\mathrm{d} t}{\,\mathrm{d}\tau}
=
\frac{t_f-t_0}{2}.
} d τ d t = 2 t f − t 0 . The endpoint states map as
x ( t 0 ) = x ( − 1 ) , x ( t f ) = x ( 1 ) . \boldsymbol{x}(t_0)=\boldsymbol{x}(-1),
\qquad
\boldsymbol{x}(t_f)=\boldsymbol{x}(1). x ( t 0 ) = x ( − 1 ) , x ( t f ) = x ( 1 ) . The cost becomes
J = Φ ( x ( − 1 ) , t 0 , x ( 1 ) , t f ) + t f − t 0 2 ∫ − 1 1 L ( x ( τ ) , u ( τ ) , t ( τ ; t 0 , t f ) ) d τ . \begin{aligned}
J
={}&
\Phi
\left(
\boldsymbol{x}(-1),t_0,
\boldsymbol{x}(1),t_f
\right)
\nonumber\\
&+
\frac{t_f-t_0}{2}
\int_{-1}^{1}
L
\left(
\boldsymbol{x}(\tau),\boldsymbol{u}(\tau),
t(\tau;t_0,t_f)
\right)
\,\mathrm{d}\tau.
\end{aligned} J = Φ ( x ( − 1 ) , t 0 , x ( 1 ) , t f ) + 2 t f − t 0 ∫ − 1 1 L ( x ( τ ) , u ( τ ) , t ( τ ; t 0 , t f ) ) d τ . The transformed dynamics are
d x d τ = t f − t 0 2 f ( x ( τ ) , u ( τ ) , t ( τ ; t 0 , t f ) ) . \boxed{
\frac{\,\mathrm{d}\boldsymbol{x}}{\,\mathrm{d}\tau}
=
\frac{t_f-t_0}{2}
\boldsymbol{f}
\left(
\boldsymbol{x}(\tau),\boldsymbol{u}(\tau),
t(\tau;t_0,t_f)
\right).
} d τ d x = 2 t f − t 0 f ( x ( τ ) , u ( τ ) , t ( τ ; t 0 , t f ) ) . The endpoint constraints become
ϕ ( x ( − 1 ) , t 0 , x ( 1 ) , t f ) = 0 . \boldsymbol{\phi}
\left(
\boldsymbol{x}(-1),t_0,
\boldsymbol{x}(1),t_f
\right)
=
\boldsymbol{0}. ϕ ( x ( − 1 ) , t 0 , x ( 1 ) , t f ) = 0 . The path constraints become
c min ≤ c ( x ( τ ) , u ( τ ) , t ( τ ; t 0 , t f ) ) ≤ c max . \boldsymbol{c}_{\min}
\leq
\boldsymbol{c}
\left(
\boldsymbol{x}(\tau),\boldsymbol{u}(\tau),
t(\tau;t_0,t_f)
\right)
\leq
\boldsymbol{c}_{\max}. c m i n ≤ c ( x ( τ ) , u ( τ ) , t ( τ ; t 0 , t f ) ) ≤ c m a x . LGR Discretization of the General Bolza Problem ¶ For a left-LGR grid, the initial point is τ 1 = − 1 \tau_1=-1 τ 1 = − 1 τ N + 1 = 1 \tau_{N+1}=1 τ N + 1 = 1
Discrete Objective ¶ The objective is approximated by
J N = Φ ( X 1 , t 0 , X N + 1 , t f ) + t f − t 0 2 ∑ i = 1 N w i L ( X i , U i , t i ) , \begin{aligned}
J_N
={}&
\Phi
\left(
\boldsymbol{X}_1,t_0,
\boldsymbol{X}_{N+1},t_f
\right)
\nonumber\\
&+
\frac{t_f-t_0}{2}
\sum_{i=1}^{N}
w_i
L
\left(
\boldsymbol{X}_i,\boldsymbol{U}_i,t_i
\right),
\end{aligned} J N = Φ ( X 1 , t 0 , X N + 1 , t f ) + 2 t f − t 0 i = 1 ∑ N w i L ( X i , U i , t i ) , where
t i = t f − t 0 2 τ i + t f + t 0 2 . t_i
=
\frac{t_f-t_0}{2}\tau_i
+
\frac{t_f+t_0}{2}. t i = 2 t f − t 0 τ i + 2 t f + t 0 . Discrete Dynamics ¶ The collocation equations are
D X 1 : N + 1 = t f − t 0 2 F 1 : N . \boxed{
\boldsymbol{D}\boldsymbol{X}_{1:N+1}
=
\frac{t_f-t_0}{2}
\boldsymbol{F}_{1:N}.
} D X 1 : N + 1 = 2 t f − t 0 F 1 : N . Equivalently,
D X 1 : N + 1 − t f − t 0 2 F 1 : N = 0 . \boldsymbol{D}\boldsymbol{X}_{1:N+1}
-
\frac{t_f-t_0}{2}
\boldsymbol{F}_{1:N}
=
\boldsymbol{0}. D X 1 : N + 1 − 2 t f − t 0 F 1 : N = 0 . Endpoint Constraints ¶ The endpoint constraints become
ϕ ( X 1 , t 0 , X N + 1 , t f ) = 0 . \boxed{
\boldsymbol{\phi}
\left(
\boldsymbol{X}_1,t_0,
\boldsymbol{X}_{N+1},t_f
\right)
=
\boldsymbol{0}.
} ϕ ( X 1 , t 0 , X N + 1 , t f ) = 0 . Path Constraints ¶ At each collocation point,
c min ≤ c ( X i , U i , t i ) ≤ c max , i = 1 , … , N . \boxed{
\boldsymbol{c}_{\min}
\leq
\boldsymbol{c}
\left(
\boldsymbol{X}_i,\boldsymbol{U}_i,t_i
\right)
\leq
\boldsymbol{c}_{\max},
\qquad
i=1,\ldots,N.
} c m i n ≤ c ( X i , U i , t i ) ≤ c m a x , i = 1 , … , N . Interpretation of the Time-Scaling Factor ¶ The only structural change in the dynamics and integral cost is the factor
t f − t 0 2 . \boxed{
\frac{t_f-t_0}{2}.
} 2 t f − t 0 . When
t 0 = − 1 , t f = 1 , t_0=-1,
\qquad
t_f=1, t 0 = − 1 , t f = 1 , this factor equals one. Thus, the simplified problem used in the derivations is a special case of the general transcription.
The initial and final times may themselves be NLP variables. In that case, they enter:
the endpoint cost;
the quadrature scaling;
the transformed dynamics;
the physical time values t i t_i t i
the endpoint constraints;
the path constraints.
Extension to LG and LGL ¶ The same time transformation applies to LG and LGL schemes.
The LG transcription additionally requires a quadrature equation for the terminal state because neither endpoint is a collocation point.
LGR ¶ The LGR transcription requires one noncollocated endpoint state and has a rectangular differentiation matrix.
LGL ¶ The LGL transcription includes both endpoint states directly and has a square singular differentiation matrix.
Differential and Integral Interpretations ¶ For LG and LGR, a suitable square portion of the rectangular differentiation matrix is invertible. Its inverse acts as an integration matrix. Consequently, their differential collocation forms have equivalent integral forms.
For LGL, the differentiation matrix is singular. Therefore, no unique ordinary inverse exists, and the differential form does not possess the same uniquely equivalent integration representation.
This distinction explains much of the difference in costate recovery quality.
Implementation Considerations ¶ A practical implementation should carefully manage:
the ordering of state, control, time, and parameter variables;
the ordering of defect, endpoint, and path constraints;
the sign convention used for equality constraints;
the association between solver multipliers and constraint blocks;
state, control, defect, and objective scaling;
exact or algorithmic derivatives;
mesh refinement;
endpoint consistency;
primal and dual residuals.
Summary ¶ LGL collocation includes both endpoints.
Its state differentiation matrix is square and singular.
The LGL costate estimate is obtained by
λ 1 : N = W − 1 Λ 1 : N . \boldsymbol{\lambda}_{1:N}
=
\boldsymbol{W}^{-1}\boldsymbol{\Lambda}_{1:N}. λ 1 : N = W − 1 Λ 1 : N . The transformed adjoint equation contains both initial and terminal boundary residuals.
These endpoint terms couple the transversality conditions to the costate dynamics.
LG and LGR generally provide cleaner costate mappings because their rectangular differentiation matrices admit equivalent integration forms.
A general Bolza problem on [ t 0 , t f ] [t_0,t_f] [ t 0 , t f ] [ − 1 , 1 ] [-1,1] [ − 1 , 1 ]
The transformed cost and dynamics acquire the factor
( t f − t 0 ) / 2. (t_f-t_0)/2. ( t f − t 0 ) /2. Endpoint and path constraints are imposed directly in the resulting NLP.
Connection. With the mathematical transcription complete, attention shifts to the data layout and function interfaces required by a general NLP solver.