Higher Order Numerical Approximation of Switching Systems

Higher order numerical schemes for affine nonlinear control systems developed elsewhere by the authors are adapted to switching systems with prescribed switching times. In addition the calculation of the required multiple switching control integrals is discussed. These schemes are particularly useful in situations involving very rapid switching.


Introduction
Switching systems have been of considerable interest for many years in the computing and scientific computing communities.They are now being intensively investigated by control and systems engineers, see e.g.[4,8,9,10] and the papers cited therein, who are interested primarily in the asymptotic stability or stabilizability of a common equilibrium solution, which is usually taken to be the zero solution.Switching systems have have also been considered from the perspective of nonautonomous dynamical systems theory, e.g.[5,6].
In principle known numerical schemes for ordinary differential equations such as Runge-Kutta schemes can be applied to switching systems, changing the vector field after each switch has occured.However, in order to maintain the usual consistency order of these schemes, the integration time steps need to be adjusted to the switching times in such a way that switching always occurs at the end of an integration interval.This is impractical in the case of fast switching, because in this case an adjustment of the scheme's integration step size to the switching times would lead to very small time steps causing an undesirably high computational load.A class of higher order numerical schemes which are suitable for systems with a very fast rate of switching is proposed in this paper.
Consider a finite collection of ordinary differential equations in R d and assume that there is a sequence of indices i j ∈ {1, • • • , m} for j ∈ Z + and an increasing sequence of times 0 We call (2) a switching system.Introducing the corresponding switching signal, i.e the piecewise constant function σ : [τ 0 , ∞) → {1, . . ., m} which is continuous from the right and defined by σ(t) := i j if τ j ≤ t < τ j+1 , j ∈ Z + , we can rewrite the switching system (2) as the single nonautonomous ordinary differential equation dx dt = f σ(t) (t, x(t)).
It is easy to see that (3) is equivalent to the affine control system with switching controls u j : R + → {0, 1} defined by Since only one of the systems is active at any one time we always have m j=1 u j (t) ≡ 1.
In this paper we adapt the higher order Taylor schemes for control affine systems proposed by the authors in [2] to switching systems (3) reformulated as (4).The switching times appear in these schemes only within the multiple switching control integrals that they contain.We show how these Taylor schemes can be adapted to give Runge-Kutta like derivative free schemes and simplified when the vector fields commute.Examples are given for linear systems and some numerical results are also presented.Finally, the calculation of the multiple switching control integrals which appear in the schemes is dicussed, in particular for periodic and rapidly switching.

Numerical schemes for switching systems
Switching systems can be interpreted as control affine systems, so the higher order numerical schemes for such systems which were developed in [2] can also be used for a switching system (4).Essentially, the idea is to take Taylor approximations of the solutions of the differential equation ( 4) interpreted as the integral equation The situation here is simpler than in [2] since the general notation introduced there is not required here and because an additional uncontrolled term f 0 (t, x) there is not present here.Let t 0 < t 1 < . . .< t n < . . .< t N = T be a partition of the time interval [t 0 , T ] with stepsizes ∆ n = t n+1 − t n and maximal step size ∆ := max n ∆ n .Consider the multiple switching control integrals and the coefficient functions where the partial differential operators are defined by (The operators are applied initially to the identity function f (x) ≡ x in ( 7)).
2.1 Taylor schemes of orders 1, 2, and 3 for switching systems The Euler scheme is the simplest nontrivial Taylor scheme, having convergence order p = 1.For the switching system (4) it is given by where the switching control integrals represent the duration of time in discretization subinterval [t n , t n+1 ] that the jth system is active.Thus the summation term in (10) is a weighted averaging of the switching systems evaluated at the start of the discretization interval.An analgous situation occurs for the higher order terms in the higher order Taylor schemes applied to the switching system (4).The Taylor scheme of order 2 for the switching system (4) is and the Taylor scheme of order 3 for the switching system (4) The coefficients of these Taylor schemes are often tedious to derive, especially in high dimensions.However maple routines are available to undertake this task, cf.Section 5.These routines are based on a modification of the corresponding routines for stochastic differential equations from [1].
Finally, we note that implicit versions of the above schemes can be derived in a similar way and have the same order of convergence.

