Theorem 5.25. Consider a hybrid system H that satisfies Hypothesis 5.6 with hybrid time trajectory Tτ, corresponding hybrid mode trajectory Tm, and a initial value x0 ∈ Rn. Let E(Tτ) be the set of event times and let the assumptions of Theorem 5.25 be fulfilled.
Further, let Assumption 5.22 hold. Then there exists a continuous solutionx of the hybrid system H in the sense of Definition 5.15 (with El as in (5.18)) if and only if
1. the initial value x0 is consistent for the DAE in the initial mode l1 ∈Tm, 2. the transition functions Tllii−1 are the identity mappings, i.e.,
Tlli−1i (x(τi′−1),x(τ˙ i′−1)) = [x(τi′−1),x(τ˙ i′−1)] = [x(τi),x(τ˙ i)],
and for everyτi ∈ E(Tτ)the statesx(τi)are consistent with the DAE in modeli ∈Tm. The continuous solution x is unique if and only if in additionulµ = 0 for all modesl ∈Tm. Proof. The proof is analogous to the proof of Theorem 5.18.
Ll
µl. To determine a consistent initial value or to check it for consistency we must therefore solve the underdetermined system
Fµll(τi, x0, y0) = 0 (5.19) for (x0, y0). We use the Gauss-Newton method [33, 113] started with a sufficiently good initial guess (˜x0,y˜0) to solve this underdetermined systems of nonlinear equations in a least squares sense. The Gauss-Newton methodfor a nonlinear system of the form
G(z) = 0,
where G:Rn →Rm is a smooth function, generates a sequence{zk} of approximations of the form
zk+1 =zk−G+;z(zk)G(zk), (5.20) starting with an initial guess z0 ∈ Rn. Here, G+;z(z) is the Moore-Penrose pseudo-inverse of the JacobianG;z(z) = [∂G∂zi
j(z)]i,j (see Definition 2.17).
Theorem 5.26. Let G : D ⊆ Rn → Rm, m < n, with D open, convex denote a con-tinuously differentiable mapping, and assume a full row rank of the Jacobian. Consider the Gauss-Newton method (5.20) and assume that a starting point z0 ∈D, and constants α, ω ≥ 0 exist such that kG+;z(z0)G(z0)k ≤ α, and kG+;z(u)(G;z(v)− G;z(u))(v −u)k ≤ ωkv−uk2 for all u, v ∈D, v−u∈range (G+;z(u)). Moreover, let
h:=αω <2, S(z0, r)⊂D with r := 2α/(2−h).
Then:
1. The sequence {zk} of Gauss-Newton iterates is well-defined, remains in S(z0, r) and converges to some z⋆ ∈S(z0, r) with G+;z(z⋆)G(z⋆) = 0.
2. Quadratic convergence can be estimated according to kzk+1−zkk ≤ 1
2ωkzk−zk−1k2. Proof. See [33, Theorem 4.19].
Thus, if the JacobianG;z(z) has full row rank and fulfills a Lipschitz condition in an open convex set, and if G+;z(zk) is bounded in this set, then we have quadratic convergence of the iterates {zk} to a least squares solutionz⋆ of G(z) = 0 satisfying
G+;z(z⋆)G(z⋆) = 0.
Since the Jacobian of (5.19) at a solution (x0, y0) has full row rank at a solution and thus in a whole neighborhood by Hypothesis 2.37, we have local quadratic convergence of
the Gauss-Newton method, provided that the starting point z0 is sufficiently close to the solution.
During the numerical integration of a hybrid system, events cause switching from mode l to mode k at a switch point τi. The transition function Tlk as defined in (5.6) maps the state at the switch point in mode l to the state at the switch point in the new modek via
Tlk(xl(τi),x˙l(τi)) = [x⋆,x˙⋆].
But, the transfered state x⋆ is not necessarily consistent with the DAE in mode k at time τi, and to continue the integration a consistent initial value xk(τi) at τi has to be computed. In order to find a reasonable continuation of the solution of the hybrid system, we try to find a consistent initial value xk(τi) at τi from among all consistent values in the constraint manifold Lk, on the basis of the given but inconsistent initial state x⋆, in such a way that the solution xk extends the past solution xl in a physically reasonable way. Since algebraic variables need to be chosen consistently with the DAE in the current mode they have to be computed as the solution of a nonlinear system describing the algebraic constraints. On the other hand, initial values for differential variables and possibly undetermined variables can be chosen freely, such that these components of the initial value vector should be kept fixed during the computation of consistent initial values in order to find a continuation of the hybrid system solution that is as smooth as possible. Thus, even if the transition function Tlk provides continuity of the state variables over a switch point the consistent reinitialization can cause discontinuities in the solution. If possible, these discontinuities should only occur in the algebraic variables, which have to be consistent, while the differential variables and undetermined variables should proceed continuously over the switch point.
By a slight modification of the above approach it is possible to prescribe initial values for the differential variables and controls, i.e., fix certain components of an initial guess ˜x0
during the computation of consistent initial values, and only compute consistent values for the algebraic variables. This requires the classification of a component of ˜x0 to be a differential variable or a control, such that eliminating the associated columns from the Jacobian of the nonlinear system (5.19) does not lead to a rank deficiency, since we must guarantee that the remaining columns of the Jacobian still have full row rank to ensure quadratic convergence of the Gauss-Newton method. Due to Hypothesis 2.37, there exist continuous matrix functions
Z2l ∈C(I,R(µl+1)nl,alµ), T2l ∈C(I,Rnl,nl−alµ), Z1l ∈C(I,Rnl,dlµ),
with the properties described in Hypothesis 2.37 and pointwise orthonormal columns. Let the matrix functions
Z2l′ ∈C(I,R(µl+1)nl,(µl+1)nl−alµ), T2l′ ∈C(I,Rnl,alµ), Z1l′ ∈C(I,Rnl,nl−dlµ), be chosen such that
Z2l′ Z2l ,
T2l′ T2l ,
Z1l′ Z1l
are pointwise orthogonal, i.e., in particular nonsingular. Then, Hypothesis 2.37 yields Fˆ2;xl T2l = 0, rankT2l =nl−alµ, rank ˆF1; ˙lxT2l =dlµ,
and multiplication with the nonsingular matrix
T2l′ T2
yields the separation rank
Fˆ1; ˙lx Fˆ2;xl
= rank
Fˆ1; ˙lxT2l′ Fˆ1; ˙lxT2l Fˆ2;xl T2l′ 0
= rank ˆF1; ˙lxT2l+ rank ˆF2;xl T2l′ =dlµ+alµ.
Now, let ˜T2l be a fixed approximation with orthonormal columns to T2l that spans the nullspace of ˆF2;xl at the desired solution. Then, we can solve
Fµll(t0,T˜2lT˜2lTx˜0+ (I−T˜2lT˜2lT)x,x, . . . , x˙ (µl+1)) = 0 (5.21) for (x,x, . . . x˙ (µl+1)), with initial guess ˜x0, where ˜T2lT˜2lT is an orthogonal projection of rank dlµ onto kernel ˆF2;xl , while I −T˜2lT˜2lT is a projection onto cokernel ˆF2;xl . The dlµ differential components of the initial guess ˜x0 are kept fixed during the Gauss-Newton iterations as the corresponding columns of the Jacobian are set to zero. Note that this approach will lead to a rank drop in the Jacobian if any of the algebraic variables are fixed. A drawback of this approach to solve the nonlinear systems (5.19) or (5.21) with the Gauss-Newton method is the limited region of convergence. This means that the Gauss-Newton method may not converge if the initial guess is not sufficiently close to the solution. Therefore, after mode switching, the starting value (τi, x⋆) for the Gauss-Newton iteration given by the transition function Tlk should be sufficiently close to a solution in the new mode to guarantee convergence.
After a mode change from mode l to mode k, differential variables in the predecessor mode l may change to algebraic or undetermined variables in the successor mode k or vice versa. Assuming that nl = nk (otherwise undetermined variables can be inserted into the system to meet this requirement, see Section 5.3) the different possibilities are summarized in Table 5.1. Whenever algebraic variables or undetermined parts change into differential variables, no problems with consistency occur as the initial conditions fits into the differential equation (cases 9,10,11,12 in Table 5.1). Thus, ifakµ ≤alµ then it is possible to obtain a continuous solution provided that the constraint manifold has not changed. On the other hand, if differential or undetermined variables change into algebraic variables, then inconsistency can occur and reinitialization results in discontinuities in the solution (cases 2,4,6,8 in Table 5.1). Ifukµ>0, then the solution is not unique and the DAE can only be solved in a least squares sense. If variables change into undetermined variables (cases 3,5,6,7,12 in Table 5.1) and the least squares solution is obtained in such a way thatkxk2 is minimized, then the continuity condition and the minimum norm condition can contradict each other. Therefore, the minimization problem for the least squares solution should be chosen as in (5.16). In addition, the index of the differential-algebraic system might have changed form µl to µk. If µl ≥ µk, then no problems occur, but if the index increases after a mode change, then higher smoothness requirements are needed to guarantee the
Differential part Algebraic part Undetermined part Changes in char. val.
1 dlµ =dkµ alµ=akµ ulµ=ukµ — 2 dlµ =dkµ alµ< akµ ulµ> ukµ uya 3 dlµ =dkµ alµ> akµ ulµ< ukµ ayu 4 dlµ > dkµ alµ< akµ ulµ=ukµ dya 5 dlµ > dkµ alµ=akµ ulµ< ukµ dyu 6 dlµ > dkµ alµ< akµ ulµ< ukµ dya+u 7 dlµ > dkµ alµ> akµ ulµ< ukµ d+ayu 8 dlµ > dkµ alµ< akµ ulµ> ukµ d+aya 9 dlµ < dkµ alµ=akµ ulµ> ukµ uyd 10 dlµ < dkµ alµ> akµ ulµ=ukµ ayd 11 dlµ < dkµ alµ> akµ ulµ> ukµ a+uyd 12 dlµ < dkµ alµ> akµ ulµ< ukµ ayd+u 13 dlµ < dkµ alµ< akµ ulµ> ukµ uyd+a
Table 5.1: Changes in the characteristic quantities after a mode change from mode lto mode k
existence of a solution and more effort is needed to obtain the reduced system in the new mode which might alter the convergence of numerical methods.
Remark 5.27. For the computation of consistent initial value for a linear DAE (2.5) (see also [82, p. 308ff]) we can consider the reduced system of the form (2.24), where the algebraic equations are displayed directly. The condition for a givenx˜0 att0 to be consistent is given by
Aˆ2(t0)˜x0+ ˆb2(t0) = 0. (5.22) In the case that the given x˜0 is not consistent, we can use (5.22) to determine a related consistentx0. Settingx˜0 =x0+δ we determine the correctionδby solving the minimization problem
kδk2 =min!
subject to the constraint
kAˆ2(t0)δ−ˆb2(t0)−Aˆ2(t0)˜x0k2 =min!.
The solution of this least squares problem is given by
δ= ˆA+2(t0)( ˆA2(t0)˜x0+ ˆb2(t0)),
where Aˆ+2(t0) is the Moore-Penrose pseudo-inverse of Aˆ2(t0). Since Aˆ2(t0) has full row rank aµ, due to Theorem 2.41, it follows that Aˆ2(t0) ˆA+2(t0) =Ia, and therefore
Aˆ2(t0)x0+ ˆb2(t0) = ˆA2(t0)(˜x0−δ) + ˆb2(t0) = ˆA2(t0)˜x0−( ˆA2(t0)˜x0+ ˆb2(t0)) + ˆb2(t0) = 0.
Also in this case, we can prescribe initial values for the differential variables, whereas initial values for the algebraic variables are not known. This requires a separation of the unknown x into differential, algebraic and unknown parts. For a system in the form (2.24) we can compute an orthogonal matrix U = [U1, U2] of size (n, n) such that
Eˆ1(t0)
U1 U2
=
E11 0 ,
where E11 has size (dµ, dµ) and is nonsingular. Then, we determine an orthogonal matrix V = [V1, V2] of size (n−dµ, n−dµ) such that
Aˆ2(t0)U2V =
A22 0 ,
where A22 is of size (aµ, aµ) and nonsingular. This allows a reinterpretation of variables as differential, algebraic or undetermined variables using the basis transformation
x=Q
˜ x1
˜ x2
˜ x3
, Q=
U1 U2V1 U2V2 , with orthogonal matrix Q and corresponding DAE
E11 0 0
0 0 0
0 0 0
˙˜
x1
˙˜
x2
˙˜
x3
=
A11 A12 A13
A21 A22 0
0 0 0
˜ x1
˜ x2
˜ x3
+
ˆb1
ˆb2 ˆb3
.
From the second block row we get a partitioning of the consistency condition (5.22) into 0 = [A21(t0)A22(t0)]
x˜1
˜ x2
+ ˆb2(t0).
Now, let an estimate x˜0 = (˜x1,0,x˜2,0) for a consistent initial value be given. Keeping x˜1,0
fixed, we can determine a correction δ2 for the estimate x˜2,0 = x2,0 +δ2 by solving the minimization problem
kδ2k2 =min!
subject to the constraint
kAˆ22(t0)δ2−Aˆ2(t0)˜x0−ˆb2(t0)k2 =min!,
i.e.,δ2 = ˆA+22(t0)( ˆA2(t0)˜x0+ˆb2(t0)). The corrected consistent initial condition is then given by x1,0 = ˜x1,0 and x2,0 = ˜x2,0−δ2, and thus x0 =Q
x˜1,0
˜
x2,0−δ2
.