Data-Driven Predictive Control of Interconnected Systems Using the Koopman Operator

Abstract

Interconnected systems are widespread in modern technological systems. Designing a reliable control strategy requires modeling and analysis of the system, which can be a complicated, or even impossible, task in some cases. However, current technological developments in data sensing, processing, and storage make data-driven control techniques an appealing alternative solution. In this work, a design methodology of a decentralized control strategy is developed for interconnected systems based only on local and interconnection time series. Then, the optimization problem associated with the predictive control design is defined. Finally, an extension to interconnected systems coupled through their input signals is discussed. Simulations of two coupled Duffing oscillators, a bipedal locomotion model, and a four water tank system show the effectiveness of the approach.

1. Introduction

Many problems in engineering involve coupled dynamics, and designing a reliable control strategy for these systems requires an adequate understanding of the underlying dynamics. Interconnected systems are present in several areas of science and technology, from robotics to power grids and social systems [1,2,3]. In these examples, there may be several subsystems that are interconnected in one way or another, and handling the interconnections presents significant challenges in the analysis and design of control strategies. Traditionally, the derivation of state equations comes from applying physical principles to generate control-oriented models. However, due to the complexity of the system, modeling its interactions is a difficult, or even impossible, task in some cases.

Nevertheless, current technological developments in data sensing, processing, and storage, together with the large amount of information available, make data-driven control techniques an appealing alternative solution [4]. Early approaches to address the control synthesis of systems with unknown dynamical coupling relied on the framework of adaptive control. In [5], the authors presented a decentralized adaptive control for estimating dynamical couplings among subsystems. This adaptive approach dealt with a class of nonlinear systems, and it was able to compensate the interconnection with only local measurements. The assumption of a matching condition was relaxed, transforming the systems into a strict feedback form and validating the boundedness of the tracking error through classical Lyapunov techniques. Similarly, in [6], the authors proposed parameter estimation in the feedback control to compensate dynamical interactions. Adaptivity was based on radial basis neural networks, and a guarantee on the asymptotic tracking error was provided without knowledge of bounds on coupling dynamics. In [7], the authors considered time delay in a model reference adaptive control scheme.

However, most of the aforementioned approaches require satisfying the condition of persistent excitation and also detailed knowledge of the interconnection structure, which can be difficult in practice. On the other hand, taking advantage of a data-driven approach is appealing, and the Koopman operator has gained attention for its ability to evolve nonlinear dynamics in a function-space (lifted space of observables) in a linear form. The main issue with this approach is that the lifted dynamics act in an infinite-dimensional space [8]. Some form of truncation will therefore be required, implying a compromise between accuracy and dimensionality. Recently, the extended dynamic mode decomposition algorithm (EDMD) has been shown to be successful in obtaining a finite approximation of linear operators on the subspace spanned by a set of basis functions. Therefore, traditional techniques for the control synthesis have been extended to the Koopman operator framework [9]. Additionally, the Koopman theory provides solutions for dynamical systems that include coupled dynamics. For instance, the authors of [10] presented a methodology for reduced-order modeling using the Koopman theory in stochastic nonlinear coupled dynamical systems. The authors developed an empirical data-driven methodology to reduce the system. Alternatively, the authors in [3] presented a model reduction of agent-based systems using the Koopman generator. An interesting example of the data-driven Koopman operator identification and its relation with networked systems was presented in [11]. The authors used the operator generator to approximate the vector field of a network with nonlinear couplings, providing in turn the network structure. Continuing with the relationship between the Koopman operator and network topology, in [12,13], the authors used the dynamic mode decomposition (DMD) algorithm to make an approximation of the underlying topology. Others variants for high-dimensional systems are focused on efficient computation. For example, reference [14] investigated the EDMD based on Cholesky decomposition for dealing with high-dimensional systems and evaluated its performance with a network of coupled oscillators and a large-scale power system.