Derivative-free schemes
It is often computationally advantageous to have numerical schemes which avoid the use of derivatives of the coefficients in much the same way that Runge-Kutta schemes do in the traditional setting.Essentially one seeks to approximate the derivatives appearing in a Taylor scheme by difference quotients of an appropiate order in order to retain the order of the Taylor scheme.For example in the second order Taylor scheme (11) for scalar switching systems one can use in the (i, j) term of the Taylor scheme.Since this is mutipled by a multiple integral of magnitude O(∆ 2 n ), the error made is of the same order O(∆ 3 n ) as the local discretization error of the original Taylor scheme.This yields the second order derivative-free scheme for switching systems given by Without switching this reduces to the second order Runge-Kutta scheme for an ordinary differential equation which is usually known as the Heun scheme.This principle can be extended to obtain higher order derivative-free schemes.See [7] for analogous higher order derivative-free schemes for the stochastic case.

Commutative switching systems
The Taylor schemes (9) simplify considerable when the coefficients f 1 , . .., f m of the switching system (4) satisfy special properties.For example, the switching system is said to be commutative when the coefficients satisfy (a maple routine to check for the validity of this condition can be found in [1]).Then, by the generalized integration-by-parts identities the sum of terms which involves more easily computed multiple switching control integrals of lower multiplicity.In this case the Taylor scheme of order 2 for commutative switching system simplifies to (4) 3 Linear switching systems It is illustrative to see how our proposed schemes simplify in case of linear dynamics.For this porpose, consider the switching system (4) with linear autonomous vector fields f j (t, x) ≡ A j x, j = 1, . . ., m.
x, so the Euler scheme becomes the Taylor scheme of order 2 and the Taylor scheme of order 3 The derivative free scheme (13) here reduces to the Taylor scheme of order 2 for linear systems.(18).If the matrices A 1 , . .., A m pairwise commute, then the Taylor scheme of order 2 for commutative switching systems (16) simplifies to (20) It is interesting to compare the Euler scheme (17) for two systems, where the first system is active from t n until t n + ∆ n,1 and the second from t n + ∆ n,1 till t n+1 , where ∆ n,1 + ∆ n,2 = ∆ n , i.e.
with the usual Euler scheme applied to the first system on the interval [t n , t n + ∆ n,1 ] i.e.
and then to the second system on the interval [t n + ∆ n,1 , t n+1 ], i.e.
which can be concatenated to give Compared to the switched Euler scheme (17) the two step scheme (22) contains an additional second order term and its local error is of magnitude ∆ instead of ∆ 2 n on this time discretization subinterval.Note, however, that for a larger number of switchings the number of terms in (22) will increase drastically while in (17) only the computation of the integral term I (j) becomes more expensive.

Multiple switching control integrals
In this section we show that the multiple control integrals which appear in the Taylor expansion and in the numerical schemes can be computed explicitly without using numerical integration.We first consider some simple illustrative cases and then discuss periodic switching in order to motivate our general result which is presented in the final subsection.We remark that the required multiple switching integrals could be calculated either off-line in advance or in parallel to the scheme.

Some simple cases
As already mentioned, the simple first order switching control integral ∆ n,j = I (j),tn,t n+1 = t n+1 tn u j (s) ds, j = 1, . . ., m, represent the total duration of time in the discretization subinterval [t n , t n+1 ] that the jth system is active, possibly switching off and on in this interval.
If there is no switch in the discretization subinterval [t n , t n+1 ] and only the j * th system is active in this interval, i.e. u j * (t) ≡ 1 for t ∈ [t n , t n+1 ], then the only the integrals involving the index j * are non zero, in fact with where k is the length of the mulit-index (j * , . . ., j * ).On the other hand, if the j * th system is never active in this interval, i.e. u j * (t) ≡ 0 for t ∈ [t n , t n+1 ], then not unexpectedly all multiple integrals involving the index j * are zero.
For illustration, suppose that system 1 is active in [t n , t n + ∆ n,1 ) and system 2 is active in [t n + ∆ n,1 , t n+1 ].Then Similarly and one can verify directly the identity

Periodic switching signals
The computational cost for the evaluation of the multiple control integrals is particularly high in the case of fast switching, i.e., when the number of switches per integration interval [t n , t n+1 ] is large.If, however, the switching signal is periodic, i.e., if there exists T per > 0 with σ(t + T per ) = σ(t) for a small period T per , then computational burden can be considerably lowered.
In the case of rapid switching, we can use a constant integration time step ∆ n ≡ ∆ which is an integer multiple of T per , i.e., ∆ = kT per for some k ∈ N. Then the multiple control integrals satisfy I (j 1 ,...,j l ),tn,t n+1 = I (j 1 ,...,j l ),t 0 ,t 1 for all n ∈ N 0 .Hence, it is sufficient to evaluate these integrals only in the first step and use the stored values in each subsequent time step.Similarly for slow periodic switching, we can also use a constant integration time step ∆ n ≡ ∆ which is an integer fraction of T per , i.e., ∆ = T per /k for some k ∈ IN .Then the multiple control integrals satisfy I (j 1 ,...,j l ),t n+k ,t n+k+1 = I (j 1 ,...,j l ),t 0 ,t 1 for all n ∈ N 0 .

