Discrete Dynamics in Nature and Society

Volume 2017 (2017), Article ID 6896153, 11 pages

https://doi.org/10.1155/2017/6896153

## Numerical Feedback Stabilization with Applications to Networks

Department of Mathematics, University of Mannheim, 68131 Mannheim, Germany

Correspondence should be addressed to Simone Göttlich

Received 22 February 2017; Accepted 13 April 2017; Published 22 May 2017

Academic Editor: Gisèle R. Goldstein

Copyright © 2017 Simone Göttlich and Peter Schillen. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

#### Abstract

The focus is on the numerical consideration of feedback boundary control problems for linear systems of conservation laws including source terms. We explain under which conditions the numerical discretization can be used to design feedback boundary values for network applications such as electric transmission lines or traffic flow systems. Several numerical examples illustrate the properties of the results for different types of networks.

#### 1. Introduction

During the last few years, a huge amount of literature has emerged that deals with theoretical stability results for boundary control of conservation laws; see, for example, [1–15]. Further work also includes hyperbolic systems with source terms (so-called balance laws) for special applications such as gas dynamics [16–18] or water flow in open canals [19–23]. To the best of our knowledge only a few results are available, where these kinds of control problems are analyzed from a numerical point of view; see [24, 25].

Many network applications are represented by one-dimensional hyperbolic systems of balance laws [4, 26–30]. For the analytical case it can be proven that feedback boundary values, designed under certain conditions, yield an exponential decay of a continuous Lyapunov function [3, 9, 31] also in the context of networks.

In this paper, we are mainly concerned with the numerical analysis of the stability of the steady-state of such networked systems. We explain how suitable boundary values must be designed such that a linear system of balance laws converges exponentially to a steady-state. Inspired by [24, 25], we restrict ourselves to the numerical approximation of such systems and study the asymptotic behavior of that solution or to be more precise the conditions under which an exponential decay of the discrete solution to hyperbolic balance law can be attained. One way to influence a system of balance laws is through the boundary values. The interpretation of the boundary values in the case of hyperbolic systems is that they describe an outflow or an inflow of a quantity. One can measure the outflow, modify it, and let it flow back again as inflow. This kind of control is then called feedback boundary control.

Our considerations are done first for arbitrary linear systems of balance laws and are then extended to networks which are in fact coupled systems of boundary problems. For simulation purposes, we present two concrete examples. One example is the so-called telegrapher’s equation, an equation described by a hyperbolic balance law and normally used to model transmission lines [29]. The other example is a model for road traffic, the so-called LWR-model [32]. It is a scalar model describing traffic flow based on physical conservation laws. The question we address here is how a feedback control can be used that small errors in some measured data are dampened.

The outline of the article is as follows: in Section 2, we answer the question under which conditions an exponential decay of the numerical solution to the hyperbolic equations can be achieved. The network applications are given in Section 3. We design feedback boundary values for different networks and analyze the stability conditions.

#### 2. Problem Description and Preliminary Work

Hyperbolic systems of balance laws are often used to describe physical applications such as gas, water, traffic, or power flow. Boundary conditions can be designed to ensure that the analytical solution for these systems converges towards a desired steady-state. To prove the convergence of the approximate solution to a desired steady-state, a technique based on Lyapunov functions [3, 4, 8, 28] is used. The required boundary conditions are called feedback boundary conditions and the method is the boundary feedback stabilization of hyperbolic systems. This is a common approach for controlling systems of conservation laws.

In two recent works [24, 25], the theoretical results mentioned above are transferred to the numerical discretization of systems of balance laws. It can be shown that the approximate solution also converges exponentially to a steady-state. This is done using a suitable discretized Lyapunov function. Furthermore, this approach allows for an explicit representation of the decay rates in contrast to many theoretical results.

##### 2.1. Theoretical Results on Feedback Stabilization

Our goal is the stability analysis of the numerical solution to systems of conservation laws with linear source terms. We assume that the hyperbolic systems under investigation are given in the so-called characteristic form. Typically, linear systems of balance laws can be expressed in a characteristic form in a straightforward way while nonlinear systems can be transformed into a characteristic form by linearization with respect to a steady-state solution; see [3, 26] for examples and more details.

Let us consider the system of balance laws in characteristic form for :where is a matrix with constant entries, that is, , and for are the eigenvalues of the matrix with . The eigenvalues are either positive or negative, so we order these values such that for and for and split the matrix into and .

For characteristic variables , we have and such that So are the components of corresponding to the positive and negative eigenvalues of .

For the above system (1), linear feedback boundary conditions [7, 9] can be prescribed bywith matrix being a matrix, that is, , with blocks and . The name feedback originates from the fact that outputs of a system are fed back as inputs. Since the inputs and outputs of a system of balance laws are the boundaries, one speaks of boundary feedback control or boundary feedback stabilization, if the control is used to reach a steady-state (or equilibrium).

Note that boundary control problems might depend on additional controls and ; see applications in Section 3. Then, the aim is to find suitable controls such that for any initial data and boundary values (3) the solution of system (1) converges to the desired steady-state.

We are now concerned to analyze the exponential stability of system (1) under boundary conditions (3) with initial data for From theory we know that the linear hyperbolic system (1) and (3) is exponentially stable around the equilibrium state if there exist and such that, for every , the solution to the Cauchy problem satisfies

In the case of linear boundary conditions (3) the exponential stability of the equilibrium can be ensured by the following theoretical result [26, 33].

Let be the set of all real diagonal matrices with strictly positive diagonal elements. We define . Then the following theorem holds.

Theorem 1 (Theorem 2, [26]). *If , there exists such that, if , then the linear hyperbolic system (1) and (3) is exponentially stable.*

The key idea of the proof is to show that the time derivative of the Lyapunov functionwith weighting matrix and positive constants , for all , is negative definite along the solution of (1) and (3).

In [25] it is shown how this theoretical result can be reinterpreted from a numerical point of view. This means there exists a discrete version of Theorem 1 for the discretized system (1) and (3). The question is then under which additional conditions exponential stability of the discretized system can be ensured. In particular, corresponding conditions for the feedback matrix and the source term matrix need to be established. This question will be addressed in the next paragraph.

##### 2.2. Numerical Feedback Stabilization

We use a splitting method [34] to discretize system (1) and (3). This means in a first step the advection part is solved using an upwind scheme either in positive or in negative direction. In a second step, an explicit Euler scheme is applied to resolve the source term equation .

Therefore, we discretize for any component , with time step size and space step size , that is, for every discrete time , component , and cell index . We consider the cells and as ghost cells to describe the boundary values. So is the cell width of an equidistant spatial grid and the number of cells in the domain . The time step size is then chosen according to the CFL condition; that is,Note that for stiff problems there might be an additional time step restriction due to the explicit Euler discretization of the source term equation.

Hence, we end up with the following numerical scheme to solve system (1) and (3) for every cell index and time :withwhere is the element in the th row and th column of the matrix and is the th row of the source term matrix. As one can recognize, (7) and (9) are the numerical approximation to the advection and source term problem and (10) and (11) describe the boundary conditions. The initialization is given by (12).

Since we are interested in boundary feedback control from a numerical point of view, we need to define a discrete Lyapunov function of (5). This is necessary to analyze the exponential stability of the discretized system under linear boundary conditions.

Therefore we consider a discrete candidate Lyapunov function of the formwith evaluated at grid point and time step and positive real numbers and

Now, we are able to present the discrete version of Theorem 1 based on the numerical discretization (7)–(12) on the spatial domain for a positive finite time horizon . The main result taken from [25] can be summarized as follows.

Theorem 2. *The discretized system (7)–(12) is exponentially stable; that is, for and , if there exist and such that**(C1) the matricesare negative definite, with defined as (C2) the source term is restricted to Additionally, the decay rate can be explicitly given byand converges to for for all .*

For the proof, the same idea is used as in Theorem 1. The major difference is to show that the discrete time derivative of the discrete Lyapunov function (14) is negative; that is, , for the numerical discretization (7)–(12).

In the next section we apply the results of Theorem 2 to networks [17, 26–29] represented by directed graphs. Networked problems, where the dynamics on edges is governed by systems of balance laws, are linked at intersections and can be therefore treated as coupled boundary value problems. Hence, the results from above can be extended in a rigorous way.

#### 3. Application to Networks

A network is modeled by a directed graph and can be viewed as a spatial one-dimensional domain. To define hyperbolic balance laws on a network one has to introduce coupling conditions, which are a special form of boundary conditions at intersections. The associated feedback control problem is then to design control laws at the coupling interfaces to stabilize the system.