To date (and to the best of our knowledge), Koopman-based lifting techniques have not been developed yet to control coupled nonlinear systems. The objective of this study is, therefore, to present a methodology to design a decentralized control for interconnected systems based on local and interconnection data time series (experimental measurements or simulation of an existing simulator in any format). The main contributions of this work are threefold: (a) a local control-oriented model is developed based on data without knowledge of local dynamics or interconnection, (b) the optimization problem associated with the predictive control design is formulated, and (c) an identification procedure for systems, which are coupled through their inputs, is proposed.

The remainder of this paper is structured as follows: Section 2 introduces the Koopman operator and its finite dimensional approximation. In Section 3, the local control-oriented model of interconnected systems is derived, while a decentralized predictive control strategy is developed in Section 4. In Section 5, the control strategy is applied to three representative case studies, and the simulations results are discussed. Finally, Section 6 concludes the work and presents some future directions.

2. Preliminaries

This section presents a brief description of the general framework of the Koopman operator for data-driven control, followed by its finite numerical approximation.

2.1. The Koopman Operator for Control Systems

Consider a discrete-time dynamical system

$$x(t+1)=T\left(x\right(t),u(t\left)\right),$$

(1)

where $x\left(t\right)\in {\mathbb{R}}^n$ and $u\left(t\right)\in {\mathbb{R}}^m$, $T:{\mathbb{R}}^n\times {\mathbb{R}}^m\to {\mathbb{R}}^n$ is assumed to be at least ${\mathcal{C}}^1$, and $t\in {\mathbb{Z}}_{\ge 0}$. The Koopman operator has gained attention for its ability to evolve nonlinear dynamics in a function-space, i.e., the lifted space of observables, in a linear form. Given the vector space of measurement functions called observables $\mathcal{F}$ with observable functions $\psi \in \mathcal{F}$, we can define the Koopman operator for the system with inputs $\mathcal{K}:\mathcal{F}\to \mathcal{F}$ as 

$$\left[\mathcal{K}\psi \right]\left(x\right(t),u(t\left)\right)=\psi \left(T\right(x\left(t\right),u\left(t\right)),u(t+1\left)\right).$$

(2)

As mentioned above, the Koopman operator is a linear representation of nonlinear dynamics in a lifted space, with the compromise that it is infinite-dimensional, making it intractable to implement directly in an engineering application. Hence, it is necessary to obtain a finite-dimensional approximation. The next section briefly presents the method to obtain the finite dimensional approximation.

2.2. Finite Dimensional Approximation

In line with the data-driven nature of the algorithm, it is assumed that enough snapshot pairs of data are available for the computation, where each pair corresponds to $y=T(x,u)$ with

$$\begin{array}{c}\hfill \mathbf{X}=[x_0,\dots ,x_M],\mathbf{Y}=[y_0,\dots ,y_M],\end{array}$$

(3)

and the vector of data inputs is denoted by $\mathbf{U}=[u_1,\cdots ,u_M]$. The objective is to build a linear model in lifted space as follows,

$${\mathbf{\Psi }}_y=A{\mathbf{\Psi }}_x+B\mathbf{U},$$

(4)

where ${\mathbf{\Psi }}_{(·)}$ are basis functions, $A\in {\mathbb{R}}^{K\times K}$, $B\in {\mathbb{R}}^{K\times m}$. It is, therefore, possible to determine the linear operators A and B that best fit the evolution of observables by solving the following optimization problem

$$\underset{A,B}{min}{\parallel {\mathbf{\Psi }}_y-A{\mathbf{\Psi }}_x-B\mathbf{U}\parallel }_F$$

(5)

where F stands for the Frobenius norm [15,16]. The minimization problem gives a linear solution in the lifted space, where the evolution of the lifted state ${\mathbf{\Psi }}_y={[{\psi }_1\left(y\right),\dots ,{\psi }_K\left(y\right)]}^{\top }$, is the lifted state, i.e., ${\mathbf{\Psi }}_x={[{\psi }_1\left(x\right),\dots ,{\psi }_K\left(x\right)]}^{\top }$, weighted by the linear operator A plus the effect of the input. To retrieve the original states from the lifted space, we solve the least square problem as follows

$$\underset{C}{min}{\parallel \mathbf{X}-C{\mathbf{\Psi }}_x\parallel }_F$$

(6)

where $C\in {\mathbb{R}}^{n\times K}$; alternatively, C may be computed explicitly as

$$C=\mathbf{X}{\mathbf{\Psi }}_x^{†},$$

(7)

where † denotes the Pseudo-inverse. The selection of a set of basis functions can be arbitrary but the optimal choice to approximate the Koopman operator associated with complex systems is still a challenge. Possible choices are sets of orthogonal polynomials [17], orthogonal polynomials with functional embeddings [18], and radial-basis functions [11], among others. In the following, radial basis functions are selected in view of their performance in predicting the states, which is crucial in the MPC approach. For practical purposes, the observable functions may also contain the state variables; therefore, the solution to matrix C can be obtained in the form $C=[\mathit{I},\mathbf{0}]$.

3. Problem Statement

It is assumed that the composite system (1) can be decomposed as

$$x_i(t+1)=f_i\left(x_i\left(t\right)\right)+B_i{u}_i\left(t\right)+\sum _{j\in {\mathcal{N}}_i}{g}_{ij}\left(x_j\left(t\right)\right),$$

(8)

where $f_i$ and $g_{ij}$ are nonlinear maps for local and interconnected dynamics, ${\mathcal{N}}_i$ are the neighboring systems that exchange energy with the local system, and B is the local input matrix. Let $x_i\in {\mathbb{R}}^{n_i}$ be the local states, $x_j\in {\mathbb{R}}^{n_j}$ the states of a particular neighbor, $u_i\in {\mathbb{R}}^{m_i}$ control input, $n={\sum }_i{n}_i$, and $m={\sum }_i{m}_i$. The models of the component systems are unknown, but sets of measurements of inputs and local and coupled states are available. Our goal is to design a control strategy that considers the nonlinear dynamics with only the available data. The Koopman operator framework is used to find a local and finite linear approximation to model the nonlinear dynamics. Additionally, it is assumed that the linearization of each subsystem is controllable, and the coupling dynamics are Lipschitz bounded. Under the previous considerations, a linear representation in the lifted space is obtained with the following structure

$$\begin{array}{cc}\hfill {\Psi }_{x_i(t+1)}& ={\mathbf{K}}_i{\Psi }_{x_i\left(t\right)}+{\mathbf{B}}_i{\mathbf{U}}_i\left(t\right)+z_j\left(t\right)\hfill \end{array}$$

(9)

$$\begin{array}{cc}\hfill z_j\left(t\right)& =\sum _{j\in {\mathcal{N}}_j}{\mathbf{K}}_{ij}{\Psi }_{x_j\left(t\right)}\hfill \end{array}$$

(10)

where $z_j$ represents the lifted state of the interconnection dynamics. In the next section, the construction of the linear representation is detailed.

3.1. EDMD for Interconnected Systems

In order to build the linear representation, it is assumed that data are available from the inputs to the local system, the measurements of the local state, and the value of the neighboring system states, which are interconnected with the local system. The notation for this data tuple is ${\{x_i^n,y_i^n,x_j^n,y_j^n,u_i^n\}}_{n=1}^M$. The data in the lifted space is given by ${\mathbf{\Psi }}_{x_i}=[{\psi }_1\left(x_i\right),\cdots ,{\psi }_k\left(x_i\right)]\top$, ${\mathbf{\Psi }}_{y_i}=[{\psi }_1\left(y_i\right),\cdots ,{\psi }_k\left(y_i\right)]\top$, ${\mathbf{\Psi }}_{x_j}=[{\psi }_1\left(x_j\right),\cdots ,{\psi }_k\left(x_j\right)]\top$, and the input matrix $\mathbf{U}$ defined in Section 2.2. The operator approximation matrices come from the solution of the following least-squares problem

$$\begin{array}{c}\hfill \underset{{\mathbf{K}}_i{\mathbf{K}}_{ij}{\mathbf{B}}_i}{min}\parallel {\mathbf{\Psi }}_{y_i}-{\mathbf{K}}_i{\mathbf{\Psi }}_{x_i}-\sum _{j\in {\mathcal{N}}_j}{\mathbf{K}}_{ij}{\mathbf{\Psi }}_{x_j}-{\mathbf{B}}_i{\mathbf{U}}_i\parallel {}_F,\end{array}$$

(11)

with the explicit solution determined by

$$\begin{array}{c}\hfill [{\mathbf{K}}_i,{\mathbf{K}}_{ij},{\mathbf{B}}_i]={\mathbf{\Psi }}_{y_i}{[{\mathbf{\Psi }}_{x_i},{\mathbf{\Psi }}_{x_j},{\mathbf{U}}_i]}^{†}.\end{array}$$

(12)

Sometimes the time-series may be affected by noise, and the least-square optimization problem is solved with regularized terms as

$$\begin{array}{c}\hfill \underset{{\mathbf{K}}_i{\mathbf{K}}_{ij}{\mathbf{B}}_i}{min}\parallel {\mathbf{\Psi }}_{y_i}-{\mathbf{K}}_i{\mathbf{\Psi }}_{x_i}-\sum _{j\in {\mathcal{N}}_j}{\mathbf{K}}_{ij}{\mathbf{\Psi }}_{x_j}-{\mathbf{B}}_i{\mathbf{U}}_i{\parallel }_F+{\alpha }_1\parallel {\mathbf{K}}_i{\parallel }_1+{\alpha }_2{\parallel {\mathbf{K}}_{ij}\parallel }_1\end{array}$$

(13)

where ${\alpha }_{1,2}$ are tuning parameters. With the main concepts of the Koopman operator and its finite-dimensional approximation, we have the basis to develop the data-driven control design in the next section.

4. Data-Driven Koopman-Based Control Design

In this section, a data-driven control design is developed based on the approximation of the Koopman operator. Hence, the linear representation is used for prediction and a decentralized model predictive control problem is solved.

Consider the following cost function

$$J_i(z_i\left(k\right),u_i\left(k\right))=\sum _{k=1}^{H_p-1}(\parallel C{z}_i{\left(k\right)\parallel }_{\mathbf{Q}}^2+\parallel u_i\left(k\right){\parallel }_R^2),$$

(14)

where $H_p$ is the time horizon. The operator ${\parallel ·\parallel }_{\mathbf{Q}}^2$ denotes the weighted Eucledian norm defined for any vector z and weighting matrix $\mathbf{Q}$ as ${\parallel z\parallel }_{\mathbf{Q}}^2=z^{\top }\mathbf{Q}z$, and $\mathbf{Q}$ and R denote the positive definite weighting matrices for the lifted state and input vectors, respectively. Next, the optimal control problem is defined as

$$\begin{array}{cc}\hfill \underset{{\mathbf{u}}_i,{\mathbf{z}}_{\mathbf{i}}}{min}& J_i(z_i\left(k\right),u_i\left(k\right))\hfill \end{array}$$

(15)

$$\begin{array}{cc}\hfill \mathrm{s}.\mathrm{t}& z_i(k+1)={\mathbf{K}}_i{z}_i\left(k\right)+{\mathbf{B}}_i{u}_i\left(k\right)+z_j\left(k\right)\hfill \\ \hfill & z_j\left(k\right)=\sum _{{\mathcal{N}}_j}{\mathbf{K}}_{ij}{z}_j\left(k\right)\hfill \end{array}$$

(16)

$$\begin{array}{cc}\hfill & u_i^{min}\le u_i\left(k\right)\le u_i^{max}\hfill \end{array}$$

(17)

$$\begin{array}{cc}\hfill & x_i^{min}\le C{z}_i\left(k\right)\le x_i^{max}\hfill \\ \hfill & z_i\left(0\right)={\mathbf{\Psi }}_{x_i}\left(x_i\left(0\right)\right)\hfill \\ \hfill & z_j\left(0\right)={\mathbf{\Psi }}_{x_j}\left(x_j\left(0\right)\right),\hfill \end{array}$$

(18)

where the constraints (16)–(18) are defined for all $k\in \{0,\cdots ,H_p-1\}$, and the finite-dimensional approximation of the Koopman operator based on matrices (12) is used. The control signals are defined as ${\mathbf{u}}_i=(u_i\left(0\right),\cdots ,u_i(H_p-1))$ and the lifted state variables as ${\mathbf{z}}_i=(z_i\left(0\right),\cdots ,z_i(H_p-1))$.

In Figure 1, the structure of interconnection of each subsystem and the feedback control are depicted. It should be noticed that the coupling among subsystems is physical and not through the optimization cost. The control strategy is summarized in Algorithm 1, which consists of a first step with the initialization of the lifted data by the observables and parameters such as the time horizon $H_p$. Then, the finite-dimensional approximation of the Koopman operator is obtained using the method presented in Section 3.1 using (11) and (12). Finally, the optimal control signal $u_i^{*}\left(k\right)$ is applied to system (8) in the original state variables by means of the matrix C obtained in (7).

Algorithm 1: Decentralized MPC

Result: Optimal control signals $u_i^{*}$ Set $H_p$,$\mathbf{Q}$, R, ${\mathbf{\Psi }}_{x_i}$, ${\mathbf{\Psi }}_{x_j}$.

Interconnection through the Input

Another important interconnection configuration is now considered, where the coupling among subsystems is also given by the input control. Specifically, we consider the class of systems described by

$$x_i(t+1)=f_i\left(x_i\left(t\right)\right)+B_i{u}_i\left(t\right)+\sum _{j\in {\mathcal{N}}_i}{B}_{ij}\left(u_j\left(t\right)\right),$$

(19)

where it is observed that the interaction or coupling between subsystems is through the control signals of the other subsystems. The data driven-linear representation of this structure is the following

$$\begin{array}{cc}\hfill z_i(t+1)& ={\mathbf{K}}_i{z}_i\left(t\right)+{\mathbf{B}}_i{u}_i\left(t\right)+w_j\left(t\right)\hfill \end{array}$$

(20)

$$\begin{array}{cc}\hfill w_j\left(t\right)& =\sum _{j\in {\mathcal{N}}_j}{\mathbf{B}}_{ij}{\Psi }_{u_j\left(t\right)}\hfill \end{array}$$

(21)

The implementation of the methodology developed in Section 3.1 is straightforward. The buffer of neighbor’s input data ${\{u_i^n,u_j^n\}}_{n=1}^M$, and the vector of observables are required, and the optimization problem (13) can be solved.

In the next section, several simulation cases are presented to illustrate the theoretical results introduced in this section.

5. Simulations

This section presents numerical examples using the data-driven predictive control based on the Koopman operator framework and its finite-dimensional representation. All the simulations codes were run using MATLAB^® on a computer with 16 GB of RAM and a 3.8 GHz Intel Core i7 processor. The optimization problems were solved by qpOASES [19]. For the finite-dimensional approximation of the linear control-oriented model, 100 thin plate spline radial basis functions with centers generated randomly with uniform distribution at ${\mathbf{x}}_i$ were used. Each radial basis function was defined by ${\mathbf{\Psi }}_x={[{\psi }_1\left(x\right),\dots ,{\psi }_K\left(x\right)]}^{\top }$, where ${\psi }_i\left(x\right)=\Vert x-x_i{\Vert }_{2}log\left(\right|x-x_i\left|\right)$. Three different numerical examples are presented to illustrate the proposed decentralized MPC in problems of diverse complexity and nonlinearity.

5.1. Two Duffing Oscillators

As a first numerical example, consider an interconnection of two duffing oscillators, where the parameters yield a dynamic behavior with two basins of attraction. The free-body diagram of a mechanical duffing oscillator is depicted in Figure 2, and the differential equations that describe this system are given by,

$$\begin{array}{cc}\hfill {\dot{x}}_{11}=& x_{12}\hfill \end{array}$$

(22)

$$\begin{array}{cc}\hfill {\dot{x}}_{12}=& -{\delta }_1{x}_{12}-{\beta }_1{x}_{11}-{\alpha }_1{x}_{11}^3+{\eta }_{21}{x}_{21}+B_1{u}_1\hfill \end{array}$$

(23)

$$\begin{array}{cc}\hfill {\dot{x}}_{21}=& x_{22}\hfill \end{array}$$

(24)

$$\begin{array}{cc}\hfill {\dot{x}}_{22}=& -{\delta }_2{x}_{22}-{\beta }_2{x}_{21}-{\alpha }_2{x}_{21}^3+{\eta }_{12}{x}_{11}+B_2{u}_2\hfill \end{array}$$

(25)

where $x_{i1}$, $x_{i2}$ are the position and velocity of the mass i. The viscosity damping is denoted by ${\delta }_i$, and ${\beta }_i$, ${\alpha }_i$ are both constants representing the stiffness of the spring, respectively. The local input of each system is given by $u_i=F_i\left(t\right)$, more detail for an ideal interconnection can be found in [20].

In this example, ${\delta }_i=0.5$, ${\beta }_i=-1$, ${\alpha }_i=1$, and the interconnection parameter ${\eta }_i=0.5$. To obtain the local linear representation, we used one-step time-series data with ${10}^5$ random initial conditions uniformly distributed in the domain $[-0.8,0.8]\times [-0.8,0.8]$ with sample time $T_s=0.01$. The input matrix was approximated using ${10}^4$ random inputs over the interval $[-1,1]$. Notice that each subsystem required neighborhood state measurements to identify the coupling dynamics. The radial basis function centers were uniformly distributed in the domain. To evaluate the performance of the interconnected model, we compared the uncontrolled system states using an input $u_i=0.5cos\left(t\right)$. The comparison between the linear representation and the original system for the initial condition $x_0=[0.31,-0.34,-0.33,0.31]$ is shown in Figure 3. It is apparent that including the coupling term in the model improved the prediction capacity. The objective of the control problem was to stabilize the system in the origin. The cost function weights were chosen as $\mathbf{Q}=I$, $H_p=50$, and $R=0.01$, and the input constraint was $-1\le u\le 1$ for the two subsystems. The closed-loop local phase-plane and input signals are presented in Figure 4.

5.2. Bipedal Robot Locomotion Model

As a second example, consider the thigh and knee dynamics of a walking bipedal robot locomotion modeled through three coupled Van der Pol oscillators [21].

$$\begin{array}{cc}\hfill \ddot{{\theta }_1}& =0.1(1-5.25{\theta }_1^2)\dot{{\theta }_1}-{\theta }_1+u_1\hfill \end{array}$$

(26)

$$\begin{array}{cc}\hfill \ddot{{\theta }_2}& =0.01(1-6070{({\theta }_2-{\theta }_{2e})}^2)\dot{{\theta }_2}-4({\theta }_2-{\theta }_{2e})\hfill \\ \hfill & +0.057{\theta }_1\dot{{\theta }_1}+0.1({\dot{\theta }}_2-{\dot{\theta }}_3)+u_2\hfill \end{array}$$

(27)

$$\begin{array}{cc}\hfill \ddot{{\theta }_3}& =0.01(1-192{({\theta }_3-{\theta }_{3e})}^2)\dot{{\theta }_3}-4({\theta }_3-{\theta }_{3e})\hfill \\ \hfill & +0.057{\theta }_1\dot{{\theta }_1}+0.1({\dot{\theta }}_3-{\dot{\theta }}_3)+u_3\hfill \end{array}$$

(28)

where ${\theta }_i$ and ${\dot{\theta }}_i$ are the angle (rad) and angular velocity (rad/s), respectively. The states of each subsystem are denoted $x_i={\theta }_i-{\theta }_{ie}$, and $u_i$ is the control input. A diagram which represents the interconnection in the state space is depicted in Figure 5. In this example, we used locally one-step time-series data with ${10}^4$ random initial conditions uniformly distributed in the domain $[-\pi /2,\pi /2]\times [-1,1]$ with sample time $T_s=0.01$. The input matrix was approximated using ${10}^3$ random inputs over the interval $[-1,1]$. The radial basis function centers were again uniformly distributed in the domain. The control design was subject to the following constraints $-\frac{\pi }{2}\le {\theta }_i-{\theta }_{ie}\le \frac{\pi }{2}$ and $|u_i|\le 1$. The desired constant angles for $i\in \{1,2,3\}$ were ${\theta }_{1e}=0$, ${\theta }_{2e}=-\pi /12$, and ${\theta }_{3e}=\pi /6$. The parameters for each decentralized local control were $\mathbf{Q}=I$, $H_p=100$, $R=0.1$. The effectiveness of the control algorithm can be observed in Figure 6, where the angles achieved the desired values, and the input constraints were satisfied.

5.3. Four Water Tanks

Finally, the four water tank control benchmark, illustrated in Figure 7, where two subsystems are coupled by means of the input, is considered. The coupled dynamics is modeled as

$$\begin{array}{cc}\hfill {\dot{x}}_{11}& =-\frac{a_1}{S}\sqrt{2g{x}_{11}}+\frac{a_3}{S}\sqrt{2g{x}_{12}}+\frac{{\gamma }_1}{S}{q}_1\hfill \end{array}$$

(29)

$$\begin{array}{cc}\hfill {\dot{x}}_{12}& =-\frac{a_3}{S}\sqrt{2g{x}_{12}}+\frac{1-{\gamma }_2}{S}{q}_2\hfill \end{array}$$

(30)

$$\begin{array}{cc}\hfill {\dot{x}}_{21}& =-\frac{a_2}{S}\sqrt{2g{x}_{21}}+\frac{a_4}{S}\sqrt{2g{x}_{22}}+\frac{{\gamma }_2}{S}{q}_2\hfill \end{array}$$

(31)

$$\begin{array}{cc}\hfill {\dot{x}}_{22}& =-\frac{a_4}{S}\sqrt{2g{x}_{22}}+\frac{1-{\gamma }_1}{S}{q}_1\hfill \end{array}$$

(32)

where S is the cross-section area, and $a_i$ is the discharge constant of each tank. Let $q_i$ and ${\gamma }_i$ with $i\in \{1,2\}$ denote the flow and the ratio of the three-way valve of pump i, respectively, and g the gravitational acceleration. The system parameters can be found in [22]. In this example, we used one-step time-series data with ${10}^4$ random initial conditions uniformly distributed in the domain $[0,1]\times [0,1]$ with sample time $T_s=0.001$. The local input matrix was approximated using ${10}^5$ random inputs over the interval $[-1,1]$ and ${10}^3$ random points of neighbor’s input. In this control problem, level tracking had to be achieved. The objective function in the data-driven predictive control is defined by

$$J_i(z_i\left(k\right),u_i\left(k\right))=\sum _{k=1}^{H_p-1}(\parallel C{z}_i\left(k\right)-{\mathrm{Ref}}_i{\parallel }_{\mathbf{Q}}^2+\parallel u_i\left(k\right){\parallel }_R^2),$$

(33)

where the matrix C selects only the first state. The cost function weights were chosen as $\mathbf{Q}=I$, $H_p=100$, $R=0.001$, the level constraint was $0.1\le x_{i1}\le 1$, and the flow maximum was $q_i=0.5$ for the two subsystems. The results of the simulation are shown in Figure 8. The decentralized data-driven controller was able to regulate the system to the given set-points while guaranteeing that the state and input constraints were satisfied.

6. Conclusions

Based on the framework of the Koopman operator and the extended dynamic mode decomposition, we proposed a method to design a model-free and decentralized control for interconnected systems. The method uses only local and interconnection series-time of synthetic data or measurements. No prior knowledge of local dynamics or interconnection is required. We have shown how to construct the optimization problem associated with the predictive control design, and an extension to systems coupled through their inputs was developed. Some simulation results were shown to illustrate the effectiveness of the proposed approach. Future work will focus on more complex systems involving more subsystems and on testing different basis functions in order to improve the performance of the linear predictor model.