A general result
The following theorem provides a general formula for switching integrals when several switches occur within a discretization subinterval.
Here the convention 0! = 1 is used.
The proof will be given in Section 7. Note that I (j 1 ,...,j l ) [1] tn,t n+1 = 0 in the vacuous case that any one of the indexed switching controls is not active at all during the integration time interval, because in this case the condition in the third line of the sum will never be satisfied, hence the sum is empty.

Numerical examples
In this section we illustrate our method by solving switching systems using the Taylor scheme.The computations have been performed with several maple and matlab routines which are available from www.uni-bayreuth.de/departments/math/∼lgruene/papers/switching.html and which are used as follows: first, the maple routine staylorcoeff.mw is used in order to compute the coefficients of the pth-order Taylor scheme.These coefficients are written to a matlab M-file staylorcoeff_name_p.m with "name" indicating the name of the vector field and "p" denoting the order of the scheme.This file serves as an input for the actual implementation of the scheme in the matlab M-file staylor.m.This file needs the M-file mint_switch.mwhich computes the multiple control integrals according to Theorem 2. For more information see the example file teststaylor.mwhich is also available from the above web address and which contains all the examples in this section.
Our first example is a two dimensional linear switching system with m = 2 vector fields f j (t, x) ≡ A j x, j = 1, 2 and matrices given by The vector fields describe two linear oscillators without damping and with different parameters.Each individual system is (neutrally) Lyapunov stable.Using our numerical algorithm we investigate the behavior of the switching systems for two different switching sequences.The first sequence is given by We compute the solution for initial value x 0 = (1, 0) T on the interval [0, 35] using the Taylor scheme of order p = 5 with constant time step ∆ = 0.05.The results, shown in Figure (1), confirm the instability of the switching system.
Using the time step ∆ = 0.04 we are able to make use of the periodicity of the switching signal.We computed the solution with initial value x 0 = (0, 1) T on the interval [0, 20] using again the 5th order Taylor scheme.The numerical solution, shown in Figure (2), indicates that the system remains stable for this faster switching sequence.Our final example is a nonlinear 2d system with m = 2 vector fields given by which was taken from [9, Example 2.4].There it was shown that despite the fact that each individual system is globally uniformly asymptotically stable, the switching system can be unstable.In order to confirm this property numerically, we arbitrarily picked the switching sequence given by Once again using the 5th order Taylor scheme, here with time step ∆ = 0.01 we simulated the solution with initial value x 0 = (2.1, 1) T on the interval [0, 2.5].
The result, shown in Figure 3, indicates the instability of the switching system (28), (29).We need to show that the scheme (9) has global discretization error of order p when the coefficients f 1 , . .., f m of the switching system (4) are p times continuously differentiable.This follows by standard arguments [3] from the consistency order of the scheme which in turn follows from the form of the error term in the corresponding Taylor expansion [2].In this paper Taylor schemes are systematically developed for general control affine systems of the form with u 0 (t) ≡ 1.
Here we show that these schemes reduce to (9) in the case of switching systems, which we can do without introducing the abstract notation used in [2].Suppose that the system here, i.e., the x variable, is d + 1 dimensional, which for convenience we write as x = (x 0 , x 1 , . . ., x d ), i.e indexing the components from 0 instead of from1.The general Taylor expansion and schemes in [2] use (with new notation here) the differential operators If we apply this to time-independent functions of a solution of an autonomous system, i.e with f j (t, x) ≡ f j (x), such that f 0 (x) ≡ 0, then the L 0 operator annihilates whatever it is applied to, so has no contributes nothing to the scheme.At the same time, the L j operators can be written With f j,0 (t, x) ≡ 1 and x 0 ≡ t, this reduces to This is exactly the situation which arises when we autonomize a nonautonomous system by expanding the dimension of the system with a new state variable x 0 ≡ t and new state equation with f j,0 (x, t) ≡ 1 for j = 1,. ..,m.
Thus the pth order Taylor scheme for general affine systems reduces to (9) for the switching systems (4).
Hence the integral in (32) evaluates to which is again of the asserted form.