For application purposes, we use the results presented in Section 2 on two concrete network examples. One example is the so-called telegrapher’s equation, an equation described by a hyperbolic balance law and normally used to model transmission lines. The other example is a model for road traffic, the so-called LWR traffic flow network model. In this context, we will explain how a feedback control can be designed that small measurement errors are dampened.

In our numerical study, we discuss several scenarios and demonstrate the properties of Theorem 2 for small and larger networks. This means, in particular, to check the assumptions (C1) and (C2) for different network topologies.

Let us now start how the results from Section 2 fit into the context of networks.

*Definition 3. *A network is described by a directed graph , where is a set of nodes and is a set of edges. We represent an edge by an interval with uniform length .

To keep the notation simple, we assume that each edge (or interval ) can be mapped to the domain .

On each edge, we consider a system of balance laws (1). These systems are coupled by so-called coupling conditions prescribed as boundary values. For the coupling conditions we use the conservation of flow with inputs at the nodes. These inputs correspond to our control variables (cf. comments after (3)). So for a system of balance laws such as (1) the coupling has the following componentwise form: and here are external control inputs at fixed edges . So with this, the whole system in vector notation readswith

One can describe the controls as so-called feedback controls. A feedback law designs the controls in such a way that they depend on the outflow at the boundaries. We prescribe linear feedback boundary conditions for (21) as with matrices , and getThese are the full network equations including boundary conditions as coupling and as feedback.

In the following, we present two applications, in particular electrical transmission lines and vehicular traffic networks. For the first example, we describe how the network equations look like and which conditions must be ensured to obtain exponential stability while for the second example we explain how errors in the data can be dampened with the presented approach.

##### 3.1. Telegrapher’s Equation on Networks

We are interested in a feedback control law for electrical transmission lines. Those settings are normally represented as networks, where on each edge the so-called telegrapher’s equation is solved. The telegrapher’s equations are a simplification of the Maxwell equation and include signal losses on transmission lines; see [29] for the modeling details.

Normally, the variables current and voltage are the conservation quantities that are described by the telegrapher’s equation. They are expressed in terms of characteristic variables as and Due to the losses alongside the lines, we have an additional reaction term depending on the positive constants of , , and (resistance, inductance, and admittance) and (capacitance). The equation written in characteristic variables on edge is thenfor , For a directed graph with edges we are able to rewrite (24) in a compact form provided that is independent of . We collect , and introduce the matricesand with The coupling condition is rewritten in terms of , where , , with . Then, the network equations can be written in a compact form such as (21) and (24).

*Remark 4. *Note that the time step size should be additionally restricted in such a way that the numerical scheme (1) and (3) is still valid and assumption (C2) in Theorem 2 is satisfied; that is, is negative definite for all edges . This is the case for (cf. [29] for the details of the proof).

###### 3.1.1. Numerical Results

For all computations we choose the length of all edges as one over two, that is, , and also discretize the edges with an equidistant grid. We use the final time horizon . If some parameters change throughout the section, it will be stated. For simplicity, we skip the index whenever the situation is clear and all parameters are fixed to the same value for all edges.

*Example 5. *As a first example we consider the network shown in Figure 1.

We intend to introduce the network framework introduced above step by step. On each edge , we have two quantities and . We first choose and for all . If , the telegrapher’s equation reduces to the linear wave equation. Also note that and . We control the inflow of at node and the inflow of at node depending on the outflow at and , respectively; that is, and . The matrices and are composed ofThe corresponding eigenvalues of and (cf. (C1) in Theorem 2) are for , and for . If we assume that for all and insert the values for the spatial discretization , , on the strips and , we getWe remark that the determination of the eigenvalues is dependent on the network discretization, that is, the mapping of the edges to spatial intervals. For the discretization used here, we obtain two eigenvalues that are zero; that is, . However, these eigenvalues will not influence the stability result since they are independent of the feedback parameters , . So we have to ensure that the matrices and are negative definite for the feedback relevant eigenvalues , We claim which is the case ifSince it should also hold that and .

For the numerical computations we set as initial values , . We will vary the values for and and expect a rapidly dropping discrete Lyapunov function (14) for small and and vice versa. This behavior can be observed in Table 1. Since the CFL constant is equal to one and also the grid size is very small, that is, , the decay rate converges to for due to (19).