Location via proxy:   [ UP ]  
[Report a bug]   [Manage cookies]                

HTML conversions sometimes display errors due to content that did not convert correctly from the source. This paper uses the following packages that are not yet supported by the HTML conversion tool. Feedback on these issues are not necessary; they are known and are being worked on.

  • failed: arydshln

Authors: achieve the best HTML results from your LaTeX submissions by following these best practices.

License: arXiv.org perpetual non-exclusive license
arXiv:2404.06746v1 [eess.SY] 10 Apr 2024

Data-driven parallel Koopman subsystem modeling and distributed moving horizon state estimation for large-scale nonlinear processes

Xiaojie Lia𝑎{}^{a}start_FLOATSUPERSCRIPT italic_a end_FLOATSUPERSCRIPT, Song Bob𝑏{}^{b}start_FLOATSUPERSCRIPT italic_b end_FLOATSUPERSCRIPT, Xuewen Zhanga𝑎{}^{a}start_FLOATSUPERSCRIPT italic_a end_FLOATSUPERSCRIPT, Yan Qinc𝑐{}^{c}start_FLOATSUPERSCRIPT italic_c end_FLOATSUPERSCRIPT, Xunyuan Yina,𝑎{}^{a,}start_FLOATSUPERSCRIPT italic_a , end_FLOATSUPERSCRIPT

a𝑎{}^{a}start_FLOATSUPERSCRIPT italic_a end_FLOATSUPERSCRIPTSchool of Chemistry, Chemical Engineering and Biotechnology, Nanyang Technological University

62 Nanyang Drive, 637459, Singapore

b𝑏{}^{b}start_FLOATSUPERSCRIPT italic_b end_FLOATSUPERSCRIPTDepartment of Chemical & Materials Engineering, University of Alberta

Edmonton, AB T6G 1H9, Canada

c𝑐{}^{c}start_FLOATSUPERSCRIPT italic_c end_FLOATSUPERSCRIPTEngineering Product Development Pillar, The Singapore University of Technology and Design

8 Somapah Road, 487372, Singapore
Corresponding author: X. Yin. Tel: +65-6316-8746. Email: xunyuan.yin@ntu.edu.sg
Abstract

In this work, we consider a state estimation problem for large-scale nonlinear processes in the absence of first-principles process models. By exploiting process operation data, both process modeling and state estimation design are addressed within a distributed framework. By leveraging the Koopman operator concept, a parallel subsystem modeling approach is proposed to establish interactive linear subsystem process models in higher-dimensional subspaces, each of which correlates with the original nonlinear subspace of the corresponding process subsystem via a nonlinear mapping. The data-driven linear subsystem models can be used to collaboratively characterize and predict the dynamical behaviors of the entire nonlinear process. Based on the established subsystem models, local state estimators that can explicitly handle process operation constraints are designed using moving horizon estimation. The local estimators are integrated via information exchange to form a distributed estimation scheme, which provides estimates of the unmeasured/unmeasurable state variables of the original nonlinear process in a linear manner. The proposed framework is applied to a chemical process and an agro-hydrological process to illustrate its effectiveness and applicability. Good open-loop predictability of the linear subsystem models is confirmed, and accurate estimates of the process states are obtained without requiring a first-principles process model.

Keywords: Data-driven state estimation, parallel subsystem identification, large-scale nonlinear process, Koopman operator, distributed moving horizon estimation.

Introduction

Large-scale and tightly integrated processes have been widely used in various industries. Structurally flexible and computationally efficient advanced control systems are needed to appropriately monitor and intervene in the operation of these processes to ensure operating safety, and to achieve efficient and sustainable process operation/production [1, 2, 3, 4]. Despite the promise, the development and implementation of advanced control systems for large-scale complex processes have been challenging. Specifically, due to the high nonlinearity of the process dynamics and rapid increase in the scale and the level of structural complexity, conventional decision-making paradigms are no longer favorable, as they can lead to difficulties in systems organization and maintenance, low fault tolerance of a decision-making system, and intractable online computation. The distributed framework is a very promising alternative for the monitoring and control of complex processes. A distributed decision-making system consists of local decision-making agents, each of which typically only needs to handle a sub-problem instead of the plant-wide problem for the entire process (known as a subsystem of the process). The local agents collaborate with each other based on real-time information exchange to assess the current process operating status and make informed decisions on the most appropriate control actions [1, 2, 3, 5, 6, 7]. The distributed architecture holds the promise to be used for developing flexible and scalable decision-making schemes for large-scale industrial systems. A complete decision-making system typically involves distributed state estimation and distributed control. Distributed estimation is used to provide real-time estimates of key process state variables that are not measured by the existing sensing instruments, and is an enabler for comprehensive process monitoring and informed decision-making of distributed control. In this work, we focus on the distributed state estimation problem for large-scale nonlinear processes, which has received less research attention as compared to the distributed (model predictive) control counterpart [8, 9, 10].

In the existing literature, there have been some results on distributed state estimation for large-scale systems/processes. Distributed estimation algorithms based on deterministic state observers were proposed [11, 12, 13]. In existing literature [14, 15, 16, 17, 18], consensus-driven distributed Kalman filtering [14, 15, 16, 17] and extended Kalman filtering algorithms [18] were proposed on the basis of sensor networks. To provide optimal estimates in the presence of constraints, distributed moving horizon estimation approaches were proposed [19, 20]. More relevant results on distributed state estimation can be found in [21, 22, 23]. It is worth mentioning that the aforementioned methods assume high-fidelity (first-principles) system/process models are available, and good estimation performance is premised on the high accuracy of the existing dynamic model. However, from an application perspective, it can be challenging to establish a first-principles model with accurate model parameters due to the availability of very limited first-principles knowledge – this is especially critical when dealing with large-scale complex nonlinear processes. In addition, for nonlinear processes, the use of a nonlinear dynamic model will complicate the design and implementation of local estimators. Therefore, to bridge the gap between theoretical research and practice in distributed estimation for complex processes, it is necessary to model the process dynamics accurately without requiring the availability of full first-principles knowledge and accurate values of the associated parameters.

One of the alternative solutions to describing the dynamics of a nonlinear process is to establish a data-driven model by leveraging the concept of the Koopman operator. A Koopman operator is a linear representation of the dynamics of a nonlinear process that can be used to predict the future evolution of the process states. Specifically, by performing a nonlinear coordinate transformation, the temporal evolution of the nonlinear process states can be characterized by a higher-dimensional linear state-space model governed by the Koopman operator [24, 25, 26]. However, in real applications, the exact Koopman operator associated with a nonlinear process may have infinite dimensions, and can be very difficult to find. This has motivated the exploration of a finite-dimensional approximation of the exact Koopman operator. There have been algorithms that can be used to calculate an approximation of the exact Koopman operator by leveraging historical process data; for example, the Dynamic Mode Decomposition [27] and its non-trivial extension – Extended Dynamic Mode Decomposition [28], the Generalized Laplace Analysis [29], and the Galerkin Method [30].

In the context of process control, Koopman identification is an emerging enabler for developing data-driven decision-making and monitoring schemes for nonlinear processes/systems using linear control and estimation methods [31, 32, 33, 34, 35]. Based on Koopman identification, a Lyapunov-based model predictive control (MPC) for nonlinear chemical processes was proposed in Narasingam and Kwon [31]. In Korda and Mezić [32], a Koopman predictor was established as the basis of linear MPC for nonlinear systems. In Son et al. [33], a hybrid Koopman MPC approach was proposed. In this work, multiple local Koopman models were developed to account for the local dynamics of a batch pulping process under different operating conditions. There have been fewer results on Koopman state estimation as compared to the control counterpart. A data-driven Koopman Kalman filtering was developed for power systems [36]. In Yin et al. [37], a data-driven moving horizon estimation algorithm was proposed based on an identified Koopman linear model for constrained nonlinear processes. More results on Koopman-based control and estimation can be found in [26, 34, 38, 39, 40, 41, 42, 43, 44]. Considering its merits in terms of being able to bypass first-principles modeling, the Koopman concept will be leveraged in this work for data-driven distributed state estimation. However, the aforementioned methods focused on control and estimation within the centralized framework, and the associated Koopman-based identification methods are not favorable for developing distributed estimation schemes. The reasons are mainly threefold: 1) the identification of one single Koopman linear model that describes the dynamics of the entire process will become much more challenging as the scale and complexity of the processes increase; 2) the nonlinear coordinate transformation incurs a significant increase in the dimension of the state-space, and the number of observables of a single Koopman model for a complex process with many state variables can be excessively large; 3) as subsystem models are needed to design the local estimators within a distributed scheme, if a plant-wide model is to be used, then it needs to be partitioned into subsystem models appropriately before designing local estimators, which brings in additional difficulty and incurs more research efforts. For data-driven state estimation for large-scale processes within a distributed framework, it is indeed favorable to identify Koopman-based interactive subsystem models in parallel as the basis of the design of local estimators.

In this work, we address a data-driven distributed state estimation problem for large-scale constrained nonlinear processes by proposing a parallel subsystem modeling approach and a data-driven distributed moving horizon estimation algorithm. Specifically, a data-driven parallel subsystem model identification method is proposed to establish subsystem models that interact with each other based on the concept of the Koopman operator, which collaboratively describes the dynamics of the entire process. Both subsystem states and inputs are projected onto a higher-dimensional space. In the proposed parallel subsystem modeling method, for each subsystem, the states of its interactive subsystems (of which the states have direct influences on the subsystem being considered) are treated as known inputs. Based on the identified subsystem models, local estimators are designed as moving horizon estimators, each of which is capable of handling both the local dynamics and the interactive dynamics of the associated subsystem, in the presence of constraints on subsystem states, inputs, and disturbances. The proposed framework consisting of parallel subsystem identification and distributed moving horizon estimation is illustrated via applications to a chemical process and an argo-hydrological process.

Preliminaries

Notation

𝕀𝕀\mathbb{I}blackboard_I is a set of consecutive positive integers 𝕀={1,,m}𝕀1𝑚\mathbb{I}=\{1,\ldots,m\}blackboard_I = { 1 , … , italic_m }. \circ denotes the operator of the function composition. col(x1,,xm)colsubscript𝑥1subscript𝑥𝑚\text{col}\left(x_{1},\ldots,x_{m}\right)col ( italic_x start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT , … , italic_x start_POSTSUBSCRIPT italic_m end_POSTSUBSCRIPT ) is a column vector consisting of x1,x2,,xmsubscript𝑥1subscript𝑥2subscript𝑥𝑚x_{1},x_{2},\ldots,x_{m}italic_x start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT , italic_x start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT , … , italic_x start_POSTSUBSCRIPT italic_m end_POSTSUBSCRIPT. *superscript\mathbb{Z}^{*}blackboard_Z start_POSTSUPERSCRIPT * end_POSTSUPERSCRIPT is a set of positive integers. 𝐈nsubscript𝐈𝑛{\bf{I}}_{n}bold_I start_POSTSUBSCRIPT italic_n end_POSTSUBSCRIPT represents an identity matrix of dimension n𝑛nitalic_n. 𝟎a×bsubscript0𝑎𝑏{\boldsymbol{0}}_{a\times b}bold_0 start_POSTSUBSCRIPT italic_a × italic_b end_POSTSUBSCRIPT denotes an a×b𝑎𝑏a\times bitalic_a × italic_b zero matrix. {x}absubscriptsuperscript𝑥𝑏𝑎\left\{x\right\}^{b}_{a}{ italic_x } start_POSTSUPERSCRIPT italic_b end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_a end_POSTSUBSCRIPT represents a column vector comprising a sequence of variable x𝑥xitalic_x from time instant a𝑎aitalic_a to instant b𝑏bitalic_b, i.e., x(a),x(a+1),,x(b)𝑥𝑎𝑥𝑎1𝑥𝑏x(a),x(a+1),\ldots,x(b)italic_x ( italic_a ) , italic_x ( italic_a + 1 ) , … , italic_x ( italic_b ). xP2subscriptsuperscriptnorm𝑥2𝑃\left\|x\right\|^{2}_{P}∥ italic_x ∥ start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_P end_POSTSUBSCRIPT denotes the square of the weighted Euclidean norm of vector x𝑥xitalic_x, which is computed as xP2=xTPxsubscriptsuperscriptnorm𝑥2𝑃superscript𝑥T𝑃𝑥\left\|x\right\|^{2}_{P}=x^{\text{T}}Px∥ italic_x ∥ start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_P end_POSTSUBSCRIPT = italic_x start_POSTSUPERSCRIPT T end_POSTSUPERSCRIPT italic_P italic_x. diag{A1,A2}diagsubscript𝐴1subscript𝐴2\text{diag}\big{\{}{A_{1},A_{2}\ldots}\big{\}}diag { italic_A start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT , italic_A start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT … } denotes a block diagonal matrix of which the blocks on the main diagonal are A1subscript𝐴1A_{1}italic_A start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT, A2subscript𝐴2A_{2}italic_A start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT, … x^(d|k)^𝑥conditional𝑑𝑘{\hat{x}}(d|k)over^ start_ARG italic_x end_ARG ( italic_d | italic_k ) represents an estimate of x(d)𝑥𝑑x(d)italic_x ( italic_d ) obtained at sampling instant k𝑘kitalic_k (kd𝑘𝑑k\geq ditalic_k ≥ italic_d).

Process description

Let us consider a class of general nonlinear processes of medium- to large-scales that can be divided into m𝑚mitalic_m subsystems; the dynamics of each subsystem i𝑖iitalic_i, i𝕀𝑖𝕀i\in\mathbb{I}italic_i ∈ blackboard_I, are assumed to be with the following form:

xi(k+1)subscript𝑥𝑖𝑘1\displaystyle x_{i}(k+1)italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k + 1 ) =fi(xi(k),Xi(k),ui(k))absentsubscript𝑓𝑖subscript𝑥𝑖𝑘subscript𝑋𝑖𝑘subscript𝑢𝑖𝑘\displaystyle=f_{i}\big{(}x_{i}(k),X_{i}(k),u_{i}(k)\big{)}= italic_f start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) , italic_X start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) , italic_u start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) ) (1a)
yi(k+1)subscript𝑦𝑖𝑘1\displaystyle y_{i}(k+1)italic_y start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k + 1 ) =hi(xi(k))absentsubscript𝑖subscript𝑥𝑖𝑘\displaystyle=h_{i}\big{(}x_{i}(k)\big{)}= italic_h start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) ) (1b)

where xi𝕏inxisubscript𝑥𝑖subscript𝕏𝑖superscriptsubscript𝑛subscript𝑥𝑖x_{i}\in{\mathbb{X}}_{i}\subset\mathbb{R}^{n_{x_{i}}}italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ∈ blackboard_X start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ⊂ blackboard_R start_POSTSUPERSCRIPT italic_n start_POSTSUBSCRIPT italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT end_POSTSUPERSCRIPT is the state vector of subsystem i𝑖iitalic_i; ui𝕌inuisubscript𝑢𝑖subscript𝕌𝑖superscriptsubscript𝑛subscript𝑢𝑖u_{i}\in\mathbb{U}_{i}\subset\mathbb{R}^{n_{u_{i}}}italic_u start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ∈ blackboard_U start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ⊂ blackboard_R start_POSTSUPERSCRIPT italic_n start_POSTSUBSCRIPT italic_u start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT end_POSTSUPERSCRIPT is a vector of known inputs that affect the dynamics of subsystem i𝑖iitalic_i directly, including manipulated inputs and/or known disturbances to the process; Xi𝕏¯inXisubscript𝑋𝑖subscript¯𝕏𝑖superscriptsubscript𝑛subscript𝑋𝑖X_{i}\in{\bar{\mathbb{X}}}_{i}\subset\mathbb{R}^{n_{X_{i}}}italic_X start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ∈ over¯ start_ARG blackboard_X end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ⊂ blackboard_R start_POSTSUPERSCRIPT italic_n start_POSTSUBSCRIPT italic_X start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT end_POSTSUPERSCRIPT is a vector of the states of the subsystems that have direct interaction with subsystem i𝑖iitalic_i; yinyisubscript𝑦𝑖superscriptsubscript𝑛subscript𝑦𝑖y_{i}\in\mathbb{R}^{n_{y_{i}}}italic_y start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ∈ blackboard_R start_POSTSUPERSCRIPT italic_n start_POSTSUBSCRIPT italic_y start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT end_POSTSUPERSCRIPT is the sensor measurements vector of subsystem i𝑖iitalic_i; fisubscript𝑓𝑖f_{i}italic_f start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT is a vector-valued nonlinear function characterizing the dependence of the future states of subsystem i𝑖iitalic_i on the current states of the local subsystem and the interactive subsystems and inputs to the local subsystem; hisubscript𝑖h_{i}italic_h start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT is the function of sensor measurements for subsystem i𝑖iitalic_i, i𝕀𝑖𝕀i\in\mathbb{I}italic_i ∈ blackboard_I. Note that set 𝕏¯isubscript¯𝕏𝑖\bar{\mathbb{X}}_{i}over¯ start_ARG blackboard_X end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT can be determined based on 𝕏jsubscript𝕏𝑗\mathbb{X}_{j}blackboard_X start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT, j𝕀ifor-all𝑗subscript𝕀𝑖\forall j\in\mathbb{I}_{i}∀ italic_j ∈ blackboard_I start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT, where 𝕀isubscript𝕀𝑖\mathbb{I}_{i}blackboard_I start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT (i𝕀𝑖𝕀i\in\mathbb{I}italic_i ∈ blackboard_I) denotes a set of the indices of the subsystems, of which the states directly affect the dynamics of subsystem i𝑖iitalic_i. For an illustrating purpose, if the state of subsystem i𝑖iitalic_i at the next sampling instant is dependent on the current states of subsystems 2222 and 5555 (i.e., Xisubscript𝑋𝑖X_{i}italic_X start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT in (1a) contain the states of subsystems 2222 and 5555), then 𝕀i={2,5}subscript𝕀𝑖25\mathbb{I}_{i}=\left\{2,5\right\}blackboard_I start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT = { 2 , 5 }.

It is worth mentioning that by incorporating all the subsystem states into an augmented vector, the dynamics of the entire process consisting of all the subsystems in (Process description) can be described by the following general form:

x(k+1)=𝑥𝑘1absent\displaystyle x(k+1)=italic_x ( italic_k + 1 ) = f(x(k),u(k))𝑓𝑥𝑘𝑢𝑘\displaystyle~{}f\big{(}x(k),u(k)\big{)}italic_f ( italic_x ( italic_k ) , italic_u ( italic_k ) ) (2)
y(k)=𝑦𝑘absent\displaystyle y(k)=italic_y ( italic_k ) = h(x(k))𝑥𝑘\displaystyle~{}h\big{(}x(k)\big{)}italic_h ( italic_x ( italic_k ) )

where x=col(x1,,xm)𝕏nx𝑥colsubscript𝑥1subscript𝑥𝑚𝕏superscriptsubscript𝑛𝑥x=\text{col}\left(x_{1},\ldots,x_{m}\right)\in\mathbb{X}\subset\mathbb{R}^{n_{% x}}italic_x = col ( italic_x start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT , … , italic_x start_POSTSUBSCRIPT italic_m end_POSTSUBSCRIPT ) ∈ blackboard_X ⊂ blackboard_R start_POSTSUPERSCRIPT italic_n start_POSTSUBSCRIPT italic_x end_POSTSUBSCRIPT end_POSTSUPERSCRIPT (with nx=nxisubscript𝑛𝑥subscript𝑛subscript𝑥𝑖n_{x}=\sum n_{x_{i}}italic_n start_POSTSUBSCRIPT italic_x end_POSTSUBSCRIPT = ∑ italic_n start_POSTSUBSCRIPT italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT) represents the process state vector; u𝑢uitalic_u is a vector of all the known inputs that are involved in uisubscript𝑢𝑖u_{i}italic_u start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT, i𝕀𝑖𝕀i\in\mathbb{I}italic_i ∈ blackboard_I; y=col(y1,,ym)𝑦colsubscript𝑦1subscript𝑦𝑚y=\text{col}\left(y_{1},\ldots,y_{m}\right)italic_y = col ( italic_y start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT , … , italic_y start_POSTSUBSCRIPT italic_m end_POSTSUBSCRIPT ) is the vector of all the sensor measurements; f=col(f1,,fm)𝑓colsubscript𝑓1subscript𝑓𝑚f=\text{col}\left(f_{1},\ldots,f_{m}\right)italic_f = col ( italic_f start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT , … , italic_f start_POSTSUBSCRIPT italic_m end_POSTSUBSCRIPT ) is an augmented nonlinear vector field characterizing the dynamics of the entire process; h=col(h1,,hm)colsubscript1subscript𝑚h=\text{col}\left(h_{1},\ldots,h_{m}\right)italic_h = col ( italic_h start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT , … , italic_h start_POSTSUBSCRIPT italic_m end_POSTSUBSCRIPT ) describes the dependence of all the sensor measurements on the process state.

Problem formulation

In this work, we aim to address the state estimation problem for the nonlinear process of medium- to large-scales in (Process description) (or equivalently, (2)) by proposing a distributed scheme consisting of local estimators developed based on subsystem models. However, as mentioned, it may not be favorable to utilize the subsystem models in (Process description) directly for designing the local estimators due to: 1) the model structure and/or the values of the model parameters may not be completely known; 2) the development of local estimators based on nonlinear models are more challenging, and the implementation of nonlinear estimators is much more computationally complex than the case when linear models are used.

Based on the above observations, we aim first to construct data-driven linear subsystem models, which can be used to approximate the dynamics of the subsystems of the process in (Process description) in a lifted function space that is related to the original function space via a nonlinear lifting mapping. Specifically, the model to be identified for subsystem i𝑖iitalic_i, i𝕀𝑖𝕀i\in\mathbb{I}italic_i ∈ blackboard_I, is in the following form:

zi(k+1)=subscript𝑧𝑖𝑘1absent\displaystyle z_{i}(k+1)=italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k + 1 ) = Aiizi(k)+j𝕀iAijzj(k)+Biu~i(k)subscript𝐴𝑖𝑖subscript𝑧𝑖𝑘subscript𝑗subscript𝕀𝑖subscript𝐴𝑖𝑗subscript𝑧𝑗𝑘subscript𝐵𝑖subscript~𝑢𝑖𝑘\displaystyle~{}A_{ii}z_{i}(k)+\sum_{j\in\mathbb{I}_{i}}A_{ij}z_{j}(k)+B_{i}{% \tilde{u}}_{i}(k)italic_A start_POSTSUBSCRIPT italic_i italic_i end_POSTSUBSCRIPT italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) + ∑ start_POSTSUBSCRIPT italic_j ∈ blackboard_I start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT italic_A start_POSTSUBSCRIPT italic_i italic_j end_POSTSUBSCRIPT italic_z start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT ( italic_k ) + italic_B start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT over~ start_ARG italic_u end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) (3a)
yi(k)=subscript𝑦𝑖𝑘absent\displaystyle y_{i}(k)=italic_y start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) = Cizi(k)subscript𝐶𝑖subscript𝑧𝑖𝑘\displaystyle~{}C_{i}z_{i}(k)italic_C start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) (3b)
x^i(k)=subscript^𝑥𝑖𝑘absent\displaystyle{\hat{x}}_{i}(k)=over^ start_ARG italic_x end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) = Dizi(k)subscript𝐷𝑖subscript𝑧𝑖𝑘\displaystyle~{}D_{i}z_{i}(k)italic_D start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) (3c)

where zinzisubscript𝑧𝑖superscriptsubscript𝑛subscript𝑧𝑖z_{i}\in\mathbb{R}^{n_{z_{i}}}italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ∈ blackboard_R start_POSTSUPERSCRIPT italic_n start_POSTSUBSCRIPT italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT end_POSTSUPERSCRIPT denotes the vector of lifted states of the linear model for subsystem i𝑖iitalic_i; u~inu~isubscript~𝑢𝑖superscriptsubscript𝑛subscript~𝑢𝑖{\tilde{u}}_{i}\in\mathbb{R}^{n_{{\tilde{u}}_{i}}}over~ start_ARG italic_u end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ∈ blackboard_R start_POSTSUPERSCRIPT italic_n start_POSTSUBSCRIPT over~ start_ARG italic_u end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT end_POSTSUPERSCRIPT denotes the vector of lifted inputs to the linear model for subsystem i𝑖iitalic_i; x^inxisubscript^𝑥𝑖superscriptsubscript𝑛subscript𝑥𝑖{\hat{x}}_{i}\in\mathbb{R}^{n_{x_{i}}}over^ start_ARG italic_x end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ∈ blackboard_R start_POSTSUPERSCRIPT italic_n start_POSTSUBSCRIPT italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT end_POSTSUPERSCRIPT denotes a prediction of the actual state of subsystem i𝑖iitalic_i based on the evolution of the identified linear model (3a); Aiinzi×nzisubscript𝐴𝑖𝑖superscriptsubscript𝑛subscript𝑧𝑖subscript𝑛subscript𝑧𝑖A_{ii}\in\mathbb{R}^{n_{z_{i}}\times n_{z_{i}}}italic_A start_POSTSUBSCRIPT italic_i italic_i end_POSTSUBSCRIPT ∈ blackboard_R start_POSTSUPERSCRIPT italic_n start_POSTSUBSCRIPT italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT × italic_n start_POSTSUBSCRIPT italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT end_POSTSUPERSCRIPT and Aijnzi×nzjsubscript𝐴𝑖𝑗superscriptsubscript𝑛subscript𝑧𝑖subscript𝑛subscript𝑧𝑗A_{ij}\in\mathbb{R}^{n_{z_{i}}\times n_{z_{j}}}italic_A start_POSTSUBSCRIPT italic_i italic_j end_POSTSUBSCRIPT ∈ blackboard_R start_POSTSUPERSCRIPT italic_n start_POSTSUBSCRIPT italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT × italic_n start_POSTSUBSCRIPT italic_z start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT end_POSTSUBSCRIPT end_POSTSUPERSCRIPT are subsystem matrices characterizing the dependence of future subsystem states on the local subsystem states and interactive subsystem states, respectively; Binzi×nu~isubscript𝐵𝑖superscriptsubscript𝑛subscript𝑧𝑖subscript𝑛subscript~𝑢𝑖B_{i}\in\mathbb{R}^{n_{z_{i}}\times n_{{\tilde{u}}_{i}}}italic_B start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ∈ blackboard_R start_POSTSUPERSCRIPT italic_n start_POSTSUBSCRIPT italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT × italic_n start_POSTSUBSCRIPT over~ start_ARG italic_u end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT end_POSTSUPERSCRIPT is an input matrix characterizing the dependence of the subsystem states on known inputs of subsystem i𝑖iitalic_i; Cinyi×nzisubscript𝐶𝑖superscriptsubscript𝑛subscript𝑦𝑖subscript𝑛subscript𝑧𝑖C_{i}\in\mathbb{R}^{n_{y_{i}}\times n_{z_{i}}}italic_C start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ∈ blackboard_R start_POSTSUPERSCRIPT italic_n start_POSTSUBSCRIPT italic_y start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT × italic_n start_POSTSUBSCRIPT italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT end_POSTSUPERSCRIPT and Dinxi×nzisubscript𝐷𝑖superscriptsubscript𝑛subscript𝑥𝑖subscript𝑛subscript𝑧𝑖D_{i}\in\mathbb{R}^{n_{x_{i}}\times n_{z_{i}}}italic_D start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ∈ blackboard_R start_POSTSUPERSCRIPT italic_n start_POSTSUBSCRIPT italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT × italic_n start_POSTSUBSCRIPT italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT end_POSTSUPERSCRIPT are output measurement and state-reconstruction matrices that need to be identified for subsystem i𝑖iitalic_i, i𝕀𝑖𝕀i\in\mathbb{I}italic_i ∈ blackboard_I.

With the data-driven subsystem models in the form of (Problem formulation), local state estimators will be designed and a distributed estimation scheme will be formed for online estimating actual states of the considered nonlinear process.

Remark 1.

Deterministic models are presented in this section to facilitate the description of the proposed subsystem model identification method. It is worth noting that the parallel subsystem model identification and distributed moving horizon estimation algorithm proposed in this work can be used to handle unknown disturbances and measurement noise.

Parallel Koopman identification of the linear subsystem models

In this section, a parallel subsystem modeling method is proposed for establishing data-driven linear subsystem models in the form of (Problem formulation) as needed for the development of a distributed estimation scheme.

The basis of Koopman operator

First, the concept and the basis of the Koopman operator for discrete nonlinear controlled processes, which are leveraged for developing the subsystem modeling approach are reviewed. Specifically, consider the following nonlinear controlled system:

x(k+1)=F(x(k),u(k))𝑥𝑘1𝐹𝑥𝑘𝑢𝑘x(k+1)=F\big{(}x(k),u(k)\big{)}italic_x ( italic_k + 1 ) = italic_F ( italic_x ( italic_k ) , italic_u ( italic_k ) ) (4)

where xnx𝑥superscriptsubscript𝑛𝑥x\in\mathbb{R}^{n_{x}}italic_x ∈ blackboard_R start_POSTSUPERSCRIPT italic_n start_POSTSUBSCRIPT italic_x end_POSTSUBSCRIPT end_POSTSUPERSCRIPT is the vector of process states, unu𝑢superscriptsubscript𝑛𝑢u\in\mathbb{R}^{n_{u}}italic_u ∈ blackboard_R start_POSTSUPERSCRIPT italic_n start_POSTSUBSCRIPT italic_u end_POSTSUBSCRIPT end_POSTSUPERSCRIPT is the vector of input vectors, and F:nx×nunx:𝐹superscriptsubscript𝑛𝑥superscriptsubscript𝑛𝑢superscriptsubscript𝑛𝑥F:\mathbb{R}^{n_{x}}\times\mathbb{R}^{n_{u}}\rightarrow\mathbb{R}^{n_{x}}italic_F : blackboard_R start_POSTSUPERSCRIPT italic_n start_POSTSUBSCRIPT italic_x end_POSTSUBSCRIPT end_POSTSUPERSCRIPT × blackboard_R start_POSTSUPERSCRIPT italic_n start_POSTSUBSCRIPT italic_u end_POSTSUBSCRIPT end_POSTSUPERSCRIPT → blackboard_R start_POSTSUPERSCRIPT italic_n start_POSTSUBSCRIPT italic_x end_POSTSUBSCRIPT end_POSTSUPERSCRIPT is a vector-valued nonlinear function.

For system (4), a Koopman operator can be constructed in an extended Euclidean space constituted by the Cartesian product of the space of the original process states x𝑥xitalic_x and the space of the known inputs u𝑢uitalic_u. Specifically, an augmented vector containing both the system states and the known inputs is created as follows:

𝒳=[xu]𝒳delimited-[]𝑥𝑢\mathcal{X}=\left[\begin{array}[]{l}x\\ u\end{array}\right]caligraphic_X = [ start_ARRAY start_ROW start_CELL italic_x end_CELL end_ROW start_ROW start_CELL italic_u end_CELL end_ROW end_ARRAY ] (5)

which constitutes the state vector of the extended space. The dynamics of 𝒳𝒳\mathcal{X}caligraphic_X be described by:

𝒳(k+1)=(𝒳(k)):=[F(x(k),u(k))𝒮u(k)]:=[F(x(k),u(k))u(k+1)]𝒳𝑘1𝒳𝑘assigndelimited-[]𝐹𝑥𝑘𝑢𝑘𝒮𝑢𝑘assigndelimited-[]𝐹𝑥𝑘𝑢𝑘𝑢𝑘1\mathcal{X}(k+1)=\mathcal{F}\left(\mathcal{X}(k)\right):=\left[\begin{array}[]% {c}F\left(x(k),u(k)\right)\\ \mathcal{S}u(k)\end{array}\right]:=\left[\begin{array}[]{c}F\left(x(k),u(k)% \right)\\ u(k+1)\end{array}\right]caligraphic_X ( italic_k + 1 ) = caligraphic_F ( caligraphic_X ( italic_k ) ) := [ start_ARRAY start_ROW start_CELL italic_F ( italic_x ( italic_k ) , italic_u ( italic_k ) ) end_CELL end_ROW start_ROW start_CELL caligraphic_S italic_u ( italic_k ) end_CELL end_ROW end_ARRAY ] := [ start_ARRAY start_ROW start_CELL italic_F ( italic_x ( italic_k ) , italic_u ( italic_k ) ) end_CELL end_ROW start_ROW start_CELL italic_u ( italic_k + 1 ) end_CELL end_ROW end_ARRAY ] (6)

where 𝒮𝒮\mathcal{S}caligraphic_S is a left shift operator, that is, 𝒮u(k)=u(k+1)𝒮𝑢𝑘𝑢𝑘1\mathcal{S}u(k)=u(k+1)caligraphic_S italic_u ( italic_k ) = italic_u ( italic_k + 1 ).

For the augmented system in (6), a linear dynamic model based on the Koopman operator concept is defined as [32, 34, 45]:

𝒦g=g𝒦𝑔𝑔\mathcal{K}g=g\circ\mathcal{F}caligraphic_K italic_g = italic_g ∘ caligraphic_F (7)

where g𝑔gitalic_g represents the observables of a higher dimensional function space 𝒢𝒢\mathcal{G}caligraphic_G, which is associated with the original extended space (which is constituted by the Cartesian product of the process states x𝑥xitalic_x and known inputs u𝑢uitalic_u) through a certain nonlinear mapping. 𝒦:𝒢𝒢:𝒦𝒢𝒢\mathcal{K}:\mathcal{G}\rightarrow\mathcal{G}caligraphic_K : caligraphic_G → caligraphic_G is the Koopman operator governing the time evolution of the observables in the lifted space 𝒢𝒢\mathcal{G}caligraphic_G. Specifically, for each sampling time k0𝑘0k\geq 0italic_k ≥ 0, (7) indicates that:

𝒦g(𝒳(k))=g(𝒳(k))=g(𝒳(k+1))𝒦𝑔𝒳𝑘𝑔𝒳𝑘𝑔𝒳𝑘1\mathcal{K}g\big{(}\mathcal{X}(k)\big{)}=g\circ\mathcal{F}\left(\mathcal{X}(k)% \right)=g\big{(}\mathcal{X}(k+1)\big{)}caligraphic_K italic_g ( caligraphic_X ( italic_k ) ) = italic_g ∘ caligraphic_F ( caligraphic_X ( italic_k ) ) = italic_g ( caligraphic_X ( italic_k + 1 ) ) (8)

One favorable property of the Koopman operator is its linearity despite the nonlinear dynamics of the original system, that is, for any pair of observables g1,g2𝒢subscript𝑔1subscript𝑔2𝒢g_{1},g_{2}\in\mathcal{G}italic_g start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT , italic_g start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT ∈ caligraphic_G, and for any scalar coefficients α,β𝛼𝛽\alpha,\beta\in\mathbb{R}italic_α , italic_β ∈ blackboard_R, it is satisfied that:

𝒦(αg1+βg2)𝒦𝛼subscript𝑔1𝛽subscript𝑔2\displaystyle\mathcal{K}\left(\alpha g_{1}+\beta g_{2}\right)caligraphic_K ( italic_α italic_g start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT + italic_β italic_g start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT ) =(αg1+βg2)absent𝛼subscript𝑔1𝛽subscript𝑔2\displaystyle=\left(\alpha g_{1}+\beta g_{2}\right)\circ\mathcal{F}= ( italic_α italic_g start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT + italic_β italic_g start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT ) ∘ caligraphic_F
=αg1+βg2absent𝛼subscript𝑔1𝛽subscript𝑔2\displaystyle=\alpha g_{1}\circ\mathcal{F}+\beta g_{2}\circ\mathcal{F}= italic_α italic_g start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT ∘ caligraphic_F + italic_β italic_g start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT ∘ caligraphic_F
=α𝒦g1+β𝒦g2absent𝛼𝒦subscript𝑔1𝛽𝒦subscript𝑔2\displaystyle=\alpha\mathcal{K}g_{1}+\beta\mathcal{K}g_{2}= italic_α caligraphic_K italic_g start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT + italic_β caligraphic_K italic_g start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT

Consequently, the Koopman operator provides a linear representation of the time evolution of a nonlinear process in a higher-dimensional space, which is typical of infinite-dimension [46].

EDMD-based parallel Koopman subsystem identification method

In this subsection, we present a Koopman-based parallel identification method that can be used to construct finite-dimensional linear subsystem models in the form of (Problem formulation) to approximate the dynamics of the subsystems of the entire large-scale nonlinear process.

To identify linear models for the m𝑚mitalic_m subsystems instead of building a single Koopman linear model in the form of (8) for the entire nonlinear process, one challenge is how to appropriately handle the subsystem interactions. In the nonlinear subsystem models in (Process description) that describe process dynamics, the dependence of the future state of subsystem i𝑖iitalic_i on the current states of the interactive subsystems Xi(k)subscript𝑋𝑖𝑘X_{i}(k)italic_X start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) can be viewed as known inputs. Accordingly, in the i𝑖iitalic_ith data-driven subsystem model (Problem formulation) that needs to be constructed via identification, the vectors of lifted states of the interactive subsystems (i.e., zjsubscript𝑧𝑗z_{j}italic_z start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT, j𝕀i𝑗subscript𝕀𝑖j\in\mathbb{I}_{i}italic_j ∈ blackboard_I start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT) can also be treated as known inputs that affect the future state of the local subsystem lifted state zisubscript𝑧𝑖z_{i}italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT. Considering the concept of Koopman operator described in the subsection titled “The basis of Koopman operator”, for the state vector xisubscript𝑥𝑖x_{i}italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT of each subsystem i𝑖iitalic_i, i𝕀𝑖𝕀i\in\mathbb{I}italic_i ∈ blackboard_I, in the form of (Process description), there exists a (possibly infinite-dimensional) space 𝒢isubscript𝒢𝑖\mathcal{G}_{i}caligraphic_G start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT, of which the observables can be used to describe the dynamics of the same subsystem using a Koopman operator in a linear state-space form as an analog of (8), with both interactive subsystem states Xisubscript𝑋𝑖X_{i}italic_X start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT and subsystem inputs uisubscript𝑢𝑖u_{i}italic_u start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT being treated as known inputs in Koopman identification. Specifically, construct an augmented state vector 𝒳isubscript𝒳𝑖\mathcal{X}_{i}caligraphic_X start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT being 𝒳i:=[xiT,XiT,uiT]Tassignsubscript𝒳𝑖superscriptsuperscriptsubscript𝑥𝑖Tsuperscriptsubscript𝑋𝑖Tsuperscriptsubscript𝑢𝑖TT\mathcal{X}_{i}:=\left[x_{i}^{\text{T}},X_{i}^{\text{T}},u_{i}^{\text{T}}% \right]^{\text{T}}caligraphic_X start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT := [ italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT start_POSTSUPERSCRIPT T end_POSTSUPERSCRIPT , italic_X start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT start_POSTSUPERSCRIPT T end_POSTSUPERSCRIPT , italic_u start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT start_POSTSUPERSCRIPT T end_POSTSUPERSCRIPT ] start_POSTSUPERSCRIPT T end_POSTSUPERSCRIPT, and its dynamics can be characterized by

𝒳i(k+1)=i(𝒳i(k)):=[fi(xi(k),Xi(k),ui(k))iXi(k)𝒮iui(k)]=[xi(k+1)Xi(k+1)ui(k+1)]subscript𝒳𝑖𝑘1subscript𝑖subscript𝒳𝑖𝑘assigndelimited-[]subscript𝑓𝑖subscript𝑥𝑖𝑘subscript𝑋𝑖𝑘subscript𝑢𝑖𝑘subscript𝑖subscript𝑋𝑖𝑘subscript𝒮𝑖subscript𝑢𝑖𝑘delimited-[]subscript𝑥𝑖𝑘1subscript𝑋𝑖𝑘1subscript𝑢𝑖𝑘1\mathcal{X}_{i}(k+1)=\mathcal{F}_{i}\left(\mathcal{X}_{i}(k)\right):=\left[% \begin{array}[]{c}f_{i}\left(x_{i}(k),X_{i}(k),u_{i}(k)\right)\\ \mathcal{R}_{i}X_{i}(k)\\ \mathcal{S}_{i}u_{i}(k)\end{array}\right]=\left[\begin{array}[]{c}x_{i}(k+1)\\ X_{i}(k+1)\\ u_{i}(k+1)\end{array}\right]caligraphic_X start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k + 1 ) = caligraphic_F start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( caligraphic_X start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) ) := [ start_ARRAY start_ROW start_CELL italic_f start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) , italic_X start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) , italic_u start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) ) end_CELL end_ROW start_ROW start_CELL caligraphic_R start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT italic_X start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) end_CELL end_ROW start_ROW start_CELL caligraphic_S start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT italic_u start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) end_CELL end_ROW end_ARRAY ] = [ start_ARRAY start_ROW start_CELL italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k + 1 ) end_CELL end_ROW start_ROW start_CELL italic_X start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k + 1 ) end_CELL end_ROW start_ROW start_CELL italic_u start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k + 1 ) end_CELL end_ROW end_ARRAY ] (9)

where isubscript𝑖\mathcal{R}_{i}caligraphic_R start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT and 𝒮isubscript𝒮𝑖\mathcal{S}_{i}caligraphic_S start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT are left shift operators, that is, iXi(k)=Xi(k+1)subscript𝑖subscript𝑋𝑖𝑘subscript𝑋𝑖𝑘1\mathcal{R}_{i}X_{i}(k)=X_{i}(k+1)caligraphic_R start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT italic_X start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) = italic_X start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k + 1 ) and 𝒮iui(k)=ui(k+1)subscript𝒮𝑖subscript𝑢𝑖𝑘subscript𝑢𝑖𝑘1\mathcal{S}_{i}u_{i}(k)=u_{i}(k+1)caligraphic_S start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT italic_u start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) = italic_u start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k + 1 ). This way, from a theoretical viewpoint, it is possible to use a Koopman operator to characterize the underlying local dynamical behavior of each nonlinear subsystem i𝑖iitalic_i in the form of (1a) using a linear state-space model established via lifting the original extended space containing 𝒳isubscript𝒳𝑖\mathcal{X}_{i}caligraphic_X start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT.

However, the dimension of each space 𝒢isubscript𝒢𝑖\mathcal{G}_{i}caligraphic_G start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT, i𝕀𝑖𝕀i\in\mathbb{I}italic_i ∈ blackboard_I, and the corresponding exact Koopman operator are typically infinite. In addition, it can be very challenging to find exact Koopman operators for the subsystems for most nonlinear processes. These facts significantly limit the use of exact Koopman operators in the analysis of nonlinear processes using linear techniques. Instead of seeking the exact Koopman operator, we propose an alternative method to construct finite-dimensional approximated Koopman operators for the subsystems of the nonlinear process. In the existing literature, approximation algorithms that may be used to construct an approximated Koopman operator by leveraging batch process data include the Dynamic Mode Decomposition [27] and its non-trivial extension – Extended Dynamic Mode Decomposition (EDMD) [28], the Generalized Laplace Analysis [29], and the Galerkin Method [30]. The EDMD algorithm is an effective approach to approximate Koopman operators in finite-dimensional state space. This approach utilizes a nonlinear mapping to transform the original state vectors into lifted states in a higher-dimensional space, where the dynamics are described by using a Koopman-based linear model. The EDMD method is capable of handling processes exhibiting strong nonlinearity [28], and is leveraged to develop the Koopman-based parallel subsystem identification method of the present work.

Specifically, to establish a linear model for each subsystem i𝑖iitalic_i, i𝕀𝑖𝕀i\in\mathbb{I}italic_i ∈ blackboard_I in the form of (Problem formulation), the first step of the proposed method is to establish a nonlinear mapping from the original subsystem state vector xisubscript𝑥𝑖x_{i}italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT to the subsystem lifted states on finite-dimensional space 𝒢i~~subscript𝒢𝑖\tilde{\mathcal{G}_{i}}over~ start_ARG caligraphic_G start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_ARG, which is a subspace of 𝒢isubscript𝒢𝑖\mathcal{G}_{i}caligraphic_G start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT and is spanned by a number of user-specified linearly independent state basis functions. On this subspace, nzisubscript𝑛subscript𝑧𝑖n_{z_{i}}italic_n start_POSTSUBSCRIPT italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT state lifting functions, denoted by ϕpi(xi):𝕏i:subscriptsuperscriptitalic-ϕ𝑖𝑝subscript𝑥𝑖subscript𝕏𝑖\phi^{i}_{p}(x_{i}):\mathbb{X}_{i}\rightarrow\mathbb{R}italic_ϕ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_p end_POSTSUBSCRIPT ( italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) : blackboard_X start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT → blackboard_R, p=1,,nzifor-all𝑝1subscript𝑛subscript𝑧𝑖\forall p=1,\ldots,n_{z_{i}}∀ italic_p = 1 , … , italic_n start_POSTSUBSCRIPT italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT, are defined for each subsystem i𝑖iitalic_i, i𝕀𝑖𝕀i\in\mathbb{I}italic_i ∈ blackboard_I. In addition, as an analog of the nonlinear mapping for each subsystem state vector xisubscript𝑥𝑖x_{i}italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT, a nonlinear mapping that transforms the original inputs uisubscript𝑢𝑖u_{i}italic_u start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT into lifted inputs u~isubscript~𝑢𝑖{\tilde{u}}_{i}over~ start_ARG italic_u end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT of the linear subsystem model (Problem formulation) is determined, and the corresponding input lifting functions are represented by δqi(ui):𝕌i,q=1,,nu~i:subscriptsuperscript𝛿𝑖𝑞subscript𝑢𝑖formulae-sequencesubscript𝕌𝑖𝑞1subscript𝑛subscript~𝑢𝑖\delta^{i}_{q}(u_{i}):\mathbb{U}_{i}\rightarrow\mathbb{R},q=1,\ldots,n_{{% \tilde{u}}_{i}}italic_δ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_q end_POSTSUBSCRIPT ( italic_u start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) : blackboard_U start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT → blackboard_R , italic_q = 1 , … , italic_n start_POSTSUBSCRIPT over~ start_ARG italic_u end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT where nu~isubscript𝑛subscript~𝑢𝑖n_{{\tilde{u}}_{i}}italic_n start_POSTSUBSCRIPT over~ start_ARG italic_u end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT represents the dimension of the vector of lifted inputs u~isubscript~𝑢𝑖\tilde{u}_{i}over~ start_ARG italic_u end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT to subsystem i𝑖iitalic_i.

For each subsystem i𝑖iitalic_i in (1a), by leveraging the nonlinear mappings associated with the subsystems states and inputs, respectively, we define the vector of observables as follows:

ψi(𝒳i)=[ϕi(xi)\hdashline[2pt/2pt]Φi(Xi)\hdashline[2pt/2pt]δi(ui)]superscript𝜓𝑖subscript𝒳𝑖delimited-[]superscriptitalic-ϕ𝑖subscript𝑥𝑖\hdashlinedelimited-[]2𝑝𝑡2𝑝𝑡superscriptΦ𝑖subscript𝑋𝑖\hdashlinedelimited-[]2𝑝𝑡2𝑝𝑡superscript𝛿𝑖subscript𝑢𝑖\psi^{i}(\mathcal{X}_{i})=\left[\begin{array}[]{c}\phi^{i}(x_{i})\\ \hdashline[2pt/2pt]\Phi^{i}(X_{i})\\ \hdashline[2pt/2pt]\delta^{i}(u_{i})\end{array}\right]italic_ψ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT ( caligraphic_X start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) = [ start_ARRAY start_ROW start_CELL italic_ϕ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT ( italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) end_CELL end_ROW start_ROW start_CELL [ 2 italic_p italic_t / 2 italic_p italic_t ] roman_Φ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT ( italic_X start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) end_CELL end_ROW start_ROW start_CELL [ 2 italic_p italic_t / 2 italic_p italic_t ] italic_δ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT ( italic_u start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) end_CELL end_ROW end_ARRAY ] (10)

where ϕisuperscriptitalic-ϕ𝑖\phi^{i}italic_ϕ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT is the vector of state lifting functions for subsystem i𝑖iitalic_i (i.e., ϕi=[ϕ1i,,ϕnzii]Tsuperscriptitalic-ϕ𝑖superscriptsubscriptsuperscriptitalic-ϕ𝑖1subscriptsuperscriptitalic-ϕ𝑖subscript𝑛subscript𝑧𝑖T\phi^{i}=\big{[}\phi^{i}_{1},\ldots,\phi^{i}_{n_{z_{i}}}\big{]}^{\text{T}}italic_ϕ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT = [ italic_ϕ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT , … , italic_ϕ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_n start_POSTSUBSCRIPT italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT end_POSTSUBSCRIPT ] start_POSTSUPERSCRIPT T end_POSTSUPERSCRIPT); δisuperscript𝛿𝑖\delta^{i}italic_δ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT is the vector of input lifting functions for subsystem i𝑖iitalic_i (i.e., δi=[δ1i,,δnu~ii]Tsuperscript𝛿𝑖superscriptsubscriptsuperscript𝛿𝑖1subscriptsuperscript𝛿𝑖subscript𝑛subscript~𝑢𝑖T\delta^{i}=\big{[}\delta^{i}_{1},\ldots,\delta^{i}_{n_{{\tilde{u}}_{i}}}\big{]% }^{\text{T}}italic_δ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT = [ italic_δ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT , … , italic_δ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_n start_POSTSUBSCRIPT over~ start_ARG italic_u end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT end_POSTSUBSCRIPT ] start_POSTSUPERSCRIPT T end_POSTSUPERSCRIPT); ΦisuperscriptΦ𝑖\Phi^{i}roman_Φ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT is a concatenated vector consisting of state lifting functions of the interactive subsystems of subsystem i𝑖iitalic_i (i.e., ϕj(xj),j𝕀isuperscriptitalic-ϕ𝑗subscript𝑥𝑗for-all𝑗subscript𝕀𝑖\phi^{j}(x_{j}),\forall j\in\mathbb{I}_{i}italic_ϕ start_POSTSUPERSCRIPT italic_j end_POSTSUPERSCRIPT ( italic_x start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT ) , ∀ italic_j ∈ blackboard_I start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT). This way, both the states of the interactive subsystems Xisubscript𝑋𝑖X_{i}italic_X start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT and the known inputs uisubscript𝑢𝑖u_{i}italic_u start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT, which are treated as external inputs to subsystem i𝑖iitalic_i are taken into account in the observables.

Based on the defined vector of observables ψisuperscript𝜓𝑖\psi^{i}italic_ψ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT, a finite-dimensional Koopman linear model can be established for each subsystem i𝑖iitalic_i, i𝕀𝑖𝕀i\in\mathbb{I}italic_i ∈ blackboard_I, as follows:

ψii(𝒳i(k))=ψi(𝒳i(k+1))𝒦~iψi(𝒳i(k))superscript𝜓𝑖subscript𝑖subscript𝒳𝑖𝑘superscript𝜓𝑖subscript𝒳𝑖𝑘1subscript~𝒦𝑖superscript𝜓𝑖subscript𝒳𝑖𝑘\psi^{i}\circ\mathcal{F}_{i}\left(\mathcal{X}_{i}(k)\right)=\psi^{i}\left(% \mathcal{X}_{i}({k+1})\right)\approx{\tilde{\mathcal{K}}}_{i}\psi^{i}\left(% \mathcal{X}_{i}(k)\right)italic_ψ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT ∘ caligraphic_F start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( caligraphic_X start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) ) = italic_ψ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT ( caligraphic_X start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k + 1 ) ) ≈ over~ start_ARG caligraphic_K end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT italic_ψ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT ( caligraphic_X start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) ) (11)

where ψinisuperscript𝜓𝑖superscriptsubscript𝑛𝑖\psi^{i}\in\mathbb{R}^{n_{i}}italic_ψ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT ∈ blackboard_R start_POSTSUPERSCRIPT italic_n start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUPERSCRIPT (ninxi+nXi+nuimuch-greater-thansubscript𝑛𝑖subscript𝑛subscript𝑥𝑖subscript𝑛subscript𝑋𝑖subscript𝑛subscript𝑢𝑖n_{i}\gg n_{x_{i}}+n_{X_{i}}+n_{u_{i}}italic_n start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ≫ italic_n start_POSTSUBSCRIPT italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT + italic_n start_POSTSUBSCRIPT italic_X start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT + italic_n start_POSTSUBSCRIPT italic_u start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT) is the subsystem observable vector; 𝒦i~ni×ni~subscript𝒦𝑖superscriptsubscript𝑛𝑖subscript𝑛𝑖\tilde{\mathcal{K}_{i}}\in\mathbb{R}^{n_{i}\times n_{i}}over~ start_ARG caligraphic_K start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_ARG ∈ blackboard_R start_POSTSUPERSCRIPT italic_n start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT × italic_n start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUPERSCRIPT is a finite-dimensional Koopman operator governing the time evolution of the subsystem observables in (11).

The Koopman linear model (11) will be utilized to establish a subsystem model for subsystem i𝑖iitalic_i in the form of (Problem formulation), of which the lifted states are defined as ϕisuperscriptitalic-ϕ𝑖\phi^{i}italic_ϕ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT and the lifted inputs are δisuperscript𝛿𝑖\delta^{i}italic_δ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT; that is, zi=ϕisubscript𝑧𝑖superscriptitalic-ϕ𝑖z_{i}=\phi^{i}italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT = italic_ϕ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT, with initial condition zi(0)=ϕi(x(0))subscript𝑧𝑖0superscriptitalic-ϕ𝑖𝑥0z_{i}(0)=\phi^{i}\left(x(0)\right)italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( 0 ) = italic_ϕ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT ( italic_x ( 0 ) ), and u~i=δisubscript~𝑢𝑖superscript𝛿𝑖{\tilde{u}}_{i}=\delta^{i}over~ start_ARG italic_u end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT = italic_δ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT. This is based on the consideration that ϕisuperscriptitalic-ϕ𝑖\phi^{i}italic_ϕ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT is only dependent on the state of the local subsystem, such that it may be easier to reconstruct the actual states xisubscript𝑥𝑖x_{i}italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT of each subsystem i𝑖iitalic_i in (1a) from the trajectories of the lifted states of (3a), which can be obtained based on the trajectories of the observables of each Koopman linear model (11).

The next step is to calculate an approximated Koopman operator 𝒦i~~subscript𝒦𝑖\tilde{\mathcal{K}_{i}}over~ start_ARG caligraphic_K start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_ARG for each subsystem i𝑖iitalic_i, i𝕀𝑖𝕀i\in\mathbb{I}italic_i ∈ blackboard_I. It is desirable to make the approximation in (11) as accurate as possible. Therefore, the objective is to find the Koopman operator that minimizes the discrepancy between the actual values of the observables and the predicted values of the observables of each subsystem. Considering the finite number of available data samples, for each subsystem i𝑖iitalic_i, i𝕀𝑖𝕀i\in\mathbb{I}italic_i ∈ blackboard_I, we aim to find a 𝒦~isubscript~𝒦𝑖\mathcal{\tilde{K}}_{i}over~ start_ARG caligraphic_K end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT that provides the smallest averaged L2limit-fromsuperscript𝐿2L^{2}-italic_L start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT -norm prediction error over the entire time window, which is described as follows:

𝒦~i=argmin𝒦~ik=1Nψi(𝒳i(k+1))𝒦~iψi(𝒳i(k))2subscript~𝒦𝑖subscriptsubscript~𝒦𝑖superscriptsubscript𝑘1𝑁superscriptnormsuperscript𝜓𝑖subscript𝒳𝑖𝑘1subscript~𝒦𝑖superscript𝜓𝑖subscript𝒳𝑖𝑘2\mathcal{{\tilde{K}}}_{i}=\arg\mathop{\min}\limits_{\mathcal{{\tilde{K}}}_{i}}% \sum_{k=1}^{N}\left\|\psi^{i}\left(\mathcal{X}_{i}(k+1)\right)-\mathcal{{% \tilde{K}}}_{i}\psi^{i}\left(\mathcal{X}_{i}(k)\right)\right\|^{2}over~ start_ARG caligraphic_K end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT = roman_arg roman_min start_POSTSUBSCRIPT over~ start_ARG caligraphic_K end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT ∑ start_POSTSUBSCRIPT italic_k = 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_N end_POSTSUPERSCRIPT ∥ italic_ψ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT ( caligraphic_X start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k + 1 ) ) - over~ start_ARG caligraphic_K end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT italic_ψ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT ( caligraphic_X start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) ) ∥ start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT (12)

where N𝑁Nitalic_N denotes the total number of data samples available for analysis.

Considering the components of the subsystem observables, the approximated Koopman operator 𝒦~isubscript~𝒦𝑖\mathcal{{\tilde{K}}}_{i}over~ start_ARG caligraphic_K end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT can be expressed in the form of blocks as follows:

𝒦~i=[AiiA_ijB_i\hdashline[2pt/2pt]***\hdashline[2pt/2pt]\hdashline[2pt/2pt]*

*

*
]
\mathcal{{\tilde{K}}}_{i}=\left[\begin{array}[]{c;{2pt/2pt}c;{2pt/2pt}c;{2pt/2% pt}c;{2pt/2pt}c}A_{ii}&\ldots&A_{ij}&\ldots&B_i\\ \hdashline[2pt/2pt]*&*&\ldots&\ldots&*\\ \hdashline[2pt/2pt]\vdots&\vdots&&&\vdots\\ \hdashline[2pt/2pt]*&*&\ldots&\ldots&*\end{array}\right]over~ start_ARG caligraphic_K end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT = [ start_ARRAY start_ROW start_CELL italic_A start_POSTSUBSCRIPT italic_i italic_i end_POSTSUBSCRIPT end_CELL start_CELL … end_CELL start_CELL A_ij end_CELL start_CELL … end_CELL start_CELL B_i end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL end_ROW start_ROW start_CELL [ 2 italic_p italic_t / 2 italic_p italic_t ] * end_CELL start_CELL * end_CELL start_CELL … end_CELL start_CELL … end_CELL start_CELL * end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL end_ROW start_ROW start_CELL [ 2 italic_p italic_t / 2 italic_p italic_t ] ⋮ end_CELL start_CELL ⋮ end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL ⋮ end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL end_ROW start_ROW start_CELL [ 2 italic_p italic_t / 2 italic_p italic_t ] * end_CELL start_CELL end_CELL start_CELL … end_CELL start_CELL … end_CELL start_CELL * end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL end_ROW end_ARRAY ]
(13)

In the first row of the matrix in (13), the elements between Aiisubscript𝐴𝑖𝑖A_{ii}italic_A start_POSTSUBSCRIPT italic_i italic_i end_POSTSUBSCRIPT and Bisubscript𝐵𝑖B_{i}italic_B start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT are constituted by blocks Aijsubscript𝐴𝑖𝑗A_{ij}italic_A start_POSTSUBSCRIPT italic_i italic_j end_POSTSUBSCRIPT, j𝕀ifor-all𝑗subscript𝕀𝑖\forall j\in\mathbb{I}_{i}∀ italic_j ∈ blackboard_I start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT, which are associated with Φi(Xi)superscriptΦ𝑖subscript𝑋𝑖\Phi^{i}(X_{i})roman_Φ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT ( italic_X start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ). The primary goal of the proposed parallel subsystem modeling approach is to establish a linear subsystem model for each subsystem of the underlying nonlinear process, in the form of (3a). Each subsystem model characterizes the dynamic behaviors of its corresponding subsystem. To create a subsystem model for each subsystem i𝑖iitalic_i, i𝕀𝑖𝕀i\in\mathbb{I}italic_i ∈ blackboard_I, we need to identify matrices Aiisubscript𝐴𝑖𝑖A_{ii}italic_A start_POSTSUBSCRIPT italic_i italic_i end_POSTSUBSCRIPT, Aijsubscript𝐴𝑖𝑗A_{ij}italic_A start_POSTSUBSCRIPT italic_i italic_j end_POSTSUBSCRIPT for j𝕀i𝑗subscript𝕀𝑖j\in{\mathbb{I}}_{i}italic_j ∈ blackboard_I start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT, and Bisubscript𝐵𝑖B_{i}italic_B start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT, which constitute the first-row blocks of the Koopman operator 𝒦~isubscript~𝒦𝑖\mathcal{{\tilde{K}}}_{i}over~ start_ARG caligraphic_K end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT. The remaining blocks of 𝒦~isubscript~𝒦𝑖\mathcal{{\tilde{K}}}_{i}over~ start_ARG caligraphic_K end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT are not required for completing the subsystem model in the form of (3a) for subsystem i𝑖iitalic_i. Therefore, only the blocks in the first row of the matrix on the right-hand side of (13) need to be identified, while all the remaining blocks are negligible. Consequently, for each subsystem i𝑖iitalic_i, i𝕀𝑖𝕀i\in\mathbb{I}italic_i ∈ blackboard_I, the least-squares problem in (12) is reformulated by retaining only the terms directly related to ϕisubscriptitalic-ϕ𝑖\phi_{i}italic_ϕ start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT as follows:

Aii,,Aij,,Bisubscript𝐴𝑖𝑖subscript𝐴𝑖𝑗subscript𝐵𝑖\displaystyle A_{ii},\ldots,A_{ij},\ldots,B_{i}italic_A start_POSTSUBSCRIPT italic_i italic_i end_POSTSUBSCRIPT , … , italic_A start_POSTSUBSCRIPT italic_i italic_j end_POSTSUBSCRIPT , … , italic_B start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT (14)
=argminAii,,Aij,,Bik=1Nϕi(xi(k+1))Aiiϕi(xi(k))j𝕀iAijϕj(xj(k))Biδi(ui(k))2absentsubscriptsubscript𝐴𝑖𝑖subscript𝐴𝑖𝑗subscript𝐵𝑖superscriptsubscript𝑘1𝑁superscriptnormsuperscriptitalic-ϕ𝑖subscript𝑥𝑖𝑘1subscript𝐴𝑖𝑖superscriptitalic-ϕ𝑖subscript𝑥𝑖𝑘subscript𝑗subscript𝕀𝑖subscript𝐴𝑖𝑗superscriptitalic-ϕ𝑗subscript𝑥𝑗𝑘subscript𝐵𝑖superscript𝛿𝑖subscript𝑢𝑖𝑘2\displaystyle=\arg\mathop{\min}\limits_{A_{ii},\ldots,A_{ij},\ldots,B_{i}}\sum% _{k=1}^{N}\bigg{\|}\phi^{i}\left(x_{i}(k+1)\right)-A_{ii}\phi^{i}\left(x_{i}(k% )\right)-\sum_{j\in\mathbb{I}_{i}}A_{ij}\phi^{j}\left(x_{j}(k)\right)-B_{i}% \delta^{i}\left(u_{i}(k)\right)\bigg{\|}^{2}= roman_arg roman_min start_POSTSUBSCRIPT italic_A start_POSTSUBSCRIPT italic_i italic_i end_POSTSUBSCRIPT , … , italic_A start_POSTSUBSCRIPT italic_i italic_j end_POSTSUBSCRIPT , … , italic_B start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT ∑ start_POSTSUBSCRIPT italic_k = 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_N end_POSTSUPERSCRIPT ∥ italic_ϕ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT ( italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k + 1 ) ) - italic_A start_POSTSUBSCRIPT italic_i italic_i end_POSTSUBSCRIPT italic_ϕ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT ( italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) ) - ∑ start_POSTSUBSCRIPT italic_j ∈ blackboard_I start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT italic_A start_POSTSUBSCRIPT italic_i italic_j end_POSTSUBSCRIPT italic_ϕ start_POSTSUPERSCRIPT italic_j end_POSTSUPERSCRIPT ( italic_x start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT ( italic_k ) ) - italic_B start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT italic_δ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT ( italic_u start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) ) ∥ start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT
=argminAii,,Aij,,Bik=1Nzi(k+1)Aiiϕizi(k)j𝕀iAijzj(k)Biu~i(k)2absentsubscriptsubscript𝐴𝑖𝑖subscript𝐴𝑖𝑗subscript𝐵𝑖superscriptsubscript𝑘1𝑁superscriptnormsubscript𝑧𝑖𝑘1subscript𝐴𝑖𝑖superscriptitalic-ϕ𝑖subscript𝑧𝑖𝑘subscript𝑗subscript𝕀𝑖subscript𝐴𝑖𝑗subscript𝑧𝑗𝑘subscript𝐵𝑖subscript~𝑢𝑖𝑘2\displaystyle=\arg\mathop{\min}\limits_{A_{ii},\ldots,A_{ij},\ldots,B_{i}}\sum% _{k=1}^{N}\bigg{\|}z_{i}(k+1)-A_{ii}\phi^{i}z_{i}(k)-\sum_{j\in\mathbb{I}_{i}}% A_{ij}z_{j}(k)-B_{i}{\tilde{u}}_{i}(k)\bigg{\|}^{2}= roman_arg roman_min start_POSTSUBSCRIPT italic_A start_POSTSUBSCRIPT italic_i italic_i end_POSTSUBSCRIPT , … , italic_A start_POSTSUBSCRIPT italic_i italic_j end_POSTSUBSCRIPT , … , italic_B start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT ∑ start_POSTSUBSCRIPT italic_k = 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_N end_POSTSUPERSCRIPT ∥ italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k + 1 ) - italic_A start_POSTSUBSCRIPT italic_i italic_i end_POSTSUBSCRIPT italic_ϕ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) - ∑ start_POSTSUBSCRIPT italic_j ∈ blackboard_I start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT italic_A start_POSTSUBSCRIPT italic_i italic_j end_POSTSUBSCRIPT italic_z start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT ( italic_k ) - italic_B start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT over~ start_ARG italic_u end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) ∥ start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT

As an analog of (14), to identify the output measurement matrix Cisubscript𝐶𝑖C_{i}italic_C start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT in (3b) and the projection matrix Disubscript𝐷𝑖D_{i}italic_D start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT in (3c) for subsystem i𝑖iitalic_i, the following least-squares problems are formulated:

Ci=subscript𝐶𝑖absent\displaystyle C_{i}=italic_C start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT = argminCik=1Nyi(k)Ciϕi(xi(k))2subscriptsubscript𝐶𝑖superscriptsubscript𝑘1𝑁superscriptnormsubscript𝑦𝑖𝑘subscript𝐶𝑖superscriptitalic-ϕ𝑖subscript𝑥𝑖𝑘2\displaystyle\arg\mathop{\min}\limits_{C_{i}}\sum_{k=1}^{N}\left\|y_{i}(k)-C_{% i}\phi^{i}(x_{i}(k))\right\|^{2}roman_arg roman_min start_POSTSUBSCRIPT italic_C start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT ∑ start_POSTSUBSCRIPT italic_k = 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_N end_POSTSUPERSCRIPT ∥ italic_y start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) - italic_C start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT italic_ϕ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT ( italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) ) ∥ start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT (15a)
Di=subscript𝐷𝑖absent\displaystyle D_{i}=italic_D start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT = argminDik=1Nxi(k)Diϕi(xi(k))2subscriptsubscript𝐷𝑖superscriptsubscript𝑘1𝑁superscriptnormsubscript𝑥𝑖𝑘subscript𝐷𝑖superscriptitalic-ϕ𝑖subscript𝑥𝑖𝑘2\displaystyle\arg\mathop{\min}\limits_{D_{i}}\sum_{k=1}^{N}\left\|x_{i}(k)-D_{% i}\phi^{i}(x_{i}(k))\right\|^{2}roman_arg roman_min start_POSTSUBSCRIPT italic_D start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT ∑ start_POSTSUBSCRIPT italic_k = 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_N end_POSTSUPERSCRIPT ∥ italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) - italic_D start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT italic_ϕ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT ( italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) ) ∥ start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT (15b)

Since the formulated optimization problems in (14) and (EDMD-based parallel Koopman subsystem identification method) are linear least-squares problems, they can be solved straightforwardly. Next, we present an implementation algorithm describing the key steps that need to be followed to conduct paralleled subsystem model identification as follows:

Algorithm 1.

Implementation algorithm for EDMD-based paralleled subsystem model identification

  1. 1.

    Collect N+1𝑁1N+1italic_N + 1 (Nnimuch-greater-than𝑁subscript𝑛𝑖N\gg n_{i}italic_N ≫ italic_n start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT) discrete samples for the process states x𝑥xitalic_x and known inputs u𝑢uitalic_u;

  2. 2.

    For each subsystem i𝑖iitalic_i, i𝕀𝑖𝕀i\in\mathbb{I}italic_i ∈ blackboard_I, do the following for data scaling:

    1. 2.1.

      Scale each subsystem state variable xi,lsubscript𝑥𝑖𝑙x_{i,l}italic_x start_POSTSUBSCRIPT italic_i , italic_l end_POSTSUBSCRIPT (denoting the l𝑙litalic_lth element of subsystem state vector xisubscript𝑥𝑖x_{i}italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT), l=1,,nxi𝑙1subscript𝑛subscript𝑥𝑖l=1,\ldots,n_{x_{i}}italic_l = 1 , … , italic_n start_POSTSUBSCRIPT italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT following xi,l(k)xi,l(k)xi,lminxi,lmaxxi,lminsubscript𝑥𝑖𝑙𝑘subscript𝑥𝑖𝑙𝑘subscriptsuperscript𝑥𝑖𝑙subscriptsuperscript𝑥𝑖𝑙subscriptsuperscript𝑥𝑖𝑙x_{i,l}(k)\leftarrow\frac{x_{i,l}(k)-x^{\min}_{i,l}}{x^{\max}_{i,l}-x^{\min}_{% i,l}}italic_x start_POSTSUBSCRIPT italic_i , italic_l end_POSTSUBSCRIPT ( italic_k ) ← divide start_ARG italic_x start_POSTSUBSCRIPT italic_i , italic_l end_POSTSUBSCRIPT ( italic_k ) - italic_x start_POSTSUPERSCRIPT roman_min end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_i , italic_l end_POSTSUBSCRIPT end_ARG start_ARG italic_x start_POSTSUPERSCRIPT roman_max end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_i , italic_l end_POSTSUBSCRIPT - italic_x start_POSTSUPERSCRIPT roman_min end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_i , italic_l end_POSTSUBSCRIPT end_ARG where xi,lmax:=max{xi,l(k):k*}assignsubscriptsuperscript𝑥𝑖𝑙:subscript𝑥𝑖𝑙𝑘𝑘superscriptx^{\max}_{i,l}:=\max\left\{x_{i,l}(k):k\in\mathbb{Z}^{*}\right\}italic_x start_POSTSUPERSCRIPT roman_max end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_i , italic_l end_POSTSUBSCRIPT := roman_max { italic_x start_POSTSUBSCRIPT italic_i , italic_l end_POSTSUBSCRIPT ( italic_k ) : italic_k ∈ blackboard_Z start_POSTSUPERSCRIPT * end_POSTSUPERSCRIPT }; xi,lmin:=min{xi,l(k):k*}assignsubscriptsuperscript𝑥𝑖𝑙:subscript𝑥𝑖𝑙𝑘𝑘superscriptx^{\min}_{i,l}:=\min\left\{x_{i,l}(k):k\in\mathbb{Z}^{*}\right\}italic_x start_POSTSUPERSCRIPT roman_min end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_i , italic_l end_POSTSUBSCRIPT := roman_min { italic_x start_POSTSUBSCRIPT italic_i , italic_l end_POSTSUBSCRIPT ( italic_k ) : italic_k ∈ blackboard_Z start_POSTSUPERSCRIPT * end_POSTSUPERSCRIPT }.

    2. 2.2.

      Scale each subsystem input variable ui,lsubscript𝑢𝑖𝑙u_{i,l}italic_u start_POSTSUBSCRIPT italic_i , italic_l end_POSTSUBSCRIPT (the l𝑙litalic_lth element of subsystem input vector uisubscript𝑢𝑖u_{i}italic_u start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT), and scale each sensor measurement yi,lsubscript𝑦𝑖𝑙y_{i,l}italic_y start_POSTSUBSCRIPT italic_i , italic_l end_POSTSUBSCRIPT (the l𝑙litalic_lth element of vector yisubscript𝑦𝑖y_{i}italic_y start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT) in a way similar to Step 2.1.

    3. 2.3.

      Create subsystem snapshot matrices 𝒙𝒊subscript𝒙𝒊{\boldsymbol{x_{i}}}bold_italic_x start_POSTSUBSCRIPT bold_italic_i end_POSTSUBSCRIPT, 𝒙¯𝒊subscriptbold-¯𝒙𝒊{\boldsymbol{{\bar{x}}_{i}}}overbold_¯ start_ARG bold_italic_x end_ARG start_POSTSUBSCRIPT bold_italic_i end_POSTSUBSCRIPT, 𝒖𝒊subscript𝒖𝒊{\boldsymbol{u_{i}}}bold_italic_u start_POSTSUBSCRIPT bold_italic_i end_POSTSUBSCRIPT and 𝒚𝒊subscript𝒚𝒊{\boldsymbol{y_{i}}}bold_italic_y start_POSTSUBSCRIPT bold_italic_i end_POSTSUBSCRIPT as: 𝒙𝒊=[xi(1),,xi(N)]subscript𝒙𝒊subscript𝑥𝑖1subscript𝑥𝑖𝑁{\boldsymbol{x_{i}}}=\left[x_{i}(1),\ldots,x_{i}(N)\right]bold_italic_x start_POSTSUBSCRIPT bold_italic_i end_POSTSUBSCRIPT = [ italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( 1 ) , … , italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_N ) ], 𝒙¯𝒊=[xi(2),,xi(N+1)]subscriptbold-¯𝒙𝒊subscript𝑥𝑖2subscript𝑥𝑖𝑁1{\boldsymbol{{\bar{x}}_{i}}}=\left[x_{i}(2),\ldots,x_{i}(N+1)\right]overbold_¯ start_ARG bold_italic_x end_ARG start_POSTSUBSCRIPT bold_italic_i end_POSTSUBSCRIPT = [ italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( 2 ) , … , italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_N + 1 ) ], 𝒖𝒊=[ui(1),,ui(N)]subscript𝒖𝒊subscript𝑢𝑖1subscript𝑢𝑖𝑁{\boldsymbol{u_{i}}}=\left[u_{i}(1),\ldots,u_{i}(N)\right]bold_italic_u start_POSTSUBSCRIPT bold_italic_i end_POSTSUBSCRIPT = [ italic_u start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( 1 ) , … , italic_u start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_N ) ], and 𝒚𝒊=[yi(1),,yi(N)]subscript𝒚𝒊subscript𝑦𝑖1subscript𝑦𝑖𝑁{\boldsymbol{y_{i}}}=\left[y_{i}(1),\ldots,y_{i}(N)\right]bold_italic_y start_POSTSUBSCRIPT bold_italic_i end_POSTSUBSCRIPT = [ italic_y start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( 1 ) , … , italic_y start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_N ) ].

    4. 2.4.

      Specify a set of state lifting functions ϕi(xi)=[ϕ1i(xi),,ϕnzii(xi)]𝑇superscriptitalic-ϕ𝑖subscript𝑥𝑖superscriptsubscriptsuperscriptitalic-ϕ𝑖1subscript𝑥𝑖subscriptsuperscriptitalic-ϕ𝑖subscript𝑛subscript𝑧𝑖subscript𝑥𝑖𝑇\phi^{i}(x_{i})=\big{[}\phi^{i}_{1}(x_{i}),\ldots,\phi^{i}_{n_{z_{i}}}(x_{i})% \big{]}^{\text{T}}italic_ϕ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT ( italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) = [ italic_ϕ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT ( italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) , … , italic_ϕ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_n start_POSTSUBSCRIPT italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT end_POSTSUBSCRIPT ( italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) ] start_POSTSUPERSCRIPT T end_POSTSUPERSCRIPT, with the first nxisubscript𝑛subscript𝑥𝑖n_{x_{i}}italic_n start_POSTSUBSCRIPT italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT functions defined as ϕli(xi)=xi,lsubscriptsuperscriptitalic-ϕ𝑖𝑙subscript𝑥𝑖subscript𝑥𝑖𝑙\phi^{i}_{l}(x_{i})=x_{i,l}italic_ϕ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_l end_POSTSUBSCRIPT ( italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) = italic_x start_POSTSUBSCRIPT italic_i , italic_l end_POSTSUBSCRIPT, l=1,,nxi𝑙1subscript𝑛subscript𝑥𝑖l=1,\ldots,n_{x_{i}}italic_l = 1 , … , italic_n start_POSTSUBSCRIPT italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT.

    5. 2.5.

      Specify a set of input lifting functions δi(ui)=[δ1i(ui),,δnu~ii(ui)]𝑇superscript𝛿𝑖subscript𝑢𝑖superscriptsubscriptsuperscript𝛿𝑖1subscript𝑢𝑖subscriptsuperscript𝛿𝑖subscript𝑛subscript~𝑢𝑖subscript𝑢𝑖𝑇\delta^{i}(u_{i})=\big{[}\delta^{i}_{1}(u_{i}),\ldots,\delta^{i}_{n_{{\tilde{u% }}_{i}}}(u_{i})\big{]}^{\text{T}}italic_δ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT ( italic_u start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) = [ italic_δ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT ( italic_u start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) , … , italic_δ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_n start_POSTSUBSCRIPT over~ start_ARG italic_u end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT end_POSTSUBSCRIPT ( italic_u start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) ] start_POSTSUPERSCRIPT T end_POSTSUPERSCRIPT.

  3. 3.

    For each subsystem i𝑖iitalic_i, i𝕀𝑖𝕀i\in\mathbb{I}italic_i ∈ blackboard_I, do the following in parallel for subsystem model identification:

    1. 3.1

      The identification agent for subsystem i𝑖iitalic_i receives ϕj(𝒙𝒊)superscriptitalic-ϕ𝑗subscript𝒙𝒊\phi^{j}\left(\boldsymbol{x_{i}}\right)italic_ϕ start_POSTSUPERSCRIPT italic_j end_POSTSUPERSCRIPT ( bold_italic_x start_POSTSUBSCRIPT bold_italic_i end_POSTSUBSCRIPT ) from the agents for each subsystem j𝑗jitalic_j, j𝕀i𝑗subscript𝕀𝑖j\in\mathbb{I}_{i}italic_j ∈ blackboard_I start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT.

    2. 3.2

      The identification agent for subsystem i𝑖iitalic_i computes subsystem matrices Aiisubscript𝐴𝑖𝑖A_{ii}italic_A start_POSTSUBSCRIPT italic_i italic_i end_POSTSUBSCRIPT, Aij(j𝕀iA_{ij}(\forall j\in\mathbb{I}_{i}italic_A start_POSTSUBSCRIPT italic_i italic_j end_POSTSUBSCRIPT ( ∀ italic_j ∈ blackboard_I start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT), and Bisubscript𝐵𝑖B_{i}italic_B start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT following: [Aii…,A_ij,…B_i]=ϕi(𝒙¯𝒊)ψi(𝓧𝒊)𝑇(ψi(𝓧𝒊)ψi(𝓧𝒊)𝑇)+delimited-[]subscript𝐴𝑖𝑖…,A_ij,…B_imissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionsuperscriptitalic-ϕ𝑖subscriptbold-¯𝒙𝒊superscript𝜓𝑖superscriptsubscript𝓧𝒊𝑇superscriptsuperscript𝜓𝑖subscript𝓧𝒊superscript𝜓𝑖superscriptsubscript𝓧𝒊𝑇\left[\begin{array}[]{c;{2pt/2pt}c;{2pt/2pt}c}A_{ii}&\ldots,A_{ij},\ldots&B_i% \end{array}\right]={\phi^{i}}\left({\boldsymbol{{\bar{x}}_{i}}}\right)\cdot{% \psi^{i}\left({\boldsymbol{\mathcal{X}_{i}}}\right)}^{\text{T}}\cdot\left(\psi% ^{i}\left({\boldsymbol{\mathcal{X}_{i}}}\right)\cdot{\psi^{i}\left({% \boldsymbol{\mathcal{X}_{i}}}\right)}^{\text{T}}\right)^{+}[ start_ARRAY start_ROW start_CELL italic_A start_POSTSUBSCRIPT italic_i italic_i end_POSTSUBSCRIPT end_CELL start_CELL …,A_ij,… end_CELL start_CELL B_i end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL end_ROW end_ARRAY ] = italic_ϕ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT ( overbold_¯ start_ARG bold_italic_x end_ARG start_POSTSUBSCRIPT bold_italic_i end_POSTSUBSCRIPT ) ⋅ italic_ψ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT ( bold_caligraphic_X start_POSTSUBSCRIPT bold_italic_i end_POSTSUBSCRIPT ) start_POSTSUPERSCRIPT T end_POSTSUPERSCRIPT ⋅ ( italic_ψ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT ( bold_caligraphic_X start_POSTSUBSCRIPT bold_italic_i end_POSTSUBSCRIPT ) ⋅ italic_ψ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT ( bold_caligraphic_X start_POSTSUBSCRIPT bold_italic_i end_POSTSUBSCRIPT ) start_POSTSUPERSCRIPT T end_POSTSUPERSCRIPT ) start_POSTSUPERSCRIPT + end_POSTSUPERSCRIPT.

    3. 3.3

      The identification agent for subsystem i𝑖iitalic_i computes output measurement matrix Cisubscript𝐶𝑖C_{i}italic_C start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT as: Ci=𝒚𝒊ϕi(𝒙𝒊)𝑇(ϕi(𝒙𝒊)ϕi(𝒙𝒊)𝑇)+subscript𝐶𝑖subscript𝒚𝒊superscriptitalic-ϕ𝑖superscriptsubscript𝒙𝒊𝑇superscriptsuperscriptitalic-ϕ𝑖subscript𝒙𝒊superscriptitalic-ϕ𝑖superscriptsubscript𝒙𝒊𝑇C_{i}={\boldsymbol{y_{i}}}\cdot\phi^{i}(\boldsymbol{x_{i}})^{\text{T}}\cdot% \left(\phi^{i}(\boldsymbol{x_{i}})\cdot\phi^{i}(\boldsymbol{x_{i}})^{\text{T}}% \right)^{+}italic_C start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT = bold_italic_y start_POSTSUBSCRIPT bold_italic_i end_POSTSUBSCRIPT ⋅ italic_ϕ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT ( bold_italic_x start_POSTSUBSCRIPT bold_italic_i end_POSTSUBSCRIPT ) start_POSTSUPERSCRIPT T end_POSTSUPERSCRIPT ⋅ ( italic_ϕ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT ( bold_italic_x start_POSTSUBSCRIPT bold_italic_i end_POSTSUBSCRIPT ) ⋅ italic_ϕ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT ( bold_italic_x start_POSTSUBSCRIPT bold_italic_i end_POSTSUBSCRIPT ) start_POSTSUPERSCRIPT T end_POSTSUPERSCRIPT ) start_POSTSUPERSCRIPT + end_POSTSUPERSCRIPT.

In Step 3.2 of Algorithm 1, ψi(𝓧𝒊)=[ϕi(𝒙𝒊)TΦ^i(X_i)^T𝛿^i(u_i)^T]Tsuperscript𝜓𝑖subscript𝓧𝒊superscriptdelimited-[]superscriptitalic-ϕ𝑖superscriptsubscript𝒙𝒊TΦ^i(X_i)^T𝛿^i(u_i)^Tmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionT\psi^{i}\left({\boldsymbol{\mathcal{X}_{i}}}\right)=\left[\begin{array}[]{c;{2% pt/2pt}c;{2pt/2pt}c}{\phi^{i}({\boldsymbol{x_{i}}})}^{\text{T}}&{\Phi^i({% \boldsymbol{X_i}})}^{\text{T}}&{\delta^i({\boldsymbol{u_i}})}^{\text{T}}\end{% array}\right]^{\text{T}}italic_ψ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT ( bold_caligraphic_X start_POSTSUBSCRIPT bold_italic_i end_POSTSUBSCRIPT ) = [ start_ARRAY start_ROW start_CELL italic_ϕ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT ( bold_italic_x start_POSTSUBSCRIPT bold_italic_i end_POSTSUBSCRIPT ) start_POSTSUPERSCRIPT T end_POSTSUPERSCRIPT end_CELL start_CELL Φ^i( bold_X_i )^ roman_T end_CELL start_CELL italic_δ ^i( bold_u_i )^ roman_T end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL end_ROW end_ARRAY ] start_POSTSUPERSCRIPT T end_POSTSUPERSCRIPT. In Step 2.4, for each subsystem i𝑖iitalic_i, i𝕀𝑖𝕀i\in\mathbb{I}italic_i ∈ blackboard_I, the first nxisubscript𝑛subscript𝑥𝑖n_{x_{i}}italic_n start_POSTSUBSCRIPT italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT state lifting functions are made identical to the states of the local subsystem, that is, ϕil(xi)=xi,lsubscriptsuperscriptitalic-ϕ𝑙𝑖subscript𝑥𝑖subscript𝑥𝑖𝑙\phi^{l}_{i}(x_{i})=x_{i,l}italic_ϕ start_POSTSUPERSCRIPT italic_l end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) = italic_x start_POSTSUBSCRIPT italic_i , italic_l end_POSTSUBSCRIPT, l=1,,nxi𝑙1subscript𝑛subscript𝑥𝑖l=1,\ldots,n_{x_{i}}italic_l = 1 , … , italic_n start_POSTSUBSCRIPT italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT. This ensures the feasibility of reconstructing the original process states based on the prediction of the lifted state variables. Consequently, the state-reconstruction matrix Disubscript𝐷𝑖D_{i}italic_D start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT is determined directly as Di=[𝐈nxi𝟎]subscript𝐷𝑖delimited-[]subscript𝐈subscript𝑛subscript𝑥𝑖𝟎missing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionD_{i}=\left[{\begin{array}[]{*{20}{c;{2pt/2pt}c}}{\bf{I}}_{n_{x_{i}}}&\bf{0}% \end{array}}\right]italic_D start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT = [ start_ARRAY start_ROW start_CELL bold_I start_POSTSUBSCRIPT italic_n start_POSTSUBSCRIPT italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT end_POSTSUBSCRIPT end_CELL start_CELL 0 end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL end_ROW end_ARRAY ], without needing to solve (15b). Also, in many case scenarios, the measurement equations for the subsystems (1b) are linear, in the form of yi=Hixisubscript𝑦𝑖subscript𝐻𝑖subscript𝑥𝑖y_{i}=H_{i}x_{i}italic_y start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT = italic_H start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT (e.g., the sensor measurements are the temperatures of certain physical units belonging to one subsystem, which are also state variables of the subsystem). In such cases, matrix C𝐶Citalic_C can be determined in a straightforward manner as Ci=HiDi=[Hi𝟎_n_y_i×(n_z_i-n_x_i)]subscript𝐶𝑖subscript𝐻𝑖subscript𝐷𝑖delimited-[]subscript𝐻𝑖𝟎_n_y_i×(n_z_i-n_x_i)missing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionC_{i}=H_{i}D_{i}=\left[{\begin{array}[]{*{20}{c;{2pt/2pt}c}}H_{i}&\boldsymbol{% 0}_{n_{y_i}\times\left(n_{z_i}-n_{x_i}\right)}\end{array}}\right]italic_C start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT = italic_H start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT italic_D start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT = [ start_ARRAY start_ROW start_CELL italic_H start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_CELL start_CELL bold_0 _n_y_i× ( n_z_i-n_x_i ) end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL end_ROW end_ARRAY ] without resorting to Step 3.3.

Remark 2.

The optimal selection of the state basis functions ϕisuperscriptitalic-ϕ𝑖\phi^{i}italic_ϕ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT and the input basis functions δisuperscript𝛿𝑖\delta^{i}italic_δ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT for each subsystem is still an open question for EDMD-based algorithms [28, 32]. From a practical point of view, attempts can be made to find an appropriate set of basis functions that are “rich enough” to account for the leading Koopman eigenfunctions, such that the Koopman-based subsystem models can be made sufficiently accurate. more detailed discussions pertaining to this point can be found in Williams et al. [28]

Remark 3.

For many processes, the key variables have vastly different magnitudes (e.g., the mass fractions and temperatures which are usually considered process state variables); this may significantly affect the accuracy of the identified models. By using Steps 2.2 and 2.3, the process state, input, and sensor measurement variables are scaled to be of similar magnitudes, such that almost equal importance is given to the variables and potential numerical issues can be avoided.

Partition-based distributed moving horizon estimation

In this section, we propose a linear partition-based distributed moving horizon estimation (MHE) method based on the linear subsystem models that are identified by applying the EDMD-based paralleled subsystem model identification approach proposed in the “EDMD-based parallel Koopman subsystem identification method” subsection. An illustrative diagram that describes the connection between the subsystem modeling and the local estimator design of the proposed framework is shown in Figure 1. In addition, a schematic of the distributed estimation scheme is given in Figure 2. The distributed state estimation scheme consists of m𝑚mitalic_m local MHE-based estimators interconnected with each other. Each local estimator accounts for the estimation of the states zisubscript𝑧𝑖z_{i}italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT (i𝕀𝑖𝕀i\in\mathbb{I}italic_i ∈ blackboard_I) of the corresponding subsystem in the higher-dimensional space. A communication network is deployed to account for the real-time information exchange among the different subsystems. At each sampling instant, each MHE-based local estimator receives sensor measurements from its corresponding subsystem and communicates with the interacting subsystems to exchange subsystem state estimates and sensor measurements through the communication network.

Subsystem models for the design of local estimators

The identified Koopman subsystem models are used to develop the local MHE-based estimators of the distributed scheme in Figure 2. To account for unknown system disturbances and sensor measurement noise, a stochastic version of the identified model for subsystem i𝑖iitalic_i is given as follows:

zi(k+1)=subscript𝑧𝑖𝑘1absent\displaystyle z_{i}(k+1)=italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k + 1 ) = Aiizi(k)+j𝕀iAijzj(k)+Biu~i(k)+wi(k)subscript𝐴𝑖𝑖subscript𝑧𝑖𝑘subscript𝑗subscript𝕀𝑖subscript𝐴𝑖𝑗subscript𝑧𝑗𝑘subscript𝐵𝑖subscript~𝑢𝑖𝑘subscript𝑤𝑖𝑘\displaystyle~{}A_{ii}z_{i}(k)+\sum_{j\in\mathbb{I}_{i}}A_{ij}z_{j}(k)+B_{i}{% \tilde{u}}_{i}(k)+w_{i}(k)italic_A start_POSTSUBSCRIPT italic_i italic_i end_POSTSUBSCRIPT italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) + ∑ start_POSTSUBSCRIPT italic_j ∈ blackboard_I start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT italic_A start_POSTSUBSCRIPT italic_i italic_j end_POSTSUBSCRIPT italic_z start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT ( italic_k ) + italic_B start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT over~ start_ARG italic_u end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) + italic_w start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) (16a)
yi(k)=subscript𝑦𝑖𝑘absent\displaystyle y_{i}(k)=italic_y start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) = Cizi(k)+vi(k)subscript𝐶𝑖subscript𝑧𝑖𝑘subscript𝑣𝑖𝑘\displaystyle~{}C_{i}z_{i}(k)+v_{i}(k)italic_C start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) + italic_v start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) (16b)
x^i(k)=subscript^𝑥𝑖𝑘absent\displaystyle{\hat{x}}_{i}(k)=over^ start_ARG italic_x end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) = Dizi(k)subscript𝐷𝑖subscript𝑧𝑖𝑘\displaystyle~{}D_{i}z_{i}(k)italic_D start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) (16c)

where wi𝕎isubscript𝑤𝑖subscript𝕎𝑖w_{i}\in\mathbb{W}_{i}italic_w start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ∈ blackboard_W start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT represents the disturbances and vi𝕍isubscript𝑣𝑖subscript𝕍𝑖v_{i}\in\mathbb{V}_{i}italic_v start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ∈ blackboard_V start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT represents the measurement noise for subsystem i𝑖iitalic_i, with 𝕎inzisubscript𝕎𝑖superscriptsubscript𝑛subscript𝑧𝑖\mathbb{W}_{i}\subset\mathbb{R}^{n_{z_{i}}}blackboard_W start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ⊂ blackboard_R start_POSTSUPERSCRIPT italic_n start_POSTSUBSCRIPT italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT end_POSTSUPERSCRIPT and 𝕍ixyisubscript𝕍𝑖superscriptsubscript𝑥subscript𝑦𝑖\mathbb{V}_{i}\subset\mathbb{R}^{x_{y_{i}}}blackboard_V start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ⊂ blackboard_R start_POSTSUPERSCRIPT italic_x start_POSTSUBSCRIPT italic_y start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT end_POSTSUPERSCRIPT being two compact sets.

Centralized MHE

Before introducing the proposed distributed MHE design, we briefly review a representative centralized MHE design reported in Rao at al. [47] Consider a general linear system model as follows:

z(k+1)𝑧𝑘1\displaystyle z(k+1)italic_z ( italic_k + 1 ) =Az(k)+Bu~(k)+w(k)absent𝐴𝑧𝑘𝐵~𝑢𝑘𝑤𝑘\displaystyle=Az(k)+B{\tilde{u}}(k)+w(k)= italic_A italic_z ( italic_k ) + italic_B over~ start_ARG italic_u end_ARG ( italic_k ) + italic_w ( italic_k ) (17)
y(k)𝑦𝑘\displaystyle y(k)italic_y ( italic_k ) =Cz(k)+v(k)absent𝐶𝑧𝑘𝑣𝑘\displaystyle=Cz(k)+v(k)= italic_C italic_z ( italic_k ) + italic_v ( italic_k )

where z𝑧zitalic_z is the system state; u~~𝑢\tilde{u}over~ start_ARG italic_u end_ARG is the known input vector; w𝑤witalic_w represents unknown disturbances; y𝑦yitalic_y is the output measurements; v𝑣vitalic_v is the measurement noise. Note that (17) maybe viewed as an aggregation of the identified linear subsystem models in (Subsystem models for the design of local estimators). For (17), the centralized MHE design is presented as follows [47]:

minz^(kN|k),{w^(d|k)}d=kNk1Θ(z^(kN|k),w^(kN|k),,w^(k1|k))subscript^𝑧𝑘conditional𝑁𝑘subscriptsuperscript^𝑤conditional𝑑𝑘𝑘1𝑑𝑘𝑁Θ^𝑧𝑘conditional𝑁𝑘^𝑤𝑘conditional𝑁𝑘^𝑤𝑘conditional1𝑘\displaystyle\mathop{\min}\limits_{{\hat{z}({k-N|k})},\left\{{\hat{w}({d|k})}% \right\}^{k-1}_{d=k-N}}\Theta\big{(}\hat{z}({k-N|k}),\hat{w}(k-N|k),\ldots,% \hat{w}(k-1|k)\big{)}roman_min start_POSTSUBSCRIPT over^ start_ARG italic_z end_ARG ( italic_k - italic_N | italic_k ) , { over^ start_ARG italic_w end_ARG ( italic_d | italic_k ) } start_POSTSUPERSCRIPT italic_k - 1 end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_d = italic_k - italic_N end_POSTSUBSCRIPT end_POSTSUBSCRIPT roman_Θ ( over^ start_ARG italic_z end_ARG ( italic_k - italic_N | italic_k ) , over^ start_ARG italic_w end_ARG ( italic_k - italic_N | italic_k ) , … , over^ start_ARG italic_w end_ARG ( italic_k - 1 | italic_k ) ) (18a)
s.t. z^(d+1|k)=Az^(d|k)+Bu~(d)+w^(d|k),d=kN,,k1formulae-sequences.t. ^𝑧𝑑conditional1𝑘𝐴^𝑧conditional𝑑𝑘𝐵~𝑢𝑑^𝑤conditional𝑑𝑘𝑑𝑘𝑁𝑘1\displaystyle\textmd{\bf{ s.t.~{}}}\quad{\hat{z}}({d+1|k})=A{\hat{z}}({d|k})+B% {\tilde{u}}(d)+{\hat{w}}({d|k}),~{}~{}d=k-N,\ldots,k-1s.t. over^ start_ARG italic_z end_ARG ( italic_d + 1 | italic_k ) = italic_A over^ start_ARG italic_z end_ARG ( italic_d | italic_k ) + italic_B over~ start_ARG italic_u end_ARG ( italic_d ) + over^ start_ARG italic_w end_ARG ( italic_d | italic_k ) , italic_d = italic_k - italic_N , … , italic_k - 1 (18b)
y(d)=Cz^(d|k)+v^(d|k),d=kN,,kformulae-sequence𝑦𝑑𝐶^𝑧conditional𝑑𝑘^𝑣conditional𝑑𝑘𝑑𝑘𝑁𝑘\displaystyle~{}~{}\quad\qquad\qquad y(d)=C{\hat{z}}({d|k})+{\hat{v}}({d|k}),~% {}~{}d=k-N,~{}\ldots,~{}kitalic_y ( italic_d ) = italic_C over^ start_ARG italic_z end_ARG ( italic_d | italic_k ) + over^ start_ARG italic_v end_ARG ( italic_d | italic_k ) , italic_d = italic_k - italic_N , … , italic_k (18c)

where the cost function is:

Θ=d=kNk1w^(d|k)Q12+d=kNkv^(d|k)R12+z^(kN|k)z¯(kN)PkN12\Theta=\sum\limits_{d=k-N}^{k-1}\left\|{\hat{w}}(d|k)\right\|^{2}_{Q^{-1}}+% \sum\limits_{d=k-N}^{k}\left\|{\hat{v}}(d|k)\right\|^{2}_{R^{-1}}+\left\|{\hat% {z}}(k-N|k)-{\bar{z}}(k-N)\right\|^{2}_{P_{k-N}^{-1}}roman_Θ = ∑ start_POSTSUBSCRIPT italic_d = italic_k - italic_N end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_k - 1 end_POSTSUPERSCRIPT ∥ over^ start_ARG italic_w end_ARG ( italic_d | italic_k ) ∥ start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_Q start_POSTSUPERSCRIPT - 1 end_POSTSUPERSCRIPT end_POSTSUBSCRIPT + ∑ start_POSTSUBSCRIPT italic_d = italic_k - italic_N end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_k end_POSTSUPERSCRIPT ∥ over^ start_ARG italic_v end_ARG ( italic_d | italic_k ) ∥ start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_R start_POSTSUPERSCRIPT - 1 end_POSTSUPERSCRIPT end_POSTSUBSCRIPT + ∥ over^ start_ARG italic_z end_ARG ( italic_k - italic_N | italic_k ) - over¯ start_ARG italic_z end_ARG ( italic_k - italic_N ) ∥ start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_P start_POSTSUBSCRIPT italic_k - italic_N end_POSTSUBSCRIPT start_POSTSUPERSCRIPT - 1 end_POSTSUPERSCRIPT end_POSTSUBSCRIPT (19a)
The weighting matrix Pksubscript𝑃𝑘P_{k}italic_P start_POSTSUBSCRIPT italic_k end_POSTSUBSCRIPT is updated using the covariance matrix update formula for the Kalman filter as:
Pk=Q+APk1ATAPk1CT(R+CPk1CT)1CPk1ATsubscript𝑃𝑘𝑄𝐴subscript𝑃𝑘1superscript𝐴T𝐴subscript𝑃𝑘1superscript𝐶Tsuperscript𝑅𝐶subscript𝑃𝑘1superscript𝐶T1𝐶subscript𝑃𝑘1superscript𝐴TP_{k}=Q+AP_{k-1}A^{\mathrm{T}}-AP_{k-1}C^{\mathrm{T}}(R+CP_{k-1}C^{\mathrm{T}}% )^{-1}CP_{k-1}A^{\mathrm{T}}italic_P start_POSTSUBSCRIPT italic_k end_POSTSUBSCRIPT = italic_Q + italic_A italic_P start_POSTSUBSCRIPT italic_k - 1 end_POSTSUBSCRIPT italic_A start_POSTSUPERSCRIPT roman_T end_POSTSUPERSCRIPT - italic_A italic_P start_POSTSUBSCRIPT italic_k - 1 end_POSTSUBSCRIPT italic_C start_POSTSUPERSCRIPT roman_T end_POSTSUPERSCRIPT ( italic_R + italic_C italic_P start_POSTSUBSCRIPT italic_k - 1 end_POSTSUBSCRIPT italic_C start_POSTSUPERSCRIPT roman_T end_POSTSUPERSCRIPT ) start_POSTSUPERSCRIPT - 1 end_POSTSUPERSCRIPT italic_C italic_P start_POSTSUBSCRIPT italic_k - 1 end_POSTSUBSCRIPT italic_A start_POSTSUPERSCRIPT roman_T end_POSTSUPERSCRIPT (19b)

In (Centralized MHE), z^^𝑧\hat{z}over^ start_ARG italic_z end_ARG is an estimate of z𝑧zitalic_z; w^^𝑤\hat{w}over^ start_ARG italic_w end_ARG and v^^𝑣\hat{v}over^ start_ARG italic_v end_ARG are the estimates of w𝑤witalic_w and v𝑣vitalic_v, respectively; Q𝑄Qitalic_Q, R𝑅Ritalic_R, and Pksubscript𝑃𝑘P_{k}italic_P start_POSTSUBSCRIPT italic_k end_POSTSUBSCRIPT are weighting matrices; N𝑁Nitalic_N is the size of the estimation window; z¯(kN)¯𝑧𝑘𝑁\bar{z}(k-N)over¯ start_ARG italic_z end_ARG ( italic_k - italic_N ) is a one-step ahead open-loop prediction of z(kN)𝑧𝑘𝑁z(k-N)italic_z ( italic_k - italic_N ) computed following z¯(kN)=Az^(kN1|k1)+Bu~(kN1)+w^(kN1|k1)¯𝑧𝑘𝑁𝐴^𝑧𝑘𝑁conditional1𝑘1𝐵~𝑢𝑘𝑁1^𝑤𝑘𝑁conditional1𝑘1\bar{z}(k-N)=A\hat{z}(k-N-1|k-1)+B\tilde{u}(k-N-1)+\hat{w}(k-N-1|k-1)over¯ start_ARG italic_z end_ARG ( italic_k - italic_N ) = italic_A over^ start_ARG italic_z end_ARG ( italic_k - italic_N - 1 | italic_k - 1 ) + italic_B over~ start_ARG italic_u end_ARG ( italic_k - italic_N - 1 ) + over^ start_ARG italic_w end_ARG ( italic_k - italic_N - 1 | italic_k - 1 ). At each sampling instant kN𝑘𝑁k\geq Nitalic_k ≥ italic_N, with the optimal estimates z^(kN|k)^𝑧𝑘conditional𝑁𝑘\hat{z}({k-N|k})over^ start_ARG italic_z end_ARG ( italic_k - italic_N | italic_k ) and {w^(d|k)}d=kNk1subscriptsuperscript^𝑤conditional𝑑𝑘𝑘1𝑑𝑘𝑁\left\{{\hat{w}({d|k})}\right\}^{k-1}_{d=k-N}{ over^ start_ARG italic_w end_ARG ( italic_d | italic_k ) } start_POSTSUPERSCRIPT italic_k - 1 end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_d = italic_k - italic_N end_POSTSUBSCRIPT as well as model constraint (18b), the centralized MHE provides a sequence of the estimates z^(kN|k),,z^(k|k)^𝑧𝑘conditional𝑁𝑘^𝑧conditional𝑘𝑘\hat{z}({k-N|k}),\ldots,\hat{z}({k|k})over^ start_ARG italic_z end_ARG ( italic_k - italic_N | italic_k ) , … , over^ start_ARG italic_z end_ARG ( italic_k | italic_k ) for the actual system states within the current estimation window. The last element of the sequence of the estimates (i.e., z^(k|k)^𝑧conditional𝑘𝑘\hat{z}({k|k})over^ start_ARG italic_z end_ARG ( italic_k | italic_k )) is used as the estimate of the state at the current instant.

Partition-based cost function for local MHE-based estimators

In this subsection, we present the cost function for each local estimator of the proposed partition-based distributed MHE method, which is obtained based on partitioning the cost function for centralized MHE in (Centralized MHE).

First, we define the individual cost function for each MHE-based estimator for subsystem i𝑖iitalic_i, i𝕀𝑖𝕀i\in\mathbb{I}italic_i ∈ blackboard_I, of the distributed MHE by partitioning the cost function in (19a). Considering the state variables of each subsystem, ΘΘ\Thetaroman_Θ in (19a) can be decomposed into m𝑚mitalic_m partitioned functions ΘisubscriptΘ𝑖\Theta_{i}roman_Θ start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT, i=1,,m𝑖1𝑚i=1,\ldots,mitalic_i = 1 , … , italic_m, such that the partitioned function associated with subsystem i𝑖iitalic_i is with the following form:

Θi=d=kNk1w^i(d|k)Qi12+d=kNkv^i(d|k)Ri12+z^i(kN|k)z¯i(kN)Pi,kN12\Theta_{i}=\sum\limits_{d=k-N}^{k-1}\left\|{\hat{w}_{i}}(d|k)\right\|^{2}_{Q_{% i}^{-1}}+\sum\limits_{d=k-N}^{k}\left\|{\hat{v}_{i}}(d|k)\right\|^{2}_{R_{i}^{% -1}}+\left\|{\hat{z}}_{i}({k-N|k})-{\bar{z}}_{i}(k-N)\right\|^{2}_{P_{i,k-N}^{% -1}}roman_Θ start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT = ∑ start_POSTSUBSCRIPT italic_d = italic_k - italic_N end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_k - 1 end_POSTSUPERSCRIPT ∥ over^ start_ARG italic_w end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_d | italic_k ) ∥ start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_Q start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT start_POSTSUPERSCRIPT - 1 end_POSTSUPERSCRIPT end_POSTSUBSCRIPT + ∑ start_POSTSUBSCRIPT italic_d = italic_k - italic_N end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_k end_POSTSUPERSCRIPT ∥ over^ start_ARG italic_v end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_d | italic_k ) ∥ start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_R start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT start_POSTSUPERSCRIPT - 1 end_POSTSUPERSCRIPT end_POSTSUBSCRIPT + ∥ over^ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k - italic_N | italic_k ) - over¯ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k - italic_N ) ∥ start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_P start_POSTSUBSCRIPT italic_i , italic_k - italic_N end_POSTSUBSCRIPT start_POSTSUPERSCRIPT - 1 end_POSTSUPERSCRIPT end_POSTSUBSCRIPT (20a)
where Pi,ksubscript𝑃𝑖𝑘P_{i,k}italic_P start_POSTSUBSCRIPT italic_i , italic_k end_POSTSUBSCRIPT is updated using the covariance matrix update formula of distributed Kalman filter proposed in Li et al. [48]
Li,ksubscript𝐿𝑖𝑘\displaystyle L_{i,k}italic_L start_POSTSUBSCRIPT italic_i , italic_k end_POSTSUBSCRIPT =(CA[:,i]Pi,k1AiiT+C[:,i]Qi)T(CA[:,i]Pi,k1A[:,i]TCT+C[:,i]QiC[:,i]T+R)1absentsuperscript𝐶subscript𝐴:𝑖subscript𝑃𝑖𝑘1superscriptsubscript𝐴𝑖𝑖Tsubscript𝐶:𝑖subscript𝑄𝑖Tsuperscript𝐶subscript𝐴:𝑖subscript𝑃𝑖𝑘1superscriptsubscript𝐴:𝑖Tsuperscript𝐶Tsubscript𝐶:𝑖subscript𝑄𝑖superscriptsubscript𝐶:𝑖T𝑅1\displaystyle=(CA_{[:,i]}P_{i,k-1}A_{ii}^{\mathrm{T}}+C_{[:,i]}Q_{i})^{\mathrm% {T}}(CA_{[:,i]}P_{i,k-1}A_{[:,i]}^{\mathrm{T}}C^{\mathrm{T}}+C_{[:,i]}Q_{i}C_{% [:,i]}^{\mathrm{T}}+R)^{-1}= ( italic_C italic_A start_POSTSUBSCRIPT [ : , italic_i ] end_POSTSUBSCRIPT italic_P start_POSTSUBSCRIPT italic_i , italic_k - 1 end_POSTSUBSCRIPT italic_A start_POSTSUBSCRIPT italic_i italic_i end_POSTSUBSCRIPT start_POSTSUPERSCRIPT roman_T end_POSTSUPERSCRIPT + italic_C start_POSTSUBSCRIPT [ : , italic_i ] end_POSTSUBSCRIPT italic_Q start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) start_POSTSUPERSCRIPT roman_T end_POSTSUPERSCRIPT ( italic_C italic_A start_POSTSUBSCRIPT [ : , italic_i ] end_POSTSUBSCRIPT italic_P start_POSTSUBSCRIPT italic_i , italic_k - 1 end_POSTSUBSCRIPT italic_A start_POSTSUBSCRIPT [ : , italic_i ] end_POSTSUBSCRIPT start_POSTSUPERSCRIPT roman_T end_POSTSUPERSCRIPT italic_C start_POSTSUPERSCRIPT roman_T end_POSTSUPERSCRIPT + italic_C start_POSTSUBSCRIPT [ : , italic_i ] end_POSTSUBSCRIPT italic_Q start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT italic_C start_POSTSUBSCRIPT [ : , italic_i ] end_POSTSUBSCRIPT start_POSTSUPERSCRIPT roman_T end_POSTSUPERSCRIPT + italic_R ) start_POSTSUPERSCRIPT - 1 end_POSTSUPERSCRIPT
Pi,ksubscript𝑃𝑖𝑘\displaystyle P_{i,k}italic_P start_POSTSUBSCRIPT italic_i , italic_k end_POSTSUBSCRIPT =Li,k(CA[:,i]Pi,k1AiiT+C[:,i]Qi)+(AiiPi,k1AiiT+Qi)absentsubscript𝐿𝑖𝑘𝐶subscript𝐴:𝑖subscript𝑃𝑖𝑘1superscriptsubscript𝐴𝑖𝑖Tsubscript𝐶:𝑖subscript𝑄𝑖subscript𝐴𝑖𝑖subscript𝑃𝑖𝑘1superscriptsubscript𝐴𝑖𝑖Tsubscript𝑄𝑖\displaystyle=-L_{i,k}(CA_{[:,i]}P_{i,k-1}A_{ii}^{\mathrm{T}}+C_{[:,i]}Q_{i})+% (A_{ii}P_{i,k-1}A_{ii}^{\mathrm{T}}+Q_{i})= - italic_L start_POSTSUBSCRIPT italic_i , italic_k end_POSTSUBSCRIPT ( italic_C italic_A start_POSTSUBSCRIPT [ : , italic_i ] end_POSTSUBSCRIPT italic_P start_POSTSUBSCRIPT italic_i , italic_k - 1 end_POSTSUBSCRIPT italic_A start_POSTSUBSCRIPT italic_i italic_i end_POSTSUBSCRIPT start_POSTSUPERSCRIPT roman_T end_POSTSUPERSCRIPT + italic_C start_POSTSUBSCRIPT [ : , italic_i ] end_POSTSUBSCRIPT italic_Q start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) + ( italic_A start_POSTSUBSCRIPT italic_i italic_i end_POSTSUBSCRIPT italic_P start_POSTSUBSCRIPT italic_i , italic_k - 1 end_POSTSUBSCRIPT italic_A start_POSTSUBSCRIPT italic_i italic_i end_POSTSUBSCRIPT start_POSTSUPERSCRIPT roman_T end_POSTSUPERSCRIPT + italic_Q start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) (20b)

with A[:,i]subscript𝐴:𝑖A_{[:,i]}italic_A start_POSTSUBSCRIPT [ : , italic_i ] end_POSTSUBSCRIPT and C[:,i]subscript𝐶:𝑖C_{[:,i]}italic_C start_POSTSUBSCRIPT [ : , italic_i ] end_POSTSUBSCRIPT being matrices consisting of the columns of A𝐴Aitalic_A and C𝐶Citalic_C that are associated with zisubscript𝑧𝑖z_{i}italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT, respectively. By neglecting the coupling between the estimation errors of different subsystems, the summation of the arrival costs in (20a) for all the subsystems can be viewed as an approximation of the arrival cost in (19a) for the centralized MHE.

Motivated by the designs of local objective functions for distributed MHE in Schneider and Marquardt [8] and Li et al. [49] for each subsystem i𝑖iitalic_i, i𝕀𝑖𝕀i\in\mathbb{I}italic_i ∈ blackboard_I, we incorporate sensor measurements from interacting subsystems into ΘisubscriptΘ𝑖\Theta_{i}roman_Θ start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT, which contains valuable information for state estimation of the local subsystem i𝑖iitalic_i, such that the cost function for each subsystem i𝑖iitalic_i, denoted by Θ~isubscript~Θ𝑖{\tilde{\Theta}}_{i}over~ start_ARG roman_Θ end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT, is established as follows:

Θ~isubscript~Θ𝑖\displaystyle{\tilde{\Theta}}_{i}over~ start_ARG roman_Θ end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT =Θi+l𝕀i{yl}kNk𝒪[:,i]lz^i(kN|k)j𝕀i𝒪[:,j]lz¯j(kN)Λi{u~}kNk1\displaystyle=\Theta_{i}+\sum_{l\in\mathbb{I}_{i}}\Big{\|}\left\{{y_{l}}\right% \}_{k-N}^{k}-\mathcal{O}^{l}_{[:,i]}{\hat{z}}_{i}(k-N|k)-\sum_{j\in\mathbb{I}_% {i}}\mathcal{O}^{l}_{[:,j]}{\bar{z}}_{j}(k-N)-\Lambda^{i}\left\{{\tilde{u}}% \right\}_{k-N}^{k-1}= roman_Θ start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT + ∑ start_POSTSUBSCRIPT italic_l ∈ blackboard_I start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT ∥ { italic_y start_POSTSUBSCRIPT italic_l end_POSTSUBSCRIPT } start_POSTSUBSCRIPT italic_k - italic_N end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_k end_POSTSUPERSCRIPT - caligraphic_O start_POSTSUPERSCRIPT italic_l end_POSTSUPERSCRIPT start_POSTSUBSCRIPT [ : , italic_i ] end_POSTSUBSCRIPT over^ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k - italic_N | italic_k ) - ∑ start_POSTSUBSCRIPT italic_j ∈ blackboard_I start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT caligraphic_O start_POSTSUPERSCRIPT italic_l end_POSTSUPERSCRIPT start_POSTSUBSCRIPT [ : , italic_j ] end_POSTSUBSCRIPT over¯ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT ( italic_k - italic_N ) - roman_Λ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT { over~ start_ARG italic_u end_ARG } start_POSTSUBSCRIPT italic_k - italic_N end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_k - 1 end_POSTSUPERSCRIPT
Γ[:,i]l{w^i(d|k)}d=kNk1𝑹12evaluated-atsubscriptsuperscriptΓ𝑙:𝑖superscriptsubscriptsubscript^𝑤𝑖conditional𝑑𝑘𝑑𝑘𝑁𝑘1superscript𝑹12\displaystyle\quad-\Gamma^{l}_{[:,i]}\left\{{\hat{w}}_{i}(d|k)\right\}_{d=k-N}% ^{k-1}\Big{\|}^{2}_{\boldsymbol{R}^{-1}}- roman_Γ start_POSTSUPERSCRIPT italic_l end_POSTSUPERSCRIPT start_POSTSUBSCRIPT [ : , italic_i ] end_POSTSUBSCRIPT { over^ start_ARG italic_w end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_d | italic_k ) } start_POSTSUBSCRIPT italic_d = italic_k - italic_N end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_k - 1 end_POSTSUPERSCRIPT ∥ start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT start_POSTSUBSCRIPT bold_italic_R start_POSTSUPERSCRIPT - 1 end_POSTSUPERSCRIPT end_POSTSUBSCRIPT
=z^i(kN|k)z¯i(kN)Pi,kN12+d=kNk1w^i(d|k)Qi12+{y}kNk𝒪[:,i]z^i(kN|k)\displaystyle=\left\|{\hat{z}}_{i}({k-N|k})-{\bar{z}}_{i}(k-N)\right\|^{2}_{P_% {i,k-N}^{-1}}+\sum\limits_{d=k-N}^{k-1}\left\|{\hat{w}_{i}}(d|k)\right\|^{2}_{% Q_{i}^{-1}}+\bigg{\|}\ \left\{{y}\right\}_{k-N}^{k}-\mathcal{O}_{[:,i]}{\hat{z% }}_{i}(k-N|k)= ∥ over^ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k - italic_N | italic_k ) - over¯ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k - italic_N ) ∥ start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_P start_POSTSUBSCRIPT italic_i , italic_k - italic_N end_POSTSUBSCRIPT start_POSTSUPERSCRIPT - 1 end_POSTSUPERSCRIPT end_POSTSUBSCRIPT + ∑ start_POSTSUBSCRIPT italic_d = italic_k - italic_N end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_k - 1 end_POSTSUPERSCRIPT ∥ over^ start_ARG italic_w end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_d | italic_k ) ∥ start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_Q start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT start_POSTSUPERSCRIPT - 1 end_POSTSUPERSCRIPT end_POSTSUBSCRIPT + ∥ { italic_y } start_POSTSUBSCRIPT italic_k - italic_N end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_k end_POSTSUPERSCRIPT - caligraphic_O start_POSTSUBSCRIPT [ : , italic_i ] end_POSTSUBSCRIPT over^ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k - italic_N | italic_k )
j𝕀i𝒪[:,j]z¯j(kN)Λ{u~}kNk1Γ[:,i]{w^i(d|k)}d=kNk1𝑹12subscript𝑗subscript𝕀𝑖subscript𝒪:𝑗subscript¯𝑧𝑗𝑘𝑁Λsuperscriptsubscript~𝑢𝑘𝑁𝑘1evaluated-atsubscriptΓ:𝑖superscriptsubscriptsubscript^𝑤𝑖conditional𝑑𝑘𝑑𝑘𝑁𝑘1superscript𝑹12\displaystyle\quad-\sum_{j\in\mathbb{I}_{i}}\mathcal{O}_{[:,j]}{\bar{z}}_{j}(k% -N)-\Lambda\left\{{\tilde{u}}\right\}_{k-N}^{k-1}-\Gamma_{[:,i]}\left\{{\hat{w% }}_{i}(d|k)\right\}_{d=k-N}^{k-1}\bigg{\|}^{2}_{\boldsymbol{R}^{-1}}- ∑ start_POSTSUBSCRIPT italic_j ∈ blackboard_I start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT caligraphic_O start_POSTSUBSCRIPT [ : , italic_j ] end_POSTSUBSCRIPT over¯ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT ( italic_k - italic_N ) - roman_Λ { over~ start_ARG italic_u end_ARG } start_POSTSUBSCRIPT italic_k - italic_N end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_k - 1 end_POSTSUPERSCRIPT - roman_Γ start_POSTSUBSCRIPT [ : , italic_i ] end_POSTSUBSCRIPT { over^ start_ARG italic_w end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_d | italic_k ) } start_POSTSUBSCRIPT italic_d = italic_k - italic_N end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_k - 1 end_POSTSUPERSCRIPT ∥ start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT start_POSTSUBSCRIPT bold_italic_R start_POSTSUPERSCRIPT - 1 end_POSTSUPERSCRIPT end_POSTSUBSCRIPT

𝑹:=diag{R,,RN+1}assign𝑹diagsubscript𝑅𝑅𝑁1\boldsymbol{R}:=\text{diag}\big{\{}\underbrace{R,\ldots,R}_{N+1}\big{\}}bold_italic_R := diag { under⏟ start_ARG italic_R , … , italic_R end_ARG start_POSTSUBSCRIPT italic_N + 1 end_POSTSUBSCRIPT } with R:=diag{R1,,Rm}assign𝑅diagsubscript𝑅1subscript𝑅𝑚R:=\text{diag}\big{\{}R_{1},\ldots,R_{m}\big{\}}italic_R := diag { italic_R start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT , … , italic_R start_POSTSUBSCRIPT italic_m end_POSTSUBSCRIPT }; Λ,Γ,𝒪ΛΓ𝒪\Lambda,~{}\Gamma,~{}\mathcal{O}roman_Λ , roman_Γ , caligraphic_O are with the following form:

ΛΛ\displaystyle{\Lambda}roman_Λ :=[𝟎𝟎𝟎𝟎CB𝟎𝟎𝟎CAN2BCAN3B𝟎𝟎CAN1BCAN2BCABCB];Γ:=[𝟎𝟎𝟎𝟎C𝟎𝟎𝟎CAN2CAN3𝟎𝟎CAN1CAN2CAC];formulae-sequenceassignabsentdelimited-[]0000missing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpression𝐶𝐵000missing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpression𝐶superscript𝐴𝑁2𝐵𝐶superscript𝐴𝑁3𝐵00missing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpression𝐶superscript𝐴𝑁1𝐵𝐶superscript𝐴𝑁2𝐵𝐶𝐴𝐵𝐶𝐵missing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionassignΓdelimited-[]0000missing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpression𝐶000missing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpression𝐶superscript𝐴𝑁2𝐶superscript𝐴𝑁300missing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpression𝐶superscript𝐴𝑁1𝐶superscript𝐴𝑁2𝐶𝐴𝐶missing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpression\displaystyle:=\left[{\begin{array}[]{*{20}{c}}\boldsymbol{0}&\boldsymbol{0}&% \cdots&\boldsymbol{0}&\boldsymbol{0}\\ CB&\boldsymbol{0}&{\cdots}&\boldsymbol{0}&\boldsymbol{0}\\ \vdots&{\vdots}&{\cdots}&{\vdots}&{\vdots}\\ {C{A^{N-2}}B}&{C{A^{N-3}}B}&{\cdots}&\boldsymbol{0}&\boldsymbol{0}\\ {C{A^{N-1}}B}&{C{A^{N-2}}B}&{\cdots}&CAB&CB\end{array}}\right];~{}~{}{\Gamma}:% =\left[{\begin{array}[]{*{20}{c}}\boldsymbol{0}&\boldsymbol{0}&\cdots&% \boldsymbol{0}&\boldsymbol{0}\\ C&\boldsymbol{0}&{\cdots}&\boldsymbol{0}&\boldsymbol{0}\\ \vdots&{\vdots}&{\cdots}&{\vdots}&{\vdots}\\ {C{A^{N-2}}}&{C{A^{N-3}}}&{\cdots}&\boldsymbol{0}&\boldsymbol{0}\\ {C{A^{N-1}}}&{C{A^{N-2}}}&{\cdots}&CA&C\end{array}}\right];:= [ start_ARRAY start_ROW start_CELL bold_0 end_CELL start_CELL bold_0 end_CELL start_CELL ⋯ end_CELL start_CELL bold_0 end_CELL start_CELL bold_0 end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL end_ROW start_ROW start_CELL italic_C italic_B end_CELL start_CELL bold_0 end_CELL start_CELL ⋯ end_CELL start_CELL bold_0 end_CELL start_CELL bold_0 end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL end_ROW start_ROW start_CELL ⋮ end_CELL start_CELL ⋮ end_CELL start_CELL ⋯ end_CELL start_CELL ⋮ end_CELL start_CELL ⋮ end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL end_ROW start_ROW start_CELL italic_C italic_A start_POSTSUPERSCRIPT italic_N - 2 end_POSTSUPERSCRIPT italic_B end_CELL start_CELL italic_C italic_A start_POSTSUPERSCRIPT italic_N - 3 end_POSTSUPERSCRIPT italic_B end_CELL start_CELL ⋯ end_CELL start_CELL bold_0 end_CELL start_CELL bold_0 end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL end_ROW start_ROW start_CELL italic_C italic_A start_POSTSUPERSCRIPT italic_N - 1 end_POSTSUPERSCRIPT italic_B end_CELL start_CELL italic_C italic_A start_POSTSUPERSCRIPT italic_N - 2 end_POSTSUPERSCRIPT italic_B end_CELL start_CELL ⋯ end_CELL start_CELL italic_C italic_A italic_B end_CELL start_CELL italic_C italic_B end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL end_ROW end_ARRAY ] ; roman_Γ := [ start_ARRAY start_ROW start_CELL bold_0 end_CELL start_CELL bold_0 end_CELL start_CELL ⋯ end_CELL start_CELL bold_0 end_CELL start_CELL bold_0 end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL end_ROW start_ROW start_CELL italic_C end_CELL start_CELL bold_0 end_CELL start_CELL ⋯ end_CELL start_CELL bold_0 end_CELL start_CELL bold_0 end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL end_ROW start_ROW start_CELL ⋮ end_CELL start_CELL ⋮ end_CELL start_CELL ⋯ end_CELL start_CELL ⋮ end_CELL start_CELL ⋮ end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL end_ROW start_ROW start_CELL italic_C italic_A start_POSTSUPERSCRIPT italic_N - 2 end_POSTSUPERSCRIPT end_CELL start_CELL italic_C italic_A start_POSTSUPERSCRIPT italic_N - 3 end_POSTSUPERSCRIPT end_CELL start_CELL ⋯ end_CELL start_CELL bold_0 end_CELL start_CELL bold_0 end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL end_ROW start_ROW start_CELL italic_C italic_A start_POSTSUPERSCRIPT italic_N - 1 end_POSTSUPERSCRIPT end_CELL start_CELL italic_C italic_A start_POSTSUPERSCRIPT italic_N - 2 end_POSTSUPERSCRIPT end_CELL start_CELL ⋯ end_CELL start_CELL italic_C italic_A end_CELL start_CELL italic_C end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL end_ROW end_ARRAY ] ;
𝒪𝒪\displaystyle\mathcal{O}caligraphic_O :=[CT(CA)T(CAN)T]T;assignabsentsuperscriptdelimited-[]superscript𝐶Tsuperscript𝐶𝐴Tsuperscript𝐶superscript𝐴𝑁TT\displaystyle:=\left[C^{\text{T}}~{}~{}\left(CA\right)^{\text{T}}~{}\ldots~{}% \left(CA^{N}\right)^{\text{T}}\right]^{\text{T}};:= [ italic_C start_POSTSUPERSCRIPT T end_POSTSUPERSCRIPT ( italic_C italic_A ) start_POSTSUPERSCRIPT T end_POSTSUPERSCRIPT … ( italic_C italic_A start_POSTSUPERSCRIPT italic_N end_POSTSUPERSCRIPT ) start_POSTSUPERSCRIPT T end_POSTSUPERSCRIPT ] start_POSTSUPERSCRIPT T end_POSTSUPERSCRIPT ;

Γ[:,i]subscriptΓ:𝑖\Gamma_{[:,i]}roman_Γ start_POSTSUBSCRIPT [ : , italic_i ] end_POSTSUBSCRIPT and 𝒪[:,i]subscript𝒪:𝑖\mathcal{O}_{[:,i]}caligraphic_O start_POSTSUBSCRIPT [ : , italic_i ] end_POSTSUBSCRIPT are matrices constituted of the columns of matrices ΓΓ\Gammaroman_Γ and 𝒪𝒪\mathcal{O}caligraphic_O that are associated with wisubscript𝑤𝑖w_{i}italic_w start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT and zisubscript𝑧𝑖z_{i}italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT of subsystem i𝑖iitalic_i, i𝕀𝑖𝕀i\in\mathbb{I}italic_i ∈ blackboard_I, respectively; Γ[:,i]jsuperscriptsubscriptΓ:𝑖𝑗\Gamma_{[:,i]}^{j}roman_Γ start_POSTSUBSCRIPT [ : , italic_i ] end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_j end_POSTSUPERSCRIPT and O[:,i]jsuperscriptsubscript𝑂:𝑖𝑗O_{[:,i]}^{j}italic_O start_POSTSUBSCRIPT [ : , italic_i ] end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_j end_POSTSUPERSCRIPT denote the rows of the matrices Γ[:,i]subscriptΓ:𝑖\Gamma_{[:,i]}roman_Γ start_POSTSUBSCRIPT [ : , italic_i ] end_POSTSUBSCRIPT and O[:,i]subscript𝑂:𝑖O_{[:,i]}italic_O start_POSTSUBSCRIPT [ : , italic_i ] end_POSTSUBSCRIPT with respect to zjsubscript𝑧𝑗z_{j}italic_z start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT, respectively; ΛisuperscriptΛ𝑖\Lambda^{i}roman_Λ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT denotes the rows of the matrix ΛΛ\Lambdaroman_Λ that are associated with zisubscript𝑧𝑖z_{i}italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT.

Formulation of the MHE-based estimators

At each sampling instant kN𝑘𝑁k\geq Nitalic_k ≥ italic_N, Estimator i𝑖iitalic_i (the local MHE-based estimator for the i𝑖iitalic_ith subsystem) will be used to provide the estimates of the subsystem states {z^i(d|k)}d=kNksubscriptsuperscriptsubscript^𝑧𝑖conditional𝑑𝑘𝑘𝑑𝑘𝑁\left\{{\hat{z}_{i}({d|k})}\right\}^{k}_{d=k-N}{ over^ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_d | italic_k ) } start_POSTSUPERSCRIPT italic_k end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_d = italic_k - italic_N end_POSTSUBSCRIPT. Accordingly, the decision variables associated with Estimator i𝑖iitalic_i include z^i(kN|k)subscript^𝑧𝑖𝑘conditional𝑁𝑘{\hat{z}}_{i}(k-N|k)over^ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k - italic_N | italic_k ) and {w^i(d|k)}d=kNk1subscriptsuperscriptsubscript^𝑤𝑖conditional𝑑𝑘𝑘1𝑑𝑘𝑁\left\{{\hat{w}_{i}({d|k})}\right\}^{k-1}_{d=k-N}{ over^ start_ARG italic_w end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_d | italic_k ) } start_POSTSUPERSCRIPT italic_k - 1 end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_d = italic_k - italic_N end_POSTSUBSCRIPT, while z^j(kN|k)subscript^𝑧𝑗𝑘conditional𝑁𝑘{\hat{z}}_{j}(k-N|k)over^ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT ( italic_k - italic_N | italic_k ) and {w^j(d|k)}d=kNk1subscriptsuperscriptsubscript^𝑤𝑗conditional𝑑𝑘𝑘1𝑑𝑘𝑁\left\{{\hat{w}_{j}({d|k})}\right\}^{k-1}_{d=k-N}{ over^ start_ARG italic_w end_ARG start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT ( italic_d | italic_k ) } start_POSTSUPERSCRIPT italic_k - 1 end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_d = italic_k - italic_N end_POSTSUBSCRIPT, j𝕀ifor-all𝑗subscript𝕀𝑖\forall j\in\mathbb{I}_{i}∀ italic_j ∈ blackboard_I start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT will be accounted for by the neighboring local estimators, and are treated as known constants to Estimator i𝑖iitalic_i each time the optimization problem associated with Estimator i𝑖iitalic_i is solved.

Consequently, by taking into account the hard constraints on the original process states x𝑥xitalic_x, the design of each estimator for subsystem i𝑖iitalic_i is formulated as follows:

Design of Local Estimators: At time instant kN𝑘𝑁k\geq Nitalic_k ≥ italic_N, given z¯i(kN)subscript¯𝑧𝑖𝑘𝑁\bar{z}_{i}(k-N)over¯ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k - italic_N ), {y¯}kNksuperscriptsubscript¯𝑦𝑘𝑁𝑘\left\{{\bar{y}}\right\}_{k-N}^{k}{ over¯ start_ARG italic_y end_ARG } start_POSTSUBSCRIPT italic_k - italic_N end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_k end_POSTSUPERSCRIPT , {u~}kNk1superscriptsubscript~𝑢𝑘𝑁𝑘1\left\{{\tilde{u}}\right\}_{k-N}^{k-1}{ over~ start_ARG italic_u end_ARG } start_POSTSUBSCRIPT italic_k - italic_N end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_k - 1 end_POSTSUPERSCRIPT, and z¯j(kN)subscript¯𝑧𝑗𝑘𝑁{\bar{z}}_{j}(k-N)over¯ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT ( italic_k - italic_N ), j𝕀i𝑗subscript𝕀𝑖j\in\mathbb{I}_{i}italic_j ∈ blackboard_I start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT, Estimator i𝑖iitalic_i, i𝕀𝑖𝕀i\in\mathbb{I}italic_i ∈ blackboard_I, solves the following optimization problem:

z^i(kN|k),{w^i(d|k)}d=kNk1=argminz^i(kN|k),{w^i(d|k)}d=kNk1Θ~i(k)subscript^𝑧𝑖𝑘conditional𝑁𝑘subscriptsuperscriptsubscript^𝑤𝑖conditional𝑑𝑘𝑘1𝑑𝑘𝑁subscriptsubscript^𝑧𝑖𝑘conditional𝑁𝑘subscriptsuperscriptsubscript^𝑤𝑖conditional𝑑𝑘𝑘1𝑑𝑘𝑁subscript~Θ𝑖𝑘\displaystyle{\hat{z}_{i}({k-N|k})},\big{\{}{\hat{w}_{i}({d|k})}\big{\}}^{k-1}% _{d=k-N}=\arg\mathop{\min}\limits_{{\hat{z}_{i}({k-N|k})},\left\{{\hat{w}_{i}(% {d|k})}\right\}^{k-1}_{d=k-N}}{\tilde{\Theta}}_{i}(k)over^ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k - italic_N | italic_k ) , { over^ start_ARG italic_w end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_d | italic_k ) } start_POSTSUPERSCRIPT italic_k - 1 end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_d = italic_k - italic_N end_POSTSUBSCRIPT = roman_arg roman_min start_POSTSUBSCRIPT over^ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k - italic_N | italic_k ) , { over^ start_ARG italic_w end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_d | italic_k ) } start_POSTSUPERSCRIPT italic_k - 1 end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_d = italic_k - italic_N end_POSTSUBSCRIPT end_POSTSUBSCRIPT over~ start_ARG roman_Θ end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) (21a)
s.t. 𝐥𝐛i𝒁^(i)𝐮𝐛is.t. superscript𝐥𝐛𝑖superscriptbold-^𝒁𝑖superscript𝐮𝐛𝑖\displaystyle\textmd{\bf{ s.t.~{}}}\qquad\qquad\quad~{}\text{\bf{lb}}^{i}\leq{% \boldsymbol{\hat{Z}}}^{(i)}\leq\text{\bf{ub}}^{i}s.t. lb start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT ≤ overbold_^ start_ARG bold_italic_Z end_ARG start_POSTSUPERSCRIPT ( italic_i ) end_POSTSUPERSCRIPT ≤ ub start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT (21b)

In (21a), the cost function is:

Θ~i(k)subscript~Θ𝑖𝑘\displaystyle{\tilde{\Theta}}_{i}(k)over~ start_ARG roman_Θ end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k ) =z^i(kN|k)z¯i(kN)Pi,kN12+d=kNk1w^i(d|k)Qi12+{y}kNk𝒪[:,i]z^i(kN|k)\displaystyle=\left\|{\hat{z}}_{i}({k-N|k})-{\bar{z}}_{i}(k-N)\right\|^{2}_{P_% {i,k-N}^{-1}}+\sum\limits_{d=k-N}^{k-1}\left\|{\hat{w}_{i}}(d|k)\right\|^{2}_{% Q_{i}^{-1}}+\bigg{\|}\ \left\{{y}\right\}_{k-N}^{k}-\mathcal{O}_{[:,i]}{\hat{z% }}_{i}(k-N|k)= ∥ over^ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k - italic_N | italic_k ) - over¯ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k - italic_N ) ∥ start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_P start_POSTSUBSCRIPT italic_i , italic_k - italic_N end_POSTSUBSCRIPT start_POSTSUPERSCRIPT - 1 end_POSTSUPERSCRIPT end_POSTSUBSCRIPT + ∑ start_POSTSUBSCRIPT italic_d = italic_k - italic_N end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_k - 1 end_POSTSUPERSCRIPT ∥ over^ start_ARG italic_w end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_d | italic_k ) ∥ start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_Q start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT start_POSTSUPERSCRIPT - 1 end_POSTSUPERSCRIPT end_POSTSUBSCRIPT + ∥ { italic_y } start_POSTSUBSCRIPT italic_k - italic_N end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_k end_POSTSUPERSCRIPT - caligraphic_O start_POSTSUBSCRIPT [ : , italic_i ] end_POSTSUBSCRIPT over^ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k - italic_N | italic_k ) (22)
j𝕀i𝒪[:,j]z¯j(kN)Λ{u~}kNk1Γ[:,i]{w^i(d|k)}d=kNk1𝑹12subscript𝑗subscript𝕀𝑖subscript𝒪:𝑗subscript¯𝑧𝑗𝑘𝑁Λsuperscriptsubscript~𝑢𝑘𝑁𝑘1evaluated-atsubscriptΓ:𝑖superscriptsubscriptsubscript^𝑤𝑖conditional𝑑𝑘𝑑𝑘𝑁𝑘1superscript𝑹12\displaystyle~{}~{}~{}~{}-\sum_{j\in\mathbb{I}_{i}}\mathcal{O}_{[:,j]}{\bar{z}% }_{j}(k-N)-\Lambda\left\{{\tilde{u}}\right\}_{k-N}^{k-1}-\Gamma_{[:,i]}\left\{% {\hat{w}}_{i}(d|k)\right\}_{d=k-N}^{k-1}\|^{2}_{\boldsymbol{R}^{-1}}- ∑ start_POSTSUBSCRIPT italic_j ∈ blackboard_I start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT caligraphic_O start_POSTSUBSCRIPT [ : , italic_j ] end_POSTSUBSCRIPT over¯ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT ( italic_k - italic_N ) - roman_Λ { over~ start_ARG italic_u end_ARG } start_POSTSUBSCRIPT italic_k - italic_N end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_k - 1 end_POSTSUPERSCRIPT - roman_Γ start_POSTSUBSCRIPT [ : , italic_i ] end_POSTSUBSCRIPT { over^ start_ARG italic_w end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_d | italic_k ) } start_POSTSUBSCRIPT italic_d = italic_k - italic_N end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_k - 1 end_POSTSUPERSCRIPT ∥ start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT start_POSTSUBSCRIPT bold_italic_R start_POSTSUPERSCRIPT - 1 end_POSTSUPERSCRIPT end_POSTSUBSCRIPT

wherein the weighting matrix Pi,ksubscript𝑃𝑖𝑘P_{i,k}italic_P start_POSTSUBSCRIPT italic_i , italic_k end_POSTSUBSCRIPT is updated following (Partition-based cost function for local MHE-based estimators); z¯i(kN)subscript¯𝑧𝑖𝑘𝑁{\bar{z}}_{i}(k-N)over¯ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k - italic_N ) is an a priori prediction of zi(kN)subscript𝑧𝑖𝑘𝑁z_{i}(k-N)italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k - italic_N ), which is computed using the subsystem model as:

z¯i(kN)=Aiiz^i(kN1|k1)+j𝕀iAijz^j(kN1|k1)+Biu~i(kN1)+w^i(kN1|k1)subscript¯𝑧𝑖𝑘𝑁subscript𝐴𝑖𝑖subscript^𝑧𝑖𝑘𝑁conditional1𝑘1subscript𝑗subscript𝕀𝑖subscript𝐴𝑖𝑗subscript^𝑧𝑗𝑘𝑁conditional1𝑘1subscript𝐵𝑖subscript~𝑢𝑖𝑘𝑁1subscript^𝑤𝑖𝑘𝑁conditional1𝑘1{\bar{z}}_{i}(k-N)=A_{ii}{\hat{z}}_{i}(k-N-1|k-1)+\sum_{j\in\mathbb{I}_{i}}A_{% ij}{\hat{z}}_{j}(k-N-1|k-1)+B_{i}{\tilde{u}}_{i}(k-N-1)+\hat{w}_{i}(k-N-1|k-1)over¯ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k - italic_N ) = italic_A start_POSTSUBSCRIPT italic_i italic_i end_POSTSUBSCRIPT over^ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k - italic_N - 1 | italic_k - 1 ) + ∑ start_POSTSUBSCRIPT italic_j ∈ blackboard_I start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT italic_A start_POSTSUBSCRIPT italic_i italic_j end_POSTSUBSCRIPT over^ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT ( italic_k - italic_N - 1 | italic_k - 1 ) + italic_B start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT over~ start_ARG italic_u end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k - italic_N - 1 ) + over^ start_ARG italic_w end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k - italic_N - 1 | italic_k - 1 ) (23)

where z^i(kN1|k1)subscript^𝑧𝑖𝑘𝑁conditional1𝑘1{\hat{z}}_{i}(k-N-1|k-1)over^ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k - italic_N - 1 | italic_k - 1 ) and w^i(kN1|k1)subscript^𝑤𝑖𝑘𝑁conditional1𝑘1\hat{w}_{i}(k-N-1|k-1)over^ start_ARG italic_w end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k - italic_N - 1 | italic_k - 1 ), respectively, denote the estimates of zi(kN1)subscript𝑧𝑖𝑘𝑁1z_{i}(k-N-1)italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k - italic_N - 1 ) and wi(kN1)subscript𝑤𝑖𝑘𝑁1w_{i}(k-N-1)italic_w start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k - italic_N - 1 ) provided by Estimator i𝑖iitalic_i at time instant k1𝑘1k-1italic_k - 1; N𝑁Nitalic_N is the length of the estimation window of the local estimator.

Equation (21b) is an inequality constraint imposed on the decision variables (i.e., z^i(kN|k)subscript^𝑧𝑖𝑘conditional𝑁𝑘{\hat{z}_{i}({k-N|k})}over^ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k - italic_N | italic_k ) and {w^i(d|k)}d=kNk1subscriptsuperscriptsubscript^𝑤𝑖conditional𝑑𝑘𝑘1𝑑𝑘𝑁\big{\{}{\hat{w}_{i}({d|k})}\big{\}}^{k-1}_{d=k-N}{ over^ start_ARG italic_w end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_d | italic_k ) } start_POSTSUPERSCRIPT italic_k - 1 end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_d = italic_k - italic_N end_POSTSUBSCRIPT) such that the estimates of the lifted states of each subsystem will be constrained within the desired range determined by the constraints on the original process states x𝑥xitalic_x. Specifically, 𝒁^(i):=𝒢[:,i]z^i(kN|k)+j𝕀i𝒢[:,j]z¯j(kN)+{u~}kNk1+𝒥[:,i]{w^i(d|k)}d=kNk1assignsuperscriptbold-^𝒁𝑖subscript𝒢:𝑖subscript^𝑧𝑖𝑘conditional𝑁𝑘subscript𝑗subscript𝕀𝑖subscript𝒢:𝑗subscript¯𝑧𝑗𝑘𝑁superscriptsubscript~𝑢𝑘𝑁𝑘1subscript𝒥:𝑖superscriptsubscriptsubscript^𝑤𝑖conditional𝑑𝑘𝑑𝑘𝑁𝑘1{\boldsymbol{\hat{Z}}}^{(i)}:=\mathcal{G}_{[:,i]}{\hat{z}}_{i}(k-N|k)+\sum_{j% \in\mathbb{I}_{i}}\mathcal{G}_{[:,j]}{\bar{z}}_{j}(k-N)+\mathcal{H}\left\{{% \tilde{u}}\right\}_{k-N}^{k-1}+{\mathcal{J}}_{[:,i]}\left\{{\hat{w}}_{i}(d|k)% \right\}_{d=k-N}^{k-1}overbold_^ start_ARG bold_italic_Z end_ARG start_POSTSUPERSCRIPT ( italic_i ) end_POSTSUPERSCRIPT := caligraphic_G start_POSTSUBSCRIPT [ : , italic_i ] end_POSTSUBSCRIPT over^ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k - italic_N | italic_k ) + ∑ start_POSTSUBSCRIPT italic_j ∈ blackboard_I start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT caligraphic_G start_POSTSUBSCRIPT [ : , italic_j ] end_POSTSUBSCRIPT over¯ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT ( italic_k - italic_N ) + caligraphic_H { over~ start_ARG italic_u end_ARG } start_POSTSUBSCRIPT italic_k - italic_N end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_k - 1 end_POSTSUPERSCRIPT + caligraphic_J start_POSTSUBSCRIPT [ : , italic_i ] end_POSTSUBSCRIPT { over^ start_ARG italic_w end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_d | italic_k ) } start_POSTSUBSCRIPT italic_d = italic_k - italic_N end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_k - 1 end_POSTSUPERSCRIPT, in which:

\displaystyle{\mathcal{H}}caligraphic_H :=[𝟎𝟎𝟎𝟎B𝟎𝟎𝟎AN2BAN3B𝟎𝟎AN1BAN2BABB],𝒥:=[𝟎𝟎𝟎𝟎𝟎𝟎𝟎AN2AN3𝟎𝟎AN1AN2AI],formulae-sequenceassignabsentdelimited-[]0000missing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpression𝐵000missing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionsuperscript𝐴𝑁2𝐵superscript𝐴𝑁3𝐵00missing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionsuperscript𝐴𝑁1𝐵superscript𝐴𝑁2𝐵𝐴𝐵𝐵missing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionassign𝒥delimited-[]0000missing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpression000missing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionsuperscript𝐴𝑁2superscript𝐴𝑁300missing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionsuperscript𝐴𝑁1superscript𝐴𝑁2𝐴𝐼missing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpression\displaystyle:=\left[{\begin{array}[]{*{20}{c}}\boldsymbol{0}&\boldsymbol{0}&% \cdots&\boldsymbol{0}&\boldsymbol{0}\\ B&\boldsymbol{0}&{\cdots}&\boldsymbol{0}&\boldsymbol{0}\\ \vdots&{\vdots}&{\cdots}&{\vdots}&{\vdots}\\ {{A^{N-2}}B}&{{A^{N-3}}B}&{\cdots}&\boldsymbol{0}&\boldsymbol{0}\\ {{A^{N-1}}B}&{{A^{N-2}}B}&{\cdots}&AB&B\end{array}}\right],~{}~{}{\mathcal{J}}% :=\left[{\begin{array}[]{*{20}{c}}\boldsymbol{0}&\boldsymbol{0}&\cdots&% \boldsymbol{0}&\boldsymbol{0}\\ &\boldsymbol{0}&{\cdots}&\boldsymbol{0}&\boldsymbol{0}\\ \vdots&{\vdots}&{\cdots}&{\vdots}&{\vdots}\\ {{A^{N-2}}}&{{A^{N-3}}}&{\cdots}&\boldsymbol{0}&\boldsymbol{0}\\ {{A^{N-1}}}&{{A^{N-2}}}&{\cdots}&A&I\end{array}}\right],:= [ start_ARRAY start_ROW start_CELL bold_0 end_CELL start_CELL bold_0 end_CELL start_CELL ⋯ end_CELL start_CELL bold_0 end_CELL start_CELL bold_0 end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL end_ROW start_ROW start_CELL italic_B end_CELL start_CELL bold_0 end_CELL start_CELL ⋯ end_CELL start_CELL bold_0 end_CELL start_CELL bold_0 end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL end_ROW start_ROW start_CELL ⋮ end_CELL start_CELL ⋮ end_CELL start_CELL ⋯ end_CELL start_CELL ⋮ end_CELL start_CELL ⋮ end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL end_ROW start_ROW start_CELL italic_A start_POSTSUPERSCRIPT italic_N - 2 end_POSTSUPERSCRIPT italic_B end_CELL start_CELL italic_A start_POSTSUPERSCRIPT italic_N - 3 end_POSTSUPERSCRIPT italic_B end_CELL start_CELL ⋯ end_CELL start_CELL bold_0 end_CELL start_CELL bold_0 end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL end_ROW start_ROW start_CELL italic_A start_POSTSUPERSCRIPT italic_N - 1 end_POSTSUPERSCRIPT italic_B end_CELL start_CELL italic_A start_POSTSUPERSCRIPT italic_N - 2 end_POSTSUPERSCRIPT italic_B end_CELL start_CELL ⋯ end_CELL start_CELL italic_A italic_B end_CELL start_CELL italic_B end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL end_ROW end_ARRAY ] , caligraphic_J := [ start_ARRAY start_ROW start_CELL bold_0 end_CELL start_CELL bold_0 end_CELL start_CELL ⋯ end_CELL start_CELL bold_0 end_CELL start_CELL bold_0 end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL end_ROW start_ROW start_CELL end_CELL start_CELL bold_0 end_CELL start_CELL ⋯ end_CELL start_CELL bold_0 end_CELL start_CELL bold_0 end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL end_ROW start_ROW start_CELL ⋮ end_CELL start_CELL ⋮ end_CELL start_CELL ⋯ end_CELL start_CELL ⋮ end_CELL start_CELL ⋮ end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL end_ROW start_ROW start_CELL italic_A start_POSTSUPERSCRIPT italic_N - 2 end_POSTSUPERSCRIPT end_CELL start_CELL italic_A start_POSTSUPERSCRIPT italic_N - 3 end_POSTSUPERSCRIPT end_CELL start_CELL ⋯ end_CELL start_CELL bold_0 end_CELL start_CELL bold_0 end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL end_ROW start_ROW start_CELL italic_A start_POSTSUPERSCRIPT italic_N - 1 end_POSTSUPERSCRIPT end_CELL start_CELL italic_A start_POSTSUPERSCRIPT italic_N - 2 end_POSTSUPERSCRIPT end_CELL start_CELL ⋯ end_CELL start_CELL italic_A end_CELL start_CELL italic_I end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL end_ROW end_ARRAY ] ,
𝒢𝒢\displaystyle\mathcal{G}caligraphic_G :=[I(A)T(AN)T]T,assignabsentsuperscriptdelimited-[]𝐼superscript𝐴Tsuperscriptsuperscript𝐴𝑁TT\displaystyle:=\left[I~{}~{}\left(A\right)^{\text{T}}~{}\ldots~{}\left(A^{N}% \right)^{\text{T}}\right]^{\text{T}},:= [ italic_I ( italic_A ) start_POSTSUPERSCRIPT T end_POSTSUPERSCRIPT … ( italic_A start_POSTSUPERSCRIPT italic_N end_POSTSUPERSCRIPT ) start_POSTSUPERSCRIPT T end_POSTSUPERSCRIPT ] start_POSTSUPERSCRIPT T end_POSTSUPERSCRIPT ,

𝒥[:,i]subscript𝒥:𝑖{\mathcal{J}}_{[:,i]}caligraphic_J start_POSTSUBSCRIPT [ : , italic_i ] end_POSTSUBSCRIPT and 𝒢[:,i]subscript𝒢:𝑖\mathcal{G}_{[:,i]}caligraphic_G start_POSTSUBSCRIPT [ : , italic_i ] end_POSTSUBSCRIPT are matrices composed of the columns of 𝒥𝒥{\mathcal{J}}caligraphic_J and 𝒢𝒢\mathcal{G}caligraphic_G that are associated with wisubscript𝑤𝑖w_{i}italic_w start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT and zisubscript𝑧𝑖z_{i}italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT, i𝕀𝑖𝕀i\in\mathbb{I}italic_i ∈ blackboard_I, respectively. 𝒁^(i)superscriptbold-^𝒁𝑖{\boldsymbol{\hat{Z}}}^{(i)}overbold_^ start_ARG bold_italic_Z end_ARG start_POSTSUPERSCRIPT ( italic_i ) end_POSTSUPERSCRIPT are composed of the estimates of states of the i𝑖iitalic_ith subsystems within the current estimation window, that is, 𝒁^(i)=[z^i(kN|k)T,,z^i(k|k)T]superscriptbold-^𝒁𝑖subscript^𝑧𝑖superscript𝑘conditional𝑁𝑘Tsubscript^𝑧𝑖superscriptconditional𝑘𝑘T{\boldsymbol{\hat{Z}}}^{(i)}=\big{[}{\hat{z}}_{i}(k-N|k)^{\text{T}},~{}\ldots,% ~{}{\hat{z}}_{i}(k|k)^{\text{T}}\big{]}overbold_^ start_ARG bold_italic_Z end_ARG start_POSTSUPERSCRIPT ( italic_i ) end_POSTSUPERSCRIPT = [ over^ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k - italic_N | italic_k ) start_POSTSUPERSCRIPT T end_POSTSUPERSCRIPT , … , over^ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k | italic_k ) start_POSTSUPERSCRIPT T end_POSTSUPERSCRIPT ]. Meanwhile, for Estimator i𝑖iitalic_i, we only need to impose constraints on z^isubscript^𝑧𝑖{\hat{z}}_{i}over^ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT, as the estimates of the other subsystems are handled by the neighboring estimators. Accordingly, 𝐥𝐛isuperscript𝐥𝐛𝑖\text{\bf{lb}}^{i}lb start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT and 𝐮𝐛isuperscript𝐮𝐛𝑖\text{\bf{ub}}^{i}ub start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT in (21b), which are the lower and upper bounds on 𝒁^(i)superscriptbold-^𝒁𝑖{\boldsymbol{\hat{Z}}}^{(i)}overbold_^ start_ARG bold_italic_Z end_ARG start_POSTSUPERSCRIPT ( italic_i ) end_POSTSUPERSCRIPT, are determined in a way such that only constraints are imposed on the elements corresponding to the first nxisubscript𝑛subscript𝑥𝑖n_{x_{i}}italic_n start_POSTSUBSCRIPT italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT states of z^isubscript^𝑧𝑖{\hat{z}}_{i}over^ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT, that is, z^i(d|k)[1:nxi]subscript^𝑧𝑖subscriptconditional𝑑𝑘delimited-[]:1subscript𝑛subscript𝑥𝑖{\hat{z}}_{i}{({d|k})_{[1:n_{x_{i}}]}}over^ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_d | italic_k ) start_POSTSUBSCRIPT [ 1 : italic_n start_POSTSUBSCRIPT italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT ] end_POSTSUBSCRIPT, for d=kN,,k𝑑𝑘𝑁𝑘d=k-N,\ldots,kitalic_d = italic_k - italic_N , … , italic_k. The remaining elements in 𝐥𝐛isuperscript𝐥𝐛𝑖\text{\bf{lb}}^{i}lb start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT and 𝐮𝐛isuperscript𝐮𝐛𝑖\text{\bf{ub}}^{i}ub start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT that are corresponding to each subsystem j𝑗jitalic_j, j𝕀i𝑗subscript𝕀𝑖j\in\mathbb{I}_{i}italic_j ∈ blackboard_I start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT, are set to be -\infty- ∞ and ++\infty+ ∞, respectively.

By incorporating constraint (21b) with the defined 𝐥𝐛isuperscript𝐥𝐛𝑖\text{\bf{lb}}^{i}lb start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT and 𝐮𝐛isuperscript𝐮𝐛𝑖\text{\bf{ub}}^{i}ub start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT into each local MHE-based estimator, it is ensured that z^i(d|k)[1:nxi]𝕏isubscript^𝑧𝑖subscriptconditional𝑑𝑘delimited-[]:1subscript𝑛subscript𝑥𝑖subscript𝕏𝑖{\hat{z}}_{i}({{d|k})_{[1:n_{x_{i}}]}}\in\mathbb{X}_{i}over^ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_d | italic_k ) start_POSTSUBSCRIPT [ 1 : italic_n start_POSTSUBSCRIPT italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT ] end_POSTSUBSCRIPT ∈ blackboard_X start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT, for d=kN,,k𝑑𝑘𝑁𝑘d=k-N,\ldots,kitalic_d = italic_k - italic_N , … , italic_k. Since these nxisubscript𝑛subscript𝑥𝑖n_{x_{i}}italic_n start_POSTSUBSCRIPT italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUBSCRIPT elements of the lifted state vector zisubscript𝑧𝑖z_{i}italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT are made identical to xisubscript𝑥𝑖x_{i}italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT according to Algorithm 1), this constraint ensures that the estimates of the original states of each subsystem xisubscript𝑥𝑖x_{i}italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT satisfy the physical hard constraints on them.

Note that the estimator design is presented without considering the scaling of the data samples. Since the data used in parallel subsystem identification is scaled in implementation, the estimates provided by the local estimators of the distributed MHE scheme do not account for the actual process states. Therefore, the determination of the bounds 𝐥𝐛isuperscript𝐥𝐛𝑖\text{\bf{lb}}^{i}lb start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT and 𝐮𝐛isuperscript𝐮𝐛𝑖\text{\bf{ub}}^{i}ub start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT based on the constraints on the original subsystem state xisubscript𝑥𝑖x_{i}italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT is also subject to scaling. In addition, the state estimates provided by the MHE-based estimators for the subsystems need to be unscaled, and the unscaled values will be used to represent the estimates of the actual process states.

Implementation algorithm

We summarize the implementation steps for the proposed distributed estimation method in an algorithm as follows:

Algorithm 2.

Implementation procedure for the distributed MHE scheme.

At time kN+1𝑘𝑁1{k\geq N+1}italic_k ≥ italic_N + 1, do the following steps:

  1. 1.

    Estimator i𝑖iitalic_i, i𝕀𝑖𝕀i\in\mathbb{I}italic_i ∈ blackboard_I, receives sensor measurements from all the subsystems y(k)𝑦𝑘y(k)italic_y ( italic_k ).

  2. 2.

    Estimator i𝑖iitalic_i, i𝕀𝑖𝕀i\in\mathbb{I}italic_i ∈ blackboard_I, computes z¯i(kN)subscript¯𝑧𝑖𝑘𝑁{\bar{z}}_{i}(k-N)over¯ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k - italic_N ) following (23).

    • 2.1

      Estimator i𝑖iitalic_i, i𝕀𝑖𝕀i\in\mathbb{I}italic_i ∈ blackboard_I, requests and receives z¯j(kN|k)subscript¯𝑧𝑗𝑘conditional𝑁𝑘\bar{z}_{j}(k-N|k)over¯ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT ( italic_k - italic_N | italic_k ) (i.e., the prediction of the lifted states of subsystem j𝑗jitalic_j), from the each Estimator j𝑗jitalic_j, j𝕀i𝑗subscript𝕀𝑖j\in\mathbb{I}_{i}italic_j ∈ blackboard_I start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT.

    • 2.2.

      Estimator i𝑖iitalic_i, i𝕀𝑖𝕀i\in\mathbb{I}italic_i ∈ blackboard_I, solves (Formulation of the MHE-based estimators) to provide optimal estimates of the subsystem state and disturbances sequence, i.e., z^i(kN|k)subscript^𝑧𝑖𝑘conditional𝑁𝑘{\hat{z}_{i}({k-N|k})}over^ start_ARG italic_z end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_k - italic_N | italic_k ) and {w^i(d|k)}d=kNk1subscriptsuperscriptsubscript^𝑤𝑖conditional𝑑𝑘𝑘1𝑑𝑘𝑁\big{\{}{\hat{w}_{i}({d|k})}\big{\}}^{k-1}_{d=k-N}{ over^ start_ARG italic_w end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ( italic_d | italic_k ) } start_POSTSUPERSCRIPT italic_k - 1 end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_d = italic_k - italic_N end_POSTSUBSCRIPT.

  3. 3.

    A sequence of optimal state estimates z^(d|k)^𝑧conditional𝑑𝑘\hat{z}({d|k})over^ start_ARG italic_z end_ARG ( italic_d | italic_k ), d=kN,k𝑑𝑘𝑁𝑘d=k-N\ldots,kitalic_d = italic_k - italic_N … , italic_k, is generated for the state vector z𝑧zitalic_z of entire process in the lifted space.

  4. 4.

    Update the weighting matrix Pi,ksubscript𝑃𝑖𝑘P_{i,k}italic_P start_POSTSUBSCRIPT italic_i , italic_k end_POSTSUBSCRIPT following (Partition-based cost function for local MHE-based estimators). Set k=k+1𝑘𝑘1k=k+1italic_k = italic_k + 1, go to Step 1.

Remark 4.

With convex sets that bound the estimates of the subsystem states and disturbances, the execution of the local estimators can be implemented as a convex quadratic problem instead of solving non-convex and computationally expensive optimization problems associated with nonlinear MHE. This way, as compared to the case when estimators are developed based on first-principles nonlinear models, the computational efficiency of the proposed distributed estimation scheme is also expected to be improved despite the dimensionality increase of the subsystem models.

Application to a chemical process

In this section, the proposed framework is illustrated via a simulated application to a chemical process consisting of four continuous stirred tank reactors.

Process description

This process involves four continuous stirred tank reactors (CSTRs) that are interconnected via mass and energy flows. A schematic of the process is presented in Figure 3. A feed stream, containing pure reactant A𝐴Aitalic_A, enters the i𝑖iitalic_ith reactor with flow rate F0isubscript𝐹0𝑖F_{0i}italic_F start_POSTSUBSCRIPT 0 italic_i end_POSTSUBSCRIPT, temperature T0isubscript𝑇0𝑖T_{0i}italic_T start_POSTSUBSCRIPT 0 italic_i end_POSTSUBSCRIPT, and molar concentration CA0isubscript𝐶𝐴0𝑖C_{A0i}italic_C start_POSTSUBSCRIPT italic_A 0 italic_i end_POSTSUBSCRIPT. The effluent of the first reactor is fed into the second reactor at flow rate F1subscript𝐹1F_{1}italic_F start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT, temperature T1subscript𝑇1T_{1}italic_T start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT, and molar concentration CA1subscript𝐶𝐴1C_{A1}italic_C start_POSTSUBSCRIPT italic_A 1 end_POSTSUBSCRIPT. A fraction of the outlet stream from the second reactor is sent to the third reactor, with flow rate F2subscript𝐹2F_{2}italic_F start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT, temperature T2subscript𝑇2T_{2}italic_T start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT, and concentration CA2subscript𝐶𝐴2C_{A2}italic_C start_POSTSUBSCRIPT italic_A 2 end_POSTSUBSCRIPT. The remaining portion of the effluent of the second reactor is sent back to the first reactor at flow rate Fr1subscript𝐹𝑟1F_{r1}italic_F start_POSTSUBSCRIPT italic_r 1 end_POSTSUBSCRIPT, temperature T2subscript𝑇2T_{2}italic_T start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT, and concentration CA2subscript𝐶𝐴2C_{A2}italic_C start_POSTSUBSCRIPT italic_A 2 end_POSTSUBSCRIPT. The outlet of the third reactor enters the fourth reactor at flow rate F3subscript𝐹3F_{3}italic_F start_POSTSUBSCRIPT 3 end_POSTSUBSCRIPT, temperature T3subscript𝑇3T_{3}italic_T start_POSTSUBSCRIPT 3 end_POSTSUBSCRIPT, and molar concentration CA3subscript𝐶𝐴3C_{A3}italic_C start_POSTSUBSCRIPT italic_A 3 end_POSTSUBSCRIPT. A portion of the outlet stream from the fourth reactor is recycled back to the first reactor at flow rate Fr2subscript𝐹𝑟2F_{r2}italic_F start_POSTSUBSCRIPT italic_r 2 end_POSTSUBSCRIPT, temperature T4subscript𝑇4T_{4}italic_T start_POSTSUBSCRIPT 4 end_POSTSUBSCRIPT, and concentration CA4subscript𝐶𝐴4C_{A4}italic_C start_POSTSUBSCRIPT italic_A 4 end_POSTSUBSCRIPT. Each of the four reactors is equipped with a jacket, which enables the addition or extraction of heat to or from its respective reactor.

Based on mass and energy balances, a first-principles dynamic model that consists of eight differential equations is presented as follows [50]:

dT1dtdsubscript𝑇1d𝑡\displaystyle\frac{\mathrm{d}T_{1}}{\mathrm{d}t}divide start_ARG roman_d italic_T start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT end_ARG start_ARG roman_d italic_t end_ARG =F01V1(T01T1)+Fr1V1(T2T1)+Fr2V1(T4T1)i=13ΔHiρcpki0eEiRT1CA1+Q1ρcpV1absentsubscript𝐹01subscript𝑉1subscript𝑇01subscript𝑇1subscript𝐹𝑟1subscript𝑉1subscript𝑇2subscript𝑇1subscript𝐹𝑟2subscript𝑉1subscript𝑇4subscript𝑇1superscriptsubscript𝑖13Δsubscript𝐻𝑖𝜌subscript𝑐𝑝subscript𝑘𝑖0superscript𝑒subscript𝐸𝑖𝑅subscript𝑇1subscript𝐶𝐴1subscript𝑄1𝜌subscript𝑐𝑝subscript𝑉1\displaystyle=\frac{F_{01}}{V_{1}}(T_{01}-T_{1})+\frac{F_{r1}}{V_{1}}(T_{2}-T_% {1})+\frac{F_{r2}}{V_{1}}(T_{4}-T_{1})-\sum_{i=1}^{3}\frac{\Delta H_{i}}{\rho c% _{p}}k_{i0}e^{\frac{-E_{i}}{RT_{1}}}C_{A1}+\frac{Q_{1}}{\rho c_{p}V_{1}}= divide start_ARG italic_F start_POSTSUBSCRIPT 01 end_POSTSUBSCRIPT end_ARG start_ARG italic_V start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT end_ARG ( italic_T start_POSTSUBSCRIPT 01 end_POSTSUBSCRIPT - italic_T start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT ) + divide start_ARG italic_F start_POSTSUBSCRIPT italic_r 1 end_POSTSUBSCRIPT end_ARG start_ARG italic_V start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT end_ARG ( italic_T start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT - italic_T start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT ) + divide start_ARG italic_F start_POSTSUBSCRIPT italic_r 2 end_POSTSUBSCRIPT end_ARG start_ARG italic_V start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT end_ARG ( italic_T start_POSTSUBSCRIPT 4 end_POSTSUBSCRIPT - italic_T start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT ) - ∑ start_POSTSUBSCRIPT italic_i = 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT 3 end_POSTSUPERSCRIPT divide start_ARG roman_Δ italic_H start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_ARG start_ARG italic_ρ italic_c start_POSTSUBSCRIPT italic_p end_POSTSUBSCRIPT end_ARG italic_k start_POSTSUBSCRIPT italic_i 0 end_POSTSUBSCRIPT italic_e start_POSTSUPERSCRIPT divide start_ARG - italic_E start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_ARG start_ARG italic_R italic_T start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT end_ARG end_POSTSUPERSCRIPT italic_C start_POSTSUBSCRIPT italic_A 1 end_POSTSUBSCRIPT + divide start_ARG italic_Q start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT end_ARG start_ARG italic_ρ italic_c start_POSTSUBSCRIPT italic_p end_POSTSUBSCRIPT italic_V start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT end_ARG (24a)
dCA1dtdsubscript𝐶𝐴1d𝑡\displaystyle\frac{\mathrm{d}C_{A1}}{\mathrm{d}t}divide start_ARG roman_d italic_C start_POSTSUBSCRIPT italic_A 1 end_POSTSUBSCRIPT end_ARG start_ARG roman_d italic_t end_ARG =F01V1(CA01CA1)+Fr1V1(CA2CA1)+Fr2V1(CA4CA1)i=13ki0eEiRT1CA1absentsubscript𝐹01subscript𝑉1subscript𝐶𝐴01subscript𝐶𝐴1subscript𝐹𝑟1subscript𝑉1subscript𝐶𝐴2subscript𝐶𝐴1subscript𝐹𝑟2subscript𝑉1subscript𝐶𝐴4subscript𝐶𝐴1superscriptsubscript𝑖13subscript𝑘𝑖0superscript𝑒subscript𝐸𝑖𝑅subscript𝑇1subscript𝐶𝐴1\displaystyle=\frac{F_{01}}{V_{1}}(C_{A01}-C_{A1})+\frac{F_{r1}}{V_{1}}(C_{A2}% -C_{A1})+\frac{F_{r2}}{V_{1}}(C_{A4}-C_{A1})-\sum_{i=1}^{3}k_{i0}e^{\frac{-E_{% i}}{RT_{1}}}C_{A1}= divide start_ARG italic_F start_POSTSUBSCRIPT 01 end_POSTSUBSCRIPT end_ARG start_ARG italic_V start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT end_ARG ( italic_C start_POSTSUBSCRIPT italic_A 01 end_POSTSUBSCRIPT - italic_C start_POSTSUBSCRIPT italic_A 1 end_POSTSUBSCRIPT ) + divide start_ARG italic_F start_POSTSUBSCRIPT italic_r 1 end_POSTSUBSCRIPT end_ARG start_ARG italic_V start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT end_ARG ( italic_C start_POSTSUBSCRIPT italic_A 2 end_POSTSUBSCRIPT - italic_C start_POSTSUBSCRIPT italic_A 1 end_POSTSUBSCRIPT ) + divide start_ARG italic_F start_POSTSUBSCRIPT italic_r 2 end_POSTSUBSCRIPT end_ARG start_ARG italic_V start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT end_ARG ( italic_C start_POSTSUBSCRIPT italic_A 4 end_POSTSUBSCRIPT - italic_C start_POSTSUBSCRIPT italic_A 1 end_POSTSUBSCRIPT ) - ∑ start_POSTSUBSCRIPT italic_i = 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT 3 end_POSTSUPERSCRIPT italic_k start_POSTSUBSCRIPT italic_i 0 end_POSTSUBSCRIPT italic_e start_POSTSUPERSCRIPT divide start_ARG - italic_E start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_ARG start_ARG italic_R italic_T start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT end_ARG end_POSTSUPERSCRIPT italic_C start_POSTSUBSCRIPT italic_A 1 end_POSTSUBSCRIPT (24b)
dT2dtdsubscript𝑇2d𝑡\displaystyle\frac{\mathrm{d}T_{2}}{\mathrm{d}t}divide start_ARG roman_d italic_T start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT end_ARG start_ARG roman_d italic_t end_ARG =F1V2(T1T2)+F02V2(T02T2)i=13ΔHiρcpki0eEiRT2CA2+Q2ρcpV2absentsubscript𝐹1subscript𝑉2subscript𝑇1subscript𝑇2subscript𝐹02subscript𝑉2subscript𝑇02subscript𝑇2superscriptsubscript𝑖13Δsubscript𝐻𝑖𝜌subscript𝑐𝑝subscript𝑘𝑖0superscript𝑒subscript𝐸𝑖𝑅subscript𝑇2subscript𝐶𝐴2subscript𝑄2𝜌subscript𝑐𝑝subscript𝑉2\displaystyle=\frac{F_{1}}{V_{2}}(T_{1}-T_{2})+\frac{F_{02}}{V_{2}}(T_{02}-T_{% 2})-\sum_{i=1}^{3}\frac{\Delta H_{i}}{\rho c_{p}}k_{i0}e^{\frac{-E_{i}}{RT_{2}% }}C_{A2}+\frac{Q_{2}}{\rho c_{p}V_{2}}= divide start_ARG italic_F start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT end_ARG start_ARG italic_V start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT end_ARG ( italic_T start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT - italic_T start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT ) + divide start_ARG italic_F start_POSTSUBSCRIPT 02 end_POSTSUBSCRIPT end_ARG start_ARG italic_V start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT end_ARG ( italic_T start_POSTSUBSCRIPT 02 end_POSTSUBSCRIPT - italic_T start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT ) - ∑ start_POSTSUBSCRIPT italic_i = 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT 3 end_POSTSUPERSCRIPT divide start_ARG roman_Δ italic_H start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_ARG start_ARG italic_ρ italic_c start_POSTSUBSCRIPT italic_p end_POSTSUBSCRIPT end_ARG italic_k start_POSTSUBSCRIPT italic_i 0 end_POSTSUBSCRIPT italic_e start_POSTSUPERSCRIPT divide start_ARG - italic_E start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_ARG start_ARG italic_R italic_T start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT end_ARG end_POSTSUPERSCRIPT italic_C start_POSTSUBSCRIPT italic_A 2 end_POSTSUBSCRIPT + divide start_ARG italic_Q start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT end_ARG start_ARG italic_ρ italic_c start_POSTSUBSCRIPT italic_p end_POSTSUBSCRIPT italic_V start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT end_ARG (24c)
dCA2dtdsubscript𝐶𝐴2d𝑡\displaystyle\frac{\mathrm{d}C_{A2}}{\mathrm{d}t}divide start_ARG roman_d italic_C start_POSTSUBSCRIPT italic_A 2 end_POSTSUBSCRIPT end_ARG start_ARG roman_d italic_t end_ARG =F1V2(CA1CA2)+F02V2(CA02CA2)i=13ki0eEiRT2CA2absentsubscript𝐹1subscript𝑉2subscript𝐶𝐴1subscript𝐶𝐴2subscript𝐹02subscript𝑉2subscript𝐶𝐴02subscript𝐶𝐴2superscriptsubscript𝑖13subscript𝑘𝑖0superscript𝑒subscript𝐸𝑖𝑅subscript𝑇2subscript𝐶𝐴2\displaystyle=\frac{F_{1}}{V_{2}}(C_{A1}-C_{A2})+\frac{F_{02}}{V_{2}}(C_{A02}-% C_{A2})-\sum_{i=1}^{3}k_{i0}e^{\frac{-E_{i}}{RT_{2}}}C_{A2}= divide start_ARG italic_F start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT end_ARG start_ARG italic_V start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT end_ARG ( italic_C start_POSTSUBSCRIPT italic_A 1 end_POSTSUBSCRIPT - italic_C start_POSTSUBSCRIPT italic_A 2 end_POSTSUBSCRIPT ) + divide start_ARG italic_F start_POSTSUBSCRIPT 02 end_POSTSUBSCRIPT end_ARG start_ARG italic_V start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT end_ARG ( italic_C start_POSTSUBSCRIPT italic_A 02 end_POSTSUBSCRIPT - italic_C start_POSTSUBSCRIPT italic_A 2 end_POSTSUBSCRIPT ) - ∑ start_POSTSUBSCRIPT italic_i = 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT 3 end_POSTSUPERSCRIPT italic_k start_POSTSUBSCRIPT italic_i 0 end_POSTSUBSCRIPT italic_e start_POSTSUPERSCRIPT divide start_ARG - italic_E start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_ARG start_ARG italic_R italic_T start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT end_ARG end_POSTSUPERSCRIPT italic_C start_POSTSUBSCRIPT italic_A 2 end_POSTSUBSCRIPT (24d)
dT3dtdsubscript𝑇3d𝑡\displaystyle\frac{\mathrm{d}T_{3}}{\mathrm{d}t}divide start_ARG roman_d italic_T start_POSTSUBSCRIPT 3 end_POSTSUBSCRIPT end_ARG start_ARG roman_d italic_t end_ARG =F2Fr1V3(T2T3)+F03V3(T03T3)i=13ΔHiρcpki0eEiRT3CA3+Q3ρcpV3absentsubscript𝐹2subscript𝐹𝑟1subscript𝑉3subscript𝑇2subscript𝑇3subscript𝐹03subscript𝑉3subscript𝑇03subscript𝑇3superscriptsubscript𝑖13Δsubscript𝐻𝑖𝜌subscript𝑐𝑝subscript𝑘𝑖0superscript𝑒subscript𝐸𝑖𝑅subscript𝑇3subscript𝐶𝐴3subscript𝑄3𝜌subscript𝑐𝑝subscript𝑉3\displaystyle=\frac{F_{2}-F_{r1}}{V_{3}}(T_{2}-T_{3})+\frac{F_{03}}{V_{3}}(T_{% 03}-T_{3})-\sum_{i=1}^{3}\frac{\Delta H_{i}}{\rho c_{p}}k_{i0}e^{\frac{-E_{i}}% {RT_{3}}}C_{A3}+\frac{Q_{3}}{\rho c_{p}V_{3}}= divide start_ARG italic_F start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT - italic_F start_POSTSUBSCRIPT italic_r 1 end_POSTSUBSCRIPT end_ARG start_ARG italic_V start_POSTSUBSCRIPT 3 end_POSTSUBSCRIPT end_ARG ( italic_T start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT - italic_T start_POSTSUBSCRIPT 3 end_POSTSUBSCRIPT ) + divide start_ARG italic_F start_POSTSUBSCRIPT 03 end_POSTSUBSCRIPT end_ARG start_ARG italic_V start_POSTSUBSCRIPT 3 end_POSTSUBSCRIPT end_ARG ( italic_T start_POSTSUBSCRIPT 03 end_POSTSUBSCRIPT - italic_T start_POSTSUBSCRIPT 3 end_POSTSUBSCRIPT ) - ∑ start_POSTSUBSCRIPT italic_i = 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT 3 end_POSTSUPERSCRIPT divide start_ARG roman_Δ italic_H start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_ARG start_ARG italic_ρ italic_c start_POSTSUBSCRIPT italic_p end_POSTSUBSCRIPT end_ARG italic_k start_POSTSUBSCRIPT italic_i 0 end_POSTSUBSCRIPT italic_e start_POSTSUPERSCRIPT divide start_ARG - italic_E start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_ARG start_ARG italic_R italic_T start_POSTSUBSCRIPT 3 end_POSTSUBSCRIPT end_ARG end_POSTSUPERSCRIPT italic_C start_POSTSUBSCRIPT italic_A 3 end_POSTSUBSCRIPT + divide start_ARG italic_Q start_POSTSUBSCRIPT 3 end_POSTSUBSCRIPT end_ARG start_ARG italic_ρ italic_c start_POSTSUBSCRIPT italic_p end_POSTSUBSCRIPT italic_V start_POSTSUBSCRIPT 3 end_POSTSUBSCRIPT end_ARG (24e)
dCA3dtdsubscript𝐶𝐴3d𝑡\displaystyle\frac{\mathrm{d}C_{A3}}{\mathrm{d}t}divide start_ARG roman_d italic_C start_POSTSUBSCRIPT italic_A 3 end_POSTSUBSCRIPT end_ARG start_ARG roman_d italic_t end_ARG =F2Fr1V3(CA2CA3)+F03V3(CA03CA3)i=13ki0eEiRT3CA3absentsubscript𝐹2subscript𝐹𝑟1subscript𝑉3subscript𝐶𝐴2subscript𝐶𝐴3subscript𝐹03subscript𝑉3subscript𝐶𝐴03subscript𝐶𝐴3superscriptsubscript𝑖13subscript𝑘𝑖0superscript𝑒subscript𝐸𝑖𝑅subscript𝑇3subscript𝐶𝐴3\displaystyle=\frac{F_{2}-F_{r1}}{V_{3}}(C_{A2}-C_{A3})+\frac{F_{03}}{V_{3}}(C% _{A03}-C_{A3})-\sum_{i=1}^{3}k_{i0}e^{\frac{-E_{i}}{RT_{3}}}C_{A3}= divide start_ARG italic_F start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT - italic_F start_POSTSUBSCRIPT italic_r 1 end_POSTSUBSCRIPT end_ARG start_ARG italic_V start_POSTSUBSCRIPT 3 end_POSTSUBSCRIPT end_ARG ( italic_C start_POSTSUBSCRIPT italic_A 2 end_POSTSUBSCRIPT - italic_C start_POSTSUBSCRIPT italic_A 3 end_POSTSUBSCRIPT ) + divide start_ARG italic_F start_POSTSUBSCRIPT 03 end_POSTSUBSCRIPT end_ARG start_ARG italic_V start_POSTSUBSCRIPT 3 end_POSTSUBSCRIPT end_ARG ( italic_C start_POSTSUBSCRIPT italic_A 03 end_POSTSUBSCRIPT - italic_C start_POSTSUBSCRIPT italic_A 3 end_POSTSUBSCRIPT ) - ∑ start_POSTSUBSCRIPT italic_i = 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT 3 end_POSTSUPERSCRIPT italic_k start_POSTSUBSCRIPT italic_i 0 end_POSTSUBSCRIPT italic_e start_POSTSUPERSCRIPT divide start_ARG - italic_E start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_ARG start_ARG italic_R italic_T start_POSTSUBSCRIPT 3 end_POSTSUBSCRIPT end_ARG end_POSTSUPERSCRIPT italic_C start_POSTSUBSCRIPT italic_A 3 end_POSTSUBSCRIPT (24f)
dT4dtdsubscript𝑇4d𝑡\displaystyle\frac{\mathrm{d}T_{4}}{\mathrm{d}t}divide start_ARG roman_d italic_T start_POSTSUBSCRIPT 4 end_POSTSUBSCRIPT end_ARG start_ARG roman_d italic_t end_ARG =F3V4(T3T4)+F04V4(T04T4)i=13ΔHiρcpki0eEiRT4CA4+Q4ρcpV4absentsubscript𝐹3subscript𝑉4subscript𝑇3subscript𝑇4subscript𝐹04subscript𝑉4subscript𝑇04subscript𝑇4superscriptsubscript𝑖13Δsubscript𝐻𝑖𝜌subscript𝑐𝑝subscript𝑘𝑖0superscript𝑒subscript𝐸𝑖𝑅subscript𝑇4subscript𝐶𝐴4subscript𝑄4𝜌subscript𝑐𝑝subscript𝑉4\displaystyle=\frac{F_{3}}{V_{4}}(T_{3}-T_{4})+\frac{F_{04}}{V_{4}}(T_{04}-T_{% 4})-\sum_{i=1}^{3}\frac{\Delta H_{i}}{\rho c_{p}}k_{i0}e^{\frac{-E_{i}}{RT_{4}% }}C_{A4}+\frac{Q_{4}}{\rho c_{p}V_{4}}= divide start_ARG italic_F start_POSTSUBSCRIPT 3 end_POSTSUBSCRIPT end_ARG start_ARG italic_V start_POSTSUBSCRIPT 4 end_POSTSUBSCRIPT end_ARG ( italic_T start_POSTSUBSCRIPT 3 end_POSTSUBSCRIPT - italic_T start_POSTSUBSCRIPT 4 end_POSTSUBSCRIPT ) + divide start_ARG italic_F start_POSTSUBSCRIPT 04 end_POSTSUBSCRIPT end_ARG start_ARG italic_V start_POSTSUBSCRIPT 4 end_POSTSUBSCRIPT end_ARG ( italic_T start_POSTSUBSCRIPT 04 end_POSTSUBSCRIPT - italic_T start_POSTSUBSCRIPT 4 end_POSTSUBSCRIPT ) - ∑ start_POSTSUBSCRIPT italic_i = 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT 3 end_POSTSUPERSCRIPT divide start_ARG roman_Δ italic_H start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_ARG start_ARG italic_ρ italic_c start_POSTSUBSCRIPT italic_p end_POSTSUBSCRIPT end_ARG italic_k start_POSTSUBSCRIPT italic_i 0 end_POSTSUBSCRIPT italic_e start_POSTSUPERSCRIPT divide start_ARG - italic_E start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_ARG start_ARG italic_R italic_T start_POSTSUBSCRIPT 4 end_POSTSUBSCRIPT end_ARG end_POSTSUPERSCRIPT italic_C start_POSTSUBSCRIPT italic_A 4 end_POSTSUBSCRIPT + divide start_ARG italic_Q start_POSTSUBSCRIPT 4 end_POSTSUBSCRIPT end_ARG start_ARG italic_ρ italic_c start_POSTSUBSCRIPT italic_p end_POSTSUBSCRIPT italic_V start_POSTSUBSCRIPT 4 end_POSTSUBSCRIPT end_ARG (24g)
dCA4dtdsubscript𝐶𝐴4d𝑡\displaystyle\frac{\mathrm{d}C_{A4}}{\mathrm{d}t}divide start_ARG roman_d italic_C start_POSTSUBSCRIPT italic_A 4 end_POSTSUBSCRIPT end_ARG start_ARG roman_d italic_t end_ARG =F3V4(CA3CA4)+F04V4(CA04CA4)i=13ki0eEiRT4CA4absentsubscript𝐹3subscript𝑉4subscript𝐶𝐴3subscript𝐶𝐴4subscript𝐹04subscript𝑉4subscript𝐶𝐴04subscript𝐶𝐴4superscriptsubscript𝑖13subscript𝑘𝑖0superscript𝑒subscript𝐸𝑖𝑅subscript𝑇4subscript𝐶𝐴4\displaystyle=\frac{F_{3}}{V_{4}}(C_{A3}-C_{A4})+\frac{F_{04}}{V_{4}}(C_{A04}-% C_{A4})-\sum_{i=1}^{3}k_{i0}e^{\frac{-E_{i}}{RT_{4}}}C_{A4}= divide start_ARG italic_F start_POSTSUBSCRIPT 3 end_POSTSUBSCRIPT end_ARG start_ARG italic_V start_POSTSUBSCRIPT 4 end_POSTSUBSCRIPT end_ARG ( italic_C start_POSTSUBSCRIPT italic_A 3 end_POSTSUBSCRIPT - italic_C start_POSTSUBSCRIPT italic_A 4 end_POSTSUBSCRIPT ) + divide start_ARG italic_F start_POSTSUBSCRIPT 04 end_POSTSUBSCRIPT end_ARG start_ARG italic_V start_POSTSUBSCRIPT 4 end_POSTSUBSCRIPT end_ARG ( italic_C start_POSTSUBSCRIPT italic_A 04 end_POSTSUBSCRIPT - italic_C start_POSTSUBSCRIPT italic_A 4 end_POSTSUBSCRIPT ) - ∑ start_POSTSUBSCRIPT italic_i = 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT 3 end_POSTSUPERSCRIPT italic_k start_POSTSUBSCRIPT italic_i 0 end_POSTSUBSCRIPT italic_e start_POSTSUPERSCRIPT divide start_ARG - italic_E start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_ARG start_ARG italic_R italic_T start_POSTSUBSCRIPT 4 end_POSTSUBSCRIPT end_ARG end_POSTSUPERSCRIPT italic_C start_POSTSUBSCRIPT italic_A 4 end_POSTSUBSCRIPT (24h)

where CAisubscript𝐶𝐴𝑖C_{Ai}italic_C start_POSTSUBSCRIPT italic_A italic_i end_POSTSUBSCRIPT represents the molar concentration of A𝐴Aitalic_A in i𝑖iitalic_ith reactor, i=1,2,3,4𝑖1234i=1,2,3,4italic_i = 1 , 2 , 3 , 4; Tisubscript𝑇𝑖T_{i}italic_T start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT is the temperature of feed stream in i𝑖iitalic_ith reactor, i=1,2,3,4𝑖1234i=1,2,3,4italic_i = 1 , 2 , 3 , 4. The definitions of the remaining variables and parameters, and the values of the parameters involved in model (Process description) can be found in Yin and Liu [13] and Rashedi et al.[50]

Problem formulation and simulation setting

In this process, the temperatures in the four reactors can be measured using sensors, while the concentrations of the reactant A𝐴Aitalic_A need to be estimated online. To achieve this objective, the process is decomposed into four subsystems according to the physical topology of the process, that is, each reactor is one subsystem. Then, the proposed method is applied: 1) to build Koopman subsystem models for the four subsystems to collaboratively characterize the dynamics of the entire process; 2) to establish a data-driven distributed state estimation scheme for estimating the full-state of the process. In light of the objective, in this work, the nonlinear model in (Process description) is not used for designing state estimators. Instead, it is only leveraged to generate data which will be used to establish data-driven Koopman subsystem models based on the proposed method.

The values of the model parameters in (Process description) are made the same as those in Yin and Liu [13]. The initial condition of the simulation for generating state trajectories is [326.3794K3.1833kmol/m3[326.3794~{}\mathrm{K}~{}~{}3.1833~{}\mathrm{kmol/m^{3}}[ 326.3794 roman_K 3.1833 roman_kmol / roman_m start_POSTSUPERSCRIPT 3 end_POSTSUPERSCRIPT 326.3745K2.9402kmol/m3328.0896K2.9863kmol/m3326.7154K3.1649kmol/m3]T326.3745~{}\mathrm{K}~{}~{}2.9402~{}\mathrm{kmol/m^{3}}~{}~{}328.0896~{}% \mathrm{K}~{}~{}2.9863~{}\mathrm{kmol/m^{3}}~{}~{}326.7154~{}\mathrm{K}~{}~{}3% .1649~{}\mathrm{kmol/m^{3}}]^{\mathrm{T}}326.3745 roman_K 2.9402 roman_kmol / roman_m start_POSTSUPERSCRIPT 3 end_POSTSUPERSCRIPT 328.0896 roman_K 2.9863 roman_kmol / roman_m start_POSTSUPERSCRIPT 3 end_POSTSUPERSCRIPT 326.7154 roman_K 3.1649 roman_kmol / roman_m start_POSTSUPERSCRIPT 3 end_POSTSUPERSCRIPT ] start_POSTSUPERSCRIPT roman_T end_POSTSUPERSCRIPT. The sequence of heating input Qisubscript𝑄𝑖Q_{i}italic_Q start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT, i=1,2,3,4𝑖1234i=1,2,3,4italic_i = 1 , 2 , 3 , 4, is generated following a uniform distribution within predefined bounds, and it varies every 1.5 hours. Specifically, the upper bounds on the heating inputs to the four reactors are determined as Q1,max=1.2×104kJ/hsubscript𝑄11.2superscript104kJhQ_{1,\max}=1.2\times 10^{4}~{}\mathrm{kJ/h}italic_Q start_POSTSUBSCRIPT 1 , roman_max end_POSTSUBSCRIPT = 1.2 × 10 start_POSTSUPERSCRIPT 4 end_POSTSUPERSCRIPT roman_kJ / roman_h, Q2,max=2.2×104kJ/hsubscript𝑄22.2superscript104kJhQ_{2,\max}=2.2\times 10^{4}~{}\mathrm{kJ/h}italic_Q start_POSTSUBSCRIPT 2 , roman_max end_POSTSUBSCRIPT = 2.2 × 10 start_POSTSUPERSCRIPT 4 end_POSTSUPERSCRIPT roman_kJ / roman_h, Q3,max=2.7×104kJ/hsubscript𝑄32.7superscript104kJhQ_{3,\max}=2.7\times 10^{4}~{}\mathrm{kJ/h}italic_Q start_POSTSUBSCRIPT 3 , roman_max end_POSTSUBSCRIPT = 2.7 × 10 start_POSTSUPERSCRIPT 4 end_POSTSUPERSCRIPT roman_kJ / roman_h, and Q4,max=1.2×104kJ/hsubscript𝑄41.2superscript104kJhQ_{4,\max}=1.2\times 10^{4}~{}\mathrm{kJ/h}italic_Q start_POSTSUBSCRIPT 4 , roman_max end_POSTSUBSCRIPT = 1.2 × 10 start_POSTSUPERSCRIPT 4 end_POSTSUPERSCRIPT roman_kJ / roman_h, while the lower bounds of the heating inputs are set as Q1,min=0.8×104kJ/hsubscript𝑄10.8superscript104kJhQ_{1,\min}=0.8\times 10^{4}~{}\mathrm{kJ/h}italic_Q start_POSTSUBSCRIPT 1 , roman_min end_POSTSUBSCRIPT = 0.8 × 10 start_POSTSUPERSCRIPT 4 end_POSTSUPERSCRIPT roman_kJ / roman_h, Q2,min=1.8×104kJ/hsubscript𝑄21.8superscript104kJhQ_{2,\min}=1.8\times 10^{4}~{}\mathrm{kJ/h}italic_Q start_POSTSUBSCRIPT 2 , roman_min end_POSTSUBSCRIPT = 1.8 × 10 start_POSTSUPERSCRIPT 4 end_POSTSUPERSCRIPT roman_kJ / roman_h, Q3,min=2.3×104kJ/hsubscript𝑄32.3superscript104kJhQ_{3,\min}=2.3\times 10^{4}~{}\mathrm{kJ/h}italic_Q start_POSTSUBSCRIPT 3 , roman_min end_POSTSUBSCRIPT = 2.3 × 10 start_POSTSUPERSCRIPT 4 end_POSTSUPERSCRIPT roman_kJ / roman_h, and Q4,min=0.8×104kJ/hsubscript𝑄40.8superscript104kJhQ_{4,\min}=0.8\times 10^{4}~{}\mathrm{kJ/h}italic_Q start_POSTSUBSCRIPT 4 , roman_min end_POSTSUBSCRIPT = 0.8 × 10 start_POSTSUPERSCRIPT 4 end_POSTSUPERSCRIPT roman_kJ / roman_h. Unknown system disturbances are generated following Gaussian distribution with zero mean and standard deviation of θwsubscript𝜃𝑤\theta_{w}italic_θ start_POSTSUBSCRIPT italic_w end_POSTSUBSCRIPT and are then made bounded within [w¯,w¯]¯𝑤¯𝑤[-\bar{w},\bar{w}][ - over¯ start_ARG italic_w end_ARG , over¯ start_ARG italic_w end_ARG ], where w¯=5θw¯𝑤5subscript𝜃𝑤\bar{w}=5\theta_{w}over¯ start_ARG italic_w end_ARG = 5 italic_θ start_POSTSUBSCRIPT italic_w end_POSTSUBSCRIPT with θw=[0.1554K0.0015kmol/m30.1554K0.0014kmol/m30.1562K0.0014kmol/m30.1556K\theta_{w}=[0.1554~{}\mathrm{K}~{}~{}0.0015~{}\mathrm{kmol/m^{3}}~{}~{}0.1554~% {}\mathrm{K}~{}~{}0.0014~{}\mathrm{kmol/m^{3}}~{}~{}0.1562~{}\mathrm{K}~{}~{}0% .0014~{}\mathrm{kmol/m^{3}}~{}~{}0.1556~{}\mathrm{K}italic_θ start_POSTSUBSCRIPT italic_w end_POSTSUBSCRIPT = [ 0.1554 roman_K 0.0015 roman_kmol / roman_m start_POSTSUPERSCRIPT 3 end_POSTSUPERSCRIPT 0.1554 roman_K 0.0014 roman_kmol / roman_m start_POSTSUPERSCRIPT 3 end_POSTSUPERSCRIPT 0.1562 roman_K 0.0014 roman_kmol / roman_m start_POSTSUPERSCRIPT 3 end_POSTSUPERSCRIPT 0.1556 roman_K 0.0015kmol/m3]T0.0015~{}\mathrm{kmol/m^{3}}]^{\mathrm{T}}0.0015 roman_kmol / roman_m start_POSTSUPERSCRIPT 3 end_POSTSUPERSCRIPT ] start_POSTSUPERSCRIPT roman_T end_POSTSUPERSCRIPT. Similarly, unknown measurement noise is generated following Gaussian distribution with zero mean and standard deviation of θvsubscript𝜃𝑣\theta_{v}italic_θ start_POSTSUBSCRIPT italic_v end_POSTSUBSCRIPT and is then made bounded within [v¯,v¯]¯𝑣¯𝑣[-\bar{v},\bar{v}][ - over¯ start_ARG italic_v end_ARG , over¯ start_ARG italic_v end_ARG ], where v¯=5θv¯𝑣5subscript𝜃𝑣\bar{v}=5\theta_{v}over¯ start_ARG italic_v end_ARG = 5 italic_θ start_POSTSUBSCRIPT italic_v end_POSTSUBSCRIPT with θv=[0.3108K0.3108K0.3125K0.3112K]Tsubscript𝜃𝑣superscriptdelimited-[]0.3108K0.3108K0.3125K0.3112KT\theta_{v}=[0.3108~{}\mathrm{K}~{}~{}0.3108~{}\mathrm{K}~{}~{}0.3125~{}\mathrm% {K}~{}~{}0.3112~{}\mathrm{K}]^{\mathrm{T}}italic_θ start_POSTSUBSCRIPT italic_v end_POSTSUBSCRIPT = [ 0.3108 roman_K 0.3108 roman_K 0.3125 roman_K 0.3112 roman_K ] start_POSTSUPERSCRIPT roman_T end_POSTSUPERSCRIPT. The sampling period is selected as 0.025 h. 2000 samples within an operating period of 50 hours are recorded. The 2000 samples are divided into three datasets. The first dataset contains the first 1000 samples within the first 25-hour period, which will be used for Koopman subsystem modeling. The remaining 1000 samples are divided evenly into two groups: 500 samples, which constitute the second dataset, will be used to validate the data-driven subsystem models; the other 500 samples constitute the third dataset, which will be used to verify the performance of the proposed data-driven distributed state estimation scheme.

In this subsection, we apply the parallel subsystem modeling approach proposed in the “EDMD-based parallel Koopman subsystem identification method” subsection to establish four Koopman models for the subsystems of the chemical process. Given that the states assigned to different subsystems have the same physical meaning, we utilize the same state and input lifting functions ϕisuperscriptitalic-ϕ𝑖\phi^{i}italic_ϕ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT and δisuperscript𝛿𝑖\delta^{i}italic_δ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT, i=1,2,3,4𝑖1234i=1,2,3,4italic_i = 1 , 2 , 3 , 4, for different subsystems. Specifically, the lifting functions are chosen as the original subsystem state xisubscript𝑥𝑖x_{i}italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT, power function xi33subscript𝑥𝑖\sqrt[3]{x_{i}}nth-root start_ARG 3 end_ARG start_ARG italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_ARG, and exponential function exisuperscript𝑒subscript𝑥𝑖e^{x_{i}}italic_e start_POSTSUPERSCRIPT italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUPERSCRIPT, that is,

ϕi=[TiCAiTi3CAi3eTieCAi]T,i=1,2,3,4.formulae-sequencesuperscriptitalic-ϕ𝑖superscriptdelimited-[]subscript𝑇𝑖subscript𝐶𝐴𝑖3subscript𝑇𝑖3subscript𝐶𝐴𝑖superscript𝑒subscript𝑇𝑖superscript𝑒subscript𝐶𝐴𝑖missing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionTfor-all𝑖1234\phi^{i}=\left[{\begin{array}[]{*{20}{c}}T_{i}&C_{Ai}&\sqrt[3]{T_{i}}&\sqrt[3]% {C_{Ai}}&{e^{T_{i}}}&{e^{C_{Ai}}}\end{array}}\right]^{\text{T}},~{}\forall i=1% ,2,3,4.italic_ϕ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT = [ start_ARRAY start_ROW start_CELL italic_T start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_CELL start_CELL italic_C start_POSTSUBSCRIPT italic_A italic_i end_POSTSUBSCRIPT end_CELL start_CELL nth-root start_ARG 3 end_ARG start_ARG italic_T start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_ARG end_CELL start_CELL nth-root start_ARG 3 end_ARG start_ARG italic_C start_POSTSUBSCRIPT italic_A italic_i end_POSTSUBSCRIPT end_ARG end_CELL start_CELL italic_e start_POSTSUPERSCRIPT italic_T start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_POSTSUPERSCRIPT end_CELL start_CELL italic_e start_POSTSUPERSCRIPT italic_C start_POSTSUBSCRIPT italic_A italic_i end_POSTSUBSCRIPT end_POSTSUPERSCRIPT end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL end_ROW end_ARRAY ] start_POSTSUPERSCRIPT T end_POSTSUPERSCRIPT , ∀ italic_i = 1 , 2 , 3 , 4 .

The input lifting functions for each subsystem are selected as δi(Qi)=[QiQi3]T,i=1,2,3,4.formulae-sequencesuperscript𝛿𝑖subscript𝑄𝑖superscriptdelimited-[]subscript𝑄𝑖3subscript𝑄𝑖Tfor-all𝑖1234\delta^{i}(Q_{i})=\left[\begin{array}[]{cc}Q_{i}&\sqrt[3]{Q_{i}}\end{array}% \right]^{\text{T}},~{}\forall i=1,2,3,4.italic_δ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT ( italic_Q start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) = [ start_ARRAY start_ROW start_CELL italic_Q start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_CELL start_CELL nth-root start_ARG 3 end_ARG start_ARG italic_Q start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT end_ARG end_CELL end_ROW end_ARRAY ] start_POSTSUPERSCRIPT T end_POSTSUPERSCRIPT , ∀ italic_i = 1 , 2 , 3 , 4 .

Koopman subsystem modeling

The 1000 samples in the first dataset are scaled following Steps 2.1 and 2.2 of Algorithm 1, and then are utilized for identifying the Koopman-based subsystem models. Based on the selected state and input lifting functions, four subsystem models in the form of (Subsystem models for the design of local estimators) are established by executing Algorithm 1. For each subsystem model, the lifted state vector zisubscript𝑧𝑖z_{i}italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT has 6 variables; the lifted input vector u~isubscript~𝑢𝑖\tilde{u}_{i}over~ start_ARG italic_u end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT for the i𝑖iitalic_ith subsystem, i=1,2,3,4𝑖1234i=1,2,3,4italic_i = 1 , 2 , 3 , 4, has 2 variables.

Cross-validation is conducted based on the second dataset which has 500 samples. The trajectories of the actual states and the open-loop predictions of the states provided by the four Koopman-based subsystem models are presented in Figure 4. Without incorporating state measurements into generating the state predictions, the Koopman subsystem models are able to provide accurate predictions for the temporal trajectories of all eight state variables.

Distributed state estimation results

Based on the established Koopman subsystem models, a distributed state estimation scheme consisting of four MHE-based local estimators is developed. Each local estimator is designed following (Formulation of the MHE-based estimators). The length of the estimation horizon for each local estimator is N=3𝑁3N=3italic_N = 3. The weighting matrices Pi,0subscript𝑃𝑖0P_{i,0}italic_P start_POSTSUBSCRIPT italic_i , 0 end_POSTSUBSCRIPT, Qisubscript𝑄𝑖Q_{i}italic_Q start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT, and Risubscript𝑅𝑖R_{i}italic_R start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT are selected as Pi,0=diag{0.01,0.01,0.01,0.01,0.01,0.01}subscript𝑃𝑖0diag0.010.010.010.010.010.01P_{i,0}=\mathrm{diag}\{0.01,0.01,0.01,0.01,0.01,0.01\}italic_P start_POSTSUBSCRIPT italic_i , 0 end_POSTSUBSCRIPT = roman_diag { 0.01 , 0.01 , 0.01 , 0.01 , 0.01 , 0.01 }, Qi=diag{0.1,0.1,0.1,0.1,0.1,0.1}subscript𝑄𝑖diag0.10.10.10.10.10.1Q_{i}=\mathrm{diag}\{0.1,0.1,0.1,0.1,0.1,0.1\}italic_Q start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT = roman_diag { 0.1 , 0.1 , 0.1 , 0.1 , 0.1 , 0.1 }, and Ri=diag{0.001,0.001,0.001}subscript𝑅𝑖diag0.0010.0010.001R_{i}=\mathrm{diag}\{0.001,0.001,0.001\}italic_R start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT = roman_diag { 0.001 , 0.001 , 0.001 }, for i=1,2,3,4𝑖1234i=1,2,3,4italic_i = 1 , 2 , 3 , 4. It is noted that the estimators provide estimates of the scaled states, which are further unscaled to reconstruct the original process states.

The initial state of the process and the initial guess values adopted by the local estimators are presented in Table 1. Hard constraints are incorporated into the local estimators of the distributed scheme, such that the estimates of the molar concentrations of the reactants are ensured to be non-negative. The trajectories of the actual states and estimates of the full-state of the process are shown in Figure 5. The results confirm that the developed data-driven distributed estimation scheme can provide accurate estimates of the eight process states in the presence of disturbances and measurement noise.

To demonstrate the superiority of the proposed method, we compare the proposed estimation design with a distributed MHE scheme where the local estimators are developed based on the linearized subsystem models obtained via Taylor series expansions. Specifically, we consider a steady-state operating point xs=[310.8376K3.0317kmol/m3x_{s}=[310.8376~{}\mathrm{K}~{}~{}3.0317~{}\mathrm{kmol/m^{3}}italic_x start_POSTSUBSCRIPT italic_s end_POSTSUBSCRIPT = [ 310.8376 roman_K 3.0317 roman_kmol / roman_m start_POSTSUPERSCRIPT 3 end_POSTSUPERSCRIPT 310.8329K2.8002kmol/m3310.8329K2.8002kmolsuperscriptm3310.8329~{}\mathrm{K}~{}~{}2.8002~{}\mathrm{kmol/m^{3}}310.8329 roman_K 2.8002 roman_kmol / roman_m start_POSTSUPERSCRIPT 3 end_POSTSUPERSCRIPT 312.4663K2.844kmol/m3312.4663K2.844kmolsuperscriptm3312.4663~{}\mathrm{K}~{}~{}2.844~{}\mathrm{kmol/m^{3}}312.4663 roman_K 2.844 roman_kmol / roman_m start_POSTSUPERSCRIPT 3 end_POSTSUPERSCRIPT 311.1576K3.0142kmol/m3]T311.1576\,\mathrm{K}~{}3.0142~{}\mathrm{kmol/m^{3}}]^{\mathrm{T}}311.1576 roman_K 3.0142 roman_kmol / roman_m start_POSTSUPERSCRIPT 3 end_POSTSUPERSCRIPT ] start_POSTSUPERSCRIPT roman_T end_POSTSUPERSCRIPT, which corresponds to constant inputs Q1=1.0×104kJ/hsubscript𝑄11.0superscript104kJhQ_{1}=1.0\times 10^{4}~{}\mathrm{kJ/h}italic_Q start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT = 1.0 × 10 start_POSTSUPERSCRIPT 4 end_POSTSUPERSCRIPT roman_kJ / roman_h, Q2=2.0×104kJ/hsubscript𝑄22.0superscript104kJhQ_{2}=2.0\times 10^{4}~{}\mathrm{kJ/h}italic_Q start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT = 2.0 × 10 start_POSTSUPERSCRIPT 4 end_POSTSUPERSCRIPT roman_kJ / roman_h, Q3=2.5×104kJ/hsubscript𝑄32.5superscript104kJhQ_{3}=2.5\times 10^{4}~{}\mathrm{kJ/h}italic_Q start_POSTSUBSCRIPT 3 end_POSTSUBSCRIPT = 2.5 × 10 start_POSTSUPERSCRIPT 4 end_POSTSUPERSCRIPT roman_kJ / roman_h, and Q4=1.0×104kJ/hsubscript𝑄41.0superscript104kJhQ_{4}=1.0\times 10^{4}~{}\mathrm{kJ/h}italic_Q start_POSTSUBSCRIPT 4 end_POSTSUBSCRIPT = 1.0 × 10 start_POSTSUPERSCRIPT 4 end_POSTSUPERSCRIPT roman_kJ / roman_h. At xssubscript𝑥𝑠x_{s}italic_x start_POSTSUBSCRIPT italic_s end_POSTSUBSCRIPT, a linearized model is obtained using Taylor-series expansion. The linearized model is further discretized with a step size h=1.5min1.5minh=1.5~{}\mathrm{min}italic_h = 1.5 roman_min, and a discrete-time model of the entire system is obtained. Then, the linearized global model is partitioned into linearized subsystem models in the form of (Problem formulation) for the subsystems of the entire process. The root mean square error (RMSE) and the average computation time for the one-step execution of one local estimator for both designs are presented in Table 2. Compared to the distributed estimation scheme based on the linearized subsystem models, the estimation error of the proposed method is significantly smaller. This is attributed to the much better predictability of the Koopman subsystem models when the operation is distant from the steady-state point at which the linearized subsystem models are obtained. It is mentioned that the RMSE is computed based on the scaled process states and the corresponding estimates, such that equitable importance is given to the states with different physical meanings. The average computation time for both designs, as presented in Table 2, confirms that both designs are computationally efficient due to the use of linear subsystem models despite nonlinear dynamics, ensuring the feasibility of online implementation. Meanwhile, the distributed estimation scheme based on the linearized subsystem models is indeed more computationally efficient, which is consistent with expectations. The increase in the computation time stems from the higher-dimensional nature of the Koopman subsystem models in comparison to the linearized subsystem models, which, however, is at the expense of a significant reduction in the estimation accuracy.

Remark 5.

The performance of the proposed method is dependent on the noise level in the dataset and the data sampling frequency. The established Koopman-based linear models can accurately capture the dynamic behaviors of the underlying nonlinear process when the noise remains within a reasonable range. In the meantime, the data frequency also can impact the performance of our data-driven method for both subsystem modeling and distributed state estimation of the proposed framework. The data needed for implementing the proposed subsystem modeling and distributed state estimation methods are expected to meet two criteria: 1) they should encompass informative transient trends of the dynamic process; 2) they should contain timely output measurements to ensure the real-time implementation of the proposed distributed moving horizon estimation approach, in accordance with the state estimate frequency requirement.

Application to an agro-hydrological process

In this subsection, the proposed methods on parallel subsystem identification and distributed MHE are applied to an agro-hydrological process through simulations for estimating the soil moisture at different locations of the soil profile. The modeling of this process and real-time estimation of the soil moisture are critical for informed closed-loop irrigation for sustainable water use and improved crop health.

Description of the process

The agro-hydrological process incorporates the vadose zone of the soil. The dynamic interactions among the soil, crops, and the atmospheric environment within the hydrological cycle are considered. A schematic of this process is shown in Figure 6. Water flow entering/leaving the soil is the known input to the process; this water flow can be caused by irrigation via sprinklers, rainfall, water extraction by the crop roots, evaporation, water runoff, and/or drainage.

Two standard assumptions are made for the agro-hydrological process: (1) soil properties, including color, texture, structure, porosity, and density, are homogeneous horizontally; (2) surface irrigation is conducted on the field uniformly.

First-principles model and discretization for data generation

The aim is to leverage the proposed parallel subsystem model identification and distributed estimation method to the agro-hydrological process for data-driven subsystem modeling and distributed soil moisture estimation. In this subsection, we introduce a first-principles model based on the Richards equation [51]; this high-fidelity model is only used for generating process data needed for implementing the proposed method.

To describe the vertical hydrological dynamical behaviors of the process, the Richards equation [51] is introduced as follows:

ht=1C(h)z[K(h)(hz+1)]𝑡1𝐶𝑧delimited-[]𝐾𝑧1\frac{\partial h}{\partial t}=\frac{1}{C(h)}\frac{\partial}{\partial z}\left[K% (h)\left(\frac{\partial h}{\partial z}+1\right)\right]divide start_ARG ∂ italic_h end_ARG start_ARG ∂ italic_t end_ARG = divide start_ARG 1 end_ARG start_ARG italic_C ( italic_h ) end_ARG divide start_ARG ∂ end_ARG start_ARG ∂ italic_z end_ARG [ italic_K ( italic_h ) ( divide start_ARG ∂ italic_h end_ARG start_ARG ∂ italic_z end_ARG + 1 ) ] (25)

where hhitalic_h (m) is the soil water pressure head; z𝑧zitalic_z (m) is the vertical position inside the soil profile; K𝐾Kitalic_K (m/s) denotes the hydraulic conductivity which is dependent on K𝐾Kitalic_K. The dependence of K𝐾Kitalic_K on hhitalic_h is characterized by the following equation [52]:

K(h)=KsatSeλ[1(1Se1m)m]2𝐾subscript𝐾𝑠𝑎𝑡superscriptsubscript𝑆𝑒𝜆superscriptdelimited-[]1superscript1superscriptsubscript𝑆𝑒1𝑚𝑚2K(h)=K_{sat}S_{e}^{\lambda}\left[1-\left(1-S_{e}^{\frac{1}{m}}\right)^{m}% \right]^{2}italic_K ( italic_h ) = italic_K start_POSTSUBSCRIPT italic_s italic_a italic_t end_POSTSUBSCRIPT italic_S start_POSTSUBSCRIPT italic_e end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_λ end_POSTSUPERSCRIPT [ 1 - ( 1 - italic_S start_POSTSUBSCRIPT italic_e end_POSTSUBSCRIPT start_POSTSUPERSCRIPT divide start_ARG 1 end_ARG start_ARG italic_m end_ARG end_POSTSUPERSCRIPT ) start_POSTSUPERSCRIPT italic_m end_POSTSUPERSCRIPT ] start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT (26)

In (26), m=11/n𝑚11𝑛m=1-1/nitalic_m = 1 - 1 / italic_n, λ𝜆\lambdaitalic_λ, α𝛼\alphaitalic_α and n𝑛nitalic_n are empirical parameters of the soil, Ksatsubscript𝐾𝑠𝑎𝑡K_{sat}italic_K start_POSTSUBSCRIPT italic_s italic_a italic_t end_POSTSUBSCRIPT denotes the saturated hydraulic conductivity, and Sesubscript𝑆𝑒S_{e}italic_S start_POSTSUBSCRIPT italic_e end_POSTSUBSCRIPT denotes the relative saturation which can be computed following:

Se=θ(h)θrθsθrsubscript𝑆𝑒𝜃subscript𝜃𝑟subscript𝜃𝑠subscript𝜃𝑟S_{e}=\frac{\theta(h)-\theta_{r}}{\theta_{s}-\theta_{r}}italic_S start_POSTSUBSCRIPT italic_e end_POSTSUBSCRIPT = divide start_ARG italic_θ ( italic_h ) - italic_θ start_POSTSUBSCRIPT italic_r end_POSTSUBSCRIPT end_ARG start_ARG italic_θ start_POSTSUBSCRIPT italic_s end_POSTSUBSCRIPT - italic_θ start_POSTSUBSCRIPT italic_r end_POSTSUBSCRIPT end_ARG

where θssubscript𝜃𝑠\theta_{s}italic_θ start_POSTSUBSCRIPT italic_s end_POSTSUBSCRIPT represents the saturated soil water content; θrsubscript𝜃𝑟\theta_{r}italic_θ start_POSTSUBSCRIPT italic_r end_POSTSUBSCRIPT represents the residual soil water content. The relationship between soil moisture content θ(h)𝜃\theta(h)italic_θ ( italic_h ) (m/m) and the pressure head can be characterized by the soil-water retention equation [53]:

θ(h)=(θsθr)(1+(αh)n)(11n)+θr𝜃subscript𝜃𝑠subscript𝜃𝑟superscript1superscript𝛼𝑛11𝑛subscript𝜃𝑟\theta(h)=\left(\theta_{s}-\theta_{r}\right)\left(1+(-\alpha h)^{n}\right)^{-% \left(1-\frac{1}{n}\right)}+\theta_{r}italic_θ ( italic_h ) = ( italic_θ start_POSTSUBSCRIPT italic_s end_POSTSUBSCRIPT - italic_θ start_POSTSUBSCRIPT italic_r end_POSTSUBSCRIPT ) ( 1 + ( - italic_α italic_h ) start_POSTSUPERSCRIPT italic_n end_POSTSUPERSCRIPT ) start_POSTSUPERSCRIPT - ( 1 - divide start_ARG 1 end_ARG start_ARG italic_n end_ARG ) end_POSTSUPERSCRIPT + italic_θ start_POSTSUBSCRIPT italic_r end_POSTSUBSCRIPT (27)

Further, the capillary capacity (i.e., C(h)𝐶C(h)italic_C ( italic_h ) in Eq. (25)) can be calculated as [53]:

C(h)=nα(θsθr)(11n)(αh)n1(1+(αh)n)(21n)𝐶𝑛𝛼subscript𝜃𝑠subscript𝜃𝑟11𝑛superscript𝛼𝑛1superscript1superscript𝛼𝑛21𝑛C(h)=n\alpha\left(\theta_{s}-\theta_{r}\right)\Big{(}1-\frac{1}{n}\Big{)}\left% (-\alpha h\right)^{n-1}\left(1+\left(-\alpha h\right)^{n}\right)^{-\left(2-% \frac{1}{n}\right)}italic_C ( italic_h ) = italic_n italic_α ( italic_θ start_POSTSUBSCRIPT italic_s end_POSTSUBSCRIPT - italic_θ start_POSTSUBSCRIPT italic_r end_POSTSUBSCRIPT ) ( 1 - divide start_ARG 1 end_ARG start_ARG italic_n end_ARG ) ( - italic_α italic_h ) start_POSTSUPERSCRIPT italic_n - 1 end_POSTSUPERSCRIPT ( 1 + ( - italic_α italic_h ) start_POSTSUPERSCRIPT italic_n end_POSTSUPERSCRIPT ) start_POSTSUPERSCRIPT - ( 2 - divide start_ARG 1 end_ARG start_ARG italic_n end_ARG ) end_POSTSUPERSCRIPT (28)

The Richards equation in (25) is a nonlinear partial differential equation (PDE) concerning both time and vertical coordinates. To obtain information on soil moisture at finite specific locations, the PDE-based model can be discretized. Specifically, finite difference is used to conduct time and space discretization to calculate numerical approximations of (25). The soil moisture in each compartment remains homogeneous horizontally. This ensures the generated data can represent an accurate numerical approximation of the Richards equation in (25). Following the discretization procedure and the boundary conditions adopted in Yin et al.[19] and Bo et al. [54], a two-point forward difference can be conducted to approximate the time derivatives, and a two-point central difference can be conducted to approximate the spatial derivatives. With discretization in both time and space, the dynamic model in (25) can be simulated to generate soil moisture samples for different locations of the soil profile. The generated samples will be used to identify the Koopman-based subsystem models in parallel.

Process setting and data generation

The values of the parameters involved in (25) are the same as those adopted in Carsel and Parrish [55]. A loam soil profile with a depth of L=120𝐿120L=120italic_L = 120 cm is considered. In terms of space discretization, the entire soil profile is divided into 96 compartments as shown in Figure 7, in each of which the soil moisture is homogeneous. z𝑧\bigtriangleup z△ italic_z represents the size of the space step between successive nodes in discretization. The capillary pressure head information in the discretized compartments constitutes a state vector of 96 state variables. Water flow enters the soil profile at the rate of 1.944×103m/hour1.944superscript103mhour1.944\times 10^{-3}~{}\text{m}/\text{hour}1.944 × 10 start_POSTSUPERSCRIPT - 3 end_POSTSUPERSCRIPT m / hour for eight hours every day, which is the same as the irrigation schedule in Yin et al. [19] Free drainage takes place on the bottom boundary of the soil profile.

In the soil profile, there are 96 compartments of the same size after discretization, and historical data of the 96 states will be used for subsystem modeling. As depicted in the illustrative diagram of the soil profile in Figure 7, we consider that 16 permanent sensors are placed at the center points of the (12×i10)12𝑖10(12\times i-10)( 12 × italic_i - 10 )th and (12×i)12𝑖(12\times i)( 12 × italic_i )th compartments, i=1,2,,8𝑖128i=1,2,\ldots,8italic_i = 1 , 2 , … , 8. Each sensor can provide real-time measurements of the soil moisture in the corresponding compartments. The 96 state variables of the process (i.e., the capillary pressure head in the 96 compartments) are divided into 8 subsystems, each of which is assigned 12 consecutive state variables. Specifically, xi=[x{12×i11},,x{12×i}]Tsubscript𝑥𝑖superscriptsuperscript𝑥12𝑖11superscript𝑥12𝑖Tx_{i}=\big{[}x^{\left\{12\times i-11\right\}},\ldots,x^{\left\{12\times i% \right\}}\big{]}^{\text{T}}italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT = [ italic_x start_POSTSUPERSCRIPT { 12 × italic_i - 11 } end_POSTSUPERSCRIPT , … , italic_x start_POSTSUPERSCRIPT { 12 × italic_i } end_POSTSUPERSCRIPT ] start_POSTSUPERSCRIPT T end_POSTSUPERSCRIPT, i=1,2,,8𝑖128i=1,2,\ldots,8italic_i = 1 , 2 , … , 8, where x{j}superscript𝑥𝑗x^{\{j\}}italic_x start_POSTSUPERSCRIPT { italic_j } end_POSTSUPERSCRIPT is the j𝑗jitalic_jth variable of state vector x𝑥xitalic_x counting from the top to the bottom. This way, each subsystem is assigned two physical sensors for measuring the capillary pressure head in the 2nd and 12th compartment of the same subsystem.

Unknown additive process disturbances and measurement noise following Gaussian distribution are generated, and then are made bounded and added to the discretized simulation model for generating data samples. We consider synchronous sampling of the soil moisture measurements every 1 min. Then, 9,600 data samples generated from simulating (25) for an operating period of 160 hours are used for identifying the subsystem models based on the proposed method. In addition, 4,800 samples for an additional 80-hour operation are used for validating the model, and another 80-hour operating period is considered for illustrating the soil moisture estimation performance.

Subsystem modeling results

Then, the proposed parallel subsystem modeling approach in proposed in the “EDMD-based parallel Koopman subsystem identification method” subsection is applied to establish 8 subsystem models for the agro-hydrological process. State and input lifting functions ϕisuperscriptitalic-ϕ𝑖\phi^{i}italic_ϕ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT and δisuperscript𝛿𝑖\delta^{i}italic_δ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT need to be specified for subsystem i𝑖iitalic_i, i=1,2,,8𝑖128i=1,2,\ldots,8italic_i = 1 , 2 , … , 8. As the states assigned to different subsystems have the same physical meaning, the lifting functions selected for different subsystems may be made identical. Inspired by the adopted basis functions in Narasingam and Kwon [31], the original subsystem state xisubscript𝑥𝑖x_{i}italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT, quadratic functions of xisubscript𝑥𝑖x_{i}italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT, and exponential functions of xisubscript𝑥𝑖x_{i}italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT are chosen as the state lifting functions for each subsystem, that is,

ϕi(xi)=[xi,1xi,12xi,12xi,122exi,1exi,12]T,i=1,2,,8.formulae-sequencesuperscriptitalic-ϕ𝑖subscript𝑥𝑖superscriptdelimited-[]subscript𝑥𝑖1subscript𝑥𝑖12superscriptsubscript𝑥𝑖12superscriptsubscript𝑥𝑖122superscript𝑒subscript𝑥𝑖1superscript𝑒subscript𝑥𝑖12missing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionTfor-all𝑖128\phi^{i}(x_{i})=\left[{\begin{array}[]{*{20}{c}}x_{i,1}&\ldots&x_{i,12}&x_{i,1% }^{2}&\ldots&x_{i,12}^{2}&{e^{x_{i,1}}}&\ldots&{e^{x_{i,12}}}\end{array}}% \right]^{\text{T}},~{}\forall i=1,2,\ldots,8.italic_ϕ start_POSTSUPERSCRIPT italic_i end_POSTSUPERSCRIPT ( italic_x start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) = [ start_ARRAY start_ROW start_CELL italic_x start_POSTSUBSCRIPT italic_i , 1 end_POSTSUBSCRIPT end_CELL start_CELL … end_CELL start_CELL italic_x start_POSTSUBSCRIPT italic_i , 12 end_POSTSUBSCRIPT end_CELL start_CELL italic_x start_POSTSUBSCRIPT italic_i , 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT end_CELL start_CELL … end_CELL start_CELL italic_x start_POSTSUBSCRIPT italic_i , 12 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT end_CELL start_CELL italic_e start_POSTSUPERSCRIPT italic_x start_POSTSUBSCRIPT italic_i , 1 end_POSTSUBSCRIPT end_POSTSUPERSCRIPT end_CELL start_CELL … end_CELL start_CELL italic_e start_POSTSUPERSCRIPT italic_x start_POSTSUBSCRIPT italic_i , 12 end_POSTSUBSCRIPT end_POSTSUPERSCRIPT end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL end_ROW end_ARRAY ] start_POSTSUPERSCRIPT T end_POSTSUPERSCRIPT , ∀ italic_i = 1 , 2 , … , 8 .

where xi,jsubscript𝑥𝑖𝑗x_{i,j}italic_x start_POSTSUBSCRIPT italic_i , italic_j end_POSTSUBSCRIPT denotes the j𝑗jitalic_jth element of the state of subsystem i𝑖iitalic_i, i=1,2,,8for-all𝑖128\forall i=1,2,\ldots,8∀ italic_i = 1 , 2 , … , 8, j=1,2,,12for-all𝑗1212\forall j=1,2,\ldots,12∀ italic_j = 1 , 2 , … , 12. The only known input to the process is the water flow entering/leaving the soil to the process. Since this input corresponds to the surface of the soil profile, this known input is assigned to the first subsystem. Therefore, we only need to construct input lifting functions for the first subsystem, which are in the form of δ1(u1)=[u1u12eu1]Tsuperscript𝛿1subscript𝑢1superscriptdelimited-[]subscript𝑢1superscriptsubscript𝑢12superscript𝑒subscript𝑢1missing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionmissing-subexpressionT\delta^{1}(u_{1})=\big{[}{\begin{array}[]{*{20}{c}}u_{1}&u_{1}^{2}&{e^{u_{1}}}% \end{array}}\big{]}^{\text{T}}italic_δ start_POSTSUPERSCRIPT 1 end_POSTSUPERSCRIPT ( italic_u start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT ) = [ start_ARRAY start_ROW start_CELL italic_u start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT end_CELL start_CELL italic_u start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT end_CELL start_CELL italic_e start_POSTSUPERSCRIPT italic_u start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT end_POSTSUPERSCRIPT end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL start_CELL end_CELL end_ROW end_ARRAY ] start_POSTSUPERSCRIPT T end_POSTSUPERSCRIPT, which accounts for the nonlinear dependence of the lifted subsystem state of the 1st subsystem on this known input.

As has been discussed, 9,600 training samples within a 160-hour operating period are used for identifying the Koopman-based subsystem models in parallel. Snapshot matrices 𝒙𝒊subscript𝒙𝒊{\boldsymbol{x_{i}}}bold_italic_x start_POSTSUBSCRIPT bold_italic_i end_POSTSUBSCRIPT, 𝒙¯𝒊subscriptbold-¯𝒙𝒊{\boldsymbol{{\bar{x}}_{i}}}overbold_¯ start_ARG bold_italic_x end_ARG start_POSTSUBSCRIPT bold_italic_i end_POSTSUBSCRIPT, 𝒖𝒊subscript𝒖𝒊{\boldsymbol{u_{i}}}bold_italic_u start_POSTSUBSCRIPT bold_italic_i end_POSTSUBSCRIPT, and 𝒚𝒊subscript𝒚𝒊{\boldsymbol{y_{i}}}bold_italic_y start_POSTSUBSCRIPT bold_italic_i end_POSTSUBSCRIPT of appropriate dimensions are created. Based on the selected state, input lifting functions, and the snapshot matrices, Algorithm 1 is implemented to identify the models for the 8 subsystems. It is noted that the snapshot matrices are constituted by scaled states and scaled manipulated inputs that are obtained following Steps 2.1 and 2.2 of Algorithm 1. Eight Koopman subsystem models are constructed. For each subsystem model in the form of (Subsystem models for the design of local estimators), the lifted state vector zisubscript𝑧𝑖z_{i}italic_z start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT has 36 variables for i=1,2,,8𝑖128i=1,2,\ldots,8italic_i = 1 , 2 , … , 8; the lifted input vector u~1subscript~𝑢1\tilde{u}_{1}over~ start_ARG italic_u end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT for the 1st subsystem has 3 variables (since there is no known input to each subsystem j𝑗jitalic_j, j=2,3,,8𝑗238j=2,3,\ldots,8italic_j = 2 , 3 , … , 8).

Cross-validation of the established linear subsystem models based on 4,800 validation samples is performed. The trajectories of the actual soil moisture in selected compartments, and the corresponding open-loop predictions provided by the eight Koopman-based subsystem models are shown in Figures 8 and 9, where the blue lines represent the actual soil pressure head in the selected compartments, and the red dashed lines are open-loop predictions of the soil pressure head information given by the identified Koopman subsystem models. Without correction/adaptation using actual soil moisture data, the linear subsystem models can collaboratively characterize the nonlinear dynamics of the agro-hydrological system accurately. Good predictions of the soil moisture in different compartments of the soil profile are obtained. As shown in Figures 8 and 9, the discrepancy between the actual states and open-loop predictions increases with the depth. The increase in the discrepancy arises due to the error propagation from the top compartment to the current compartment. It is worth mentioning that while we only present the estimation results for selected compartments of the soil profile, the Koopman-based subsystem models can provide accurate open-loop predictions of all the states (i.e., soil moisture in all the compartments of the soil profile).

Remark 6.

In general, there is no systematic procedure that can be followed to determine the state and input lifting functions. Different lifting functions can indeed affect the accuracy of the resulting subsystem models. The selection of the lifting functions can be conducted through extensive trial-and-error tests, or neural networks that may be introduced and trained to account for the manually selected lifting functions in parametric forms; some relevant results were reported in Yeung et al. [56] and Han et al.[57]

Soil moisture estimation

In this section, a partition-based distributed MHE scheme is developed based on the Koopman subsystem models, and the estimation results for the soil moisture in different compartments, especially the compartments where no permanent sensors are deployed, are presented.

A linear distributed MHE scheme consisting of 8 local estimators is developed. Each local estimator is designed following (Formulation of the MHE-based estimators). The length of the estimation horizon for each local estimator is N=4𝑁4N=4italic_N = 4. The weighting matrices are set as follows: Pi,0=0.1×I36subscript𝑃𝑖00.1subscript𝐼36P_{i,0}=0.1\times{I}_{36}italic_P start_POSTSUBSCRIPT italic_i , 0 end_POSTSUBSCRIPT = 0.1 × italic_I start_POSTSUBSCRIPT 36 end_POSTSUBSCRIPT, Qi=0.01×I36subscript𝑄𝑖0.01subscript𝐼36Q_{i}=0.01\times{I}_{36}italic_Q start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT = 0.01 × italic_I start_POSTSUBSCRIPT 36 end_POSTSUBSCRIPT, and Ri=0.6×I6subscript𝑅𝑖0.6subscript𝐼6R_{i}=0.6\times{I}_{6}italic_R start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT = 0.6 × italic_I start_POSTSUBSCRIPT 6 end_POSTSUBSCRIPT, for i=1,2,,8𝑖128i=1,2,\ldots,8italic_i = 1 , 2 , … , 8. Note that due to the scaling, the scaled subsystem states stay within the range [0,1]01[0,1][ 0 , 1 ], which simplifies the selection of the tuning parameters of the MHE-based estimators. The initial guesses of the soil moisture estimates provided to the local estimators are made different from the actual values of the respective initial subsystem states. In the MHE-based local estimators, hard constraints are imposed on the soil moisture estimates provided by the MHE-based estimators for the subsystems, following the hard constraints considered in Yin et al. [19] and Bo at al. [54] Specifically, the estimates of the soil moisture are required to stay within [1.0,1.0×106]1.01.0superscript106[-1.0,-1.0\times 10^{-6}][ - 1.0 , - 1.0 × 10 start_POSTSUPERSCRIPT - 6 end_POSTSUPERSCRIPT ] as has been employed in Yin et al. [19]. It is noted that these constraints are imposed on the first state to the twelfth state of each local estimator that represents the soil moisture in the 12 compartments of the same subsystem. No constraints need to be imposed on the 24 lifted state variables of the subsystems.

Next, we present the soil moisture estimation results for the entire agro-hydrological process provided by the distributed MHE. Figures 10 and 11 show the trajectories of the actual soil moisture and the corresponding estimates in the selected compartments of the soil profile. The estimates can track trends of the actual soil moisture in the corresponding compartments, and the discrepancy between the actual values and the corresponding estimates remains small during the entire process operation. It is worth mentioning that despite the noticeable discrepancy between the actual states and open-loop predictions of the soil moisture in the 95th compartment (not measurable) shown in Figure 9, the proposed distributed MHE method can still provide good state estimates. This is because the developed distributed MHE is a closed-loop estimation scheme that receives and takes advantage of the information of measurable states from the local subsystem and the interacting subsystems at each time instant, whereas the Koopman-based subsystem models provide open-loop predictions of the states without using any feedback information for correction/adaptation. The trajectory of the estimation error norm is presented in Figure 12. The estimation error provided by the proposed method remains bounded during the operating period despite the use of data-driven linear subsystem models in the state estimator design.

Remark 7.

The subsystem decomposition and configuration can affect the prediction performance of the Koopman-based subsystem models. In this work, in each of the two case studies, we consider that the entire process is decomposed by taking into account the physical topology of the process. Meanwhile, in the existing literature, there have been systematic approaches to process decomposition and subsystem configuration for distributed estimation and control, including the community detection-based approaches with the aim of minimizing the strength of interactions among the configured subsystems [58, 59], the decomposition of the process into subsystems according to different time scales of the state variables [60, 61], and the hierarchical clustering [62, 63] which aim to identify the optimal subsystem designs to strengthen the connection between components within each subnetwork while weakening the connection between different subnetworks. Some of these subsystem configuration approaches may be applied to decompose a large-scale process into smaller subsystems, as the basis of subsystem modeling using the proposed approach.

Concluding Remarks

In this work, we addressed a data-driven distributed state estimation problem for general large-scale nonlinear processes that consist of interactive subsystems. The proposed framework comprises two key approaches, that is, the parallel subsystem modeling approach and the distributed state estimation approach based on the identified subsystem models. In parallel subsystem modeling, Koopman identification was leveraged as the basis of identifying linear subsystem models that interact and coordinate with each other to collaboratively approximate the nonlinear dynamics of the considered process. An implementation algorithm was proposed based on Extended Dynamic Mode Decomposition to construct the linear subsystem models in parallel. Both the states and the inputs of each subsystem are projected onto the higher-dimensional spaces. Linear subsystem models are constructed in higher-dimensional spaces. Then, we proposed a partition-based distributed state estimation approach. The established linear subsystem models are used to develop linear MHE-based state estimators, which communicate with each other and are evaluated in an iterative fashion, to provide real-time estimates of the states of the large-scale nonlinear processes. Two case studies were conducted to illustrate the proposed framework. In the first case study, a benchmark chemical process example was used to demonstrate the effectiveness and superiority of the proposed approach. In the second case study, a large-scale agro-hydrolocal process was considered. Eight linear Koopman subsystem models that interact with each other were established. The open-loop predictions given by the linear subsystem models track the actual state of the underlying nonlinear process accurately. Based on the identified subsystem models and the proposed distributed estimation algorithm, a distributed soil moisture estimation scheme was formed. Accurate soil moisture estimates were achieved across all the compartments of the soil profile, without utilizing the first-principles process model in the estimation system design.

It is noted when there is any available physical knowledge about a process, it will be favorable to incorporate the existing first-principles knowledge into Koopman modeling to build physics-enabled Koopman models, which constitutes one of potential future research directions.

Data Availability and Reproducibility Statement

The numerical data used for generating Figures 3, 4, 7-11, and Table 2 is available in the Supplementary Material. The compressed file contains data related to the case studies on the chemical reaction process and the argo-hydrological process, including data used for the identification of Koopman subsystem models, maximum and minimum values of states utilized for scaling, the identified Koopman subsystem models, and the state estimates calculated based on the proposed method. The data can be used to reproduce and generate the figures and tables presented in this paper. Particularly, the RMSE in Table 2 can be obtained by utilizing the actual states , and the state estimates provided by the proposed distributed MHE based on the linearized subsystem models and the proposed Koopman subsystem models. The subsystem matrices used for generating open-loop predictions for Figures 3, 7, and 8 can be derived by using the provided data and following Algorithm 1. The state estimates used to generate Figures 4, 9, and 10 can be obtained by implementing Algorithm 2.

Acknowledgment

This research is supported by Ministry of Education, Singapore, under its Academic Research Fund Tier 1 (RG63/22), and Nanyang Technological University, Singapore (Start-Up Grant).

Literature Cited

  • [1] Daoutidis P, Lee JH, Harjunkoski I, Skogestad S, Baldea M, Georgakis C. Integrating operations and control: A perspective and roadmap for future research. Computers & Chemical Engineering, 2018;115:179–184.
  • [2] Christofides PD, Scattolini R, de la Peña DM, Liu J. Distributed model predictive control: A tutorial review and future research directions. Computers & Chemical Engineering, 2013;51:21–41.
  • [3] Daoutidis P, Zachar M, Jogwar SS. Sustainability and process control: A survey and perspective. Journal of Process Control, 2016;44:184–206.
  • [4] Scattolini R. Architectures for distributed and hierarchical model predictive control - A review. Journal of Process Control, 2009;19:723–731.
  • [5] Chen S, Wu Z, Rincon D, Christofides PD. Machine learning-based distributed model predictive control of nonlinear processes. AIChE Journal, 2020;66(11):e17013.
  • [6] Pourkargar DB, Almansoori A, Daoutidis P. Impact of decomposition on distributed model predictive control: A process network case study. Industrial & Engineering Chemistry Research, 2017;56(34):9606–9616.
  • [7] Tang W, Daoutidis P. Coordinating distributed MPC efficiently on a plantwide scale: The Lyapunov envelope algorithm. Computers & Chemical Engineering, 2021;155:107532.
  • [8] Schneider R, Marquardt W. Convergence and stability of a constrained partition-based moving horizon estimator. IEEE Transactions on Automatic Control, 2015;61(5):1316–1321.
  • [9] Farina M, Ferrari-Trecate G, Scattolini R. Moving-horizon partition-based state estimation of large-scale systems. Automatica, 2010;46(5):910–918.
  • [10] Yin X, Zeng J, Liu J. Forming distributed state estimation network from decentralized estimators. IEEE Transactions on Control Systems Technology, 2018;27(6):2430–2443.
  • [11] ikmeşe B, Mandić M, Speyer JL. Decentralized observers with consensus filters for distributed discrete-time linear systems. Automatica, 2014;50(4):1037–1052.
  • [12] Hong Y, Chen G, Bushnell L. Distributed observers design for leader-following control of multi-agent networks. Automatica, 2008;44(3):846–850.
  • [13] Yin X, Liu J. Distributed state estimation for a class of nonlinear processes based on high-gain observers. Chemical Engineering Research and Design, 2020;160:20–30.
  • [14] Khan UA, Moura JMF. Distributing the Kalman filter for large-scale systems. IEEE Transactions on Signal Processing, 2008;56(10):4919–4935.
  • [15] Olfati-Saber R, Jalalkamali P. Coupled distributed estimation and control for mobile sensor networks. IEEE Transactions on Automatic Control, 2012;57(10):2609–2614.
  • [16] Roshany-Yamchi S, Cychowski M, Negenborn RR, De Schutter B, Delaney K, Connell J. Kalman filter-based distributed predictive control of large-scale multi-rate systems: Application to power networks. IEEE Transactions on Control Systems Technology, 2013;21(1):27–39.
  • [17] Vadigepalli R, Doyle FJ. A distributed state estimation and control algorithm for plantwide processes. IEEE Transactions on Control Systems Technology, 2003;11:119–127.
  • [18] Battistelli G, Chisci L. Stability of consensus extended Kalman filter for distributed state estimation. Automatica, 2016;68:169–178.
  • [19] Yin X, Bo S, Liu J, Huang B. Consensus-based approach for parameter and state estimation of agro-hydrological systems. AIChE Journal, 2021;67(2):e17096.
  • [20] Farina M, Ferrari-Trecate G, Scattolini R. Distributed moving horizon estimation for nonlinear constrained systems. International Journal of Robust and Nonlinear Control, 2012;22:123–143.
  • [21] Pourkargar DB, Moharir M, Almansoori A, Daoutidis P. Distributed estimation and nonlinear model predictive control using community detection. Industrial & Engineering Chemistry Research, 2019;58(30):13495–13507.
  • [22] Zhang J, Liu J. Distributed moving horizon estimation for nonlinear systems with bounded uncertainties. Journal of Process Control, 2013;23(9):1281–1295.
  • [23] Yin X, Decardi-Nelson B, Liu J. Subsystem decomposition and distributed moving horizon estimation of wastewater treatment plants. Chemical Engineering Research and Design, 2018;134:405–419.
  • [24] Koopman BO. Hamiltonian systems and transformation in Hilbert space. Proceedings of the National Academy of Sciences of the United States of America, 1931;17(5):315.
  • [25] Mezić I, Banaszuk A. Comparison of systems with complex behavior. Physica D: Nonlinear Phenomena, 2004;197(1-2):101–133.
  • [26] Bruder D, Fu X, Gillespie RB, Remy CD, Vasudevan R. Data-driven control of soft robots using Koopman operator theory. IEEE Transactions on Robotics, 2021;37(3):948–961.
  • [27] Schmid PJ. Dynamic mode decomposition of numerical and experimental data. Journal of Fluid Mechanics, 2010;656:5–28.
  • [28] Williams MO, Kevrekidis IG, Rowley CW. A data–driven approximation of the Koopman operator: Extending dynamic mode decomposition. Journal of Nonlinear Science, 2015;25(6):1307–1346.
  • [29] Mauroy A, Mezić I, Moehlis J. Isostables, isochrons, and Koopman spectrum for the action–angle representation of stable fixed point dynamics. Physica D: Nonlinear Phenomena, 2013;261:19–30.
  • [30] Froyland G, Gottwald GA, Hammerlindl A. A computational method to extract macroscopic variables and their dynamics in multiscale systems. SIAM Journal on Applied Dynamical Systems, 2014;13(4):1816–1846.
  • [31] Narasingam A, Kwon JSI. Koopman Lyapunov-based model predictive control of nonlinear chemical process systems. AIChE Journal, 2019;65(11):e16743.
  • [32] Korda M, Mezić I. Linear predictors for nonlinear dynamical systems: Koopman operator meets model predictive control. Automatica, 2018;93:149–160.
  • [33] Son SH, Choi HK, Moon J, Kwon JSI. Hybrid Koopman model predictive control of nonlinear systems using multiple EDMD models: An application to a batch pulp digester with feed fluctuation. Control Engineering Practice, 2022;118:104956.
  • [34] Narasingam A, Kwon JSI. Application of Koopman operator for model-based control of fracture propagation and proppant transport in hydraulic fracturing operation. Journal of Process Control, 2020;91:25–36.
  • [35] Kaiser E, Kutz JN, Brunton SL. Data-driven discovery of Koopman eigenfunctions for control. Machine Learning: Science and Technology, 2021;2(3):035023.
  • [36] Netto M, Mili L. A robust data-driven Koopman Kalman filter for power systems dynamic state estimation. IEEE Transactions on Power Systems, 2018;33(6):7228–7237.
  • [37] Yin X, Qin Y, Liu J, Huang B. Data-driven moving horizon state estimation of nonlinear processes using Koopman operator. In press, doi:10.1016/j.cherd.2023.10.033.
  • [38] Proctor JL, Brunton SL, Kutz JN. Generalizing Koopman theory to allow for inputs and control. SIAM Journal on Applied Dynamical Systems, 2018;17(1):909–930.
  • [39] Brunton SL, Brunton BW, Proctor JL, Kutz JN. Koopman invariant subspaces and finite linear representations of nonlinear dynamical systems for control. PloS one, 2016;11(2):e0150171.
  • [40] Peitz S, Otto SE, Rowley CW. Data-driven model predictive control using interpolated Koopman generators. SIAM Journal on Applied Dynamical Systems, 2020;19(3):2162–2193.
  • [41] Narasingam A, Kwon JSI. Closed-loop stabilization of nonlinear systems using Koopman Lyapunov-based model predictive control. In Conference on Decision and Control (CDC), 2020;704–709, Jeju Island, South Korea.
  • [42] Han M, Li Z, Yin X, Yin X. Robust learning and control of time-delay nonlinear systems with deep recurrent Koopman operators. IEEE Transactions on Industrial Informatics, in press, doi: 10.1109/TII.2023.3328432.
  • [43] Zhang X, Han M, Yin X. Reduced-order Koopman modeling and predictive control of nonlinear processes. Computers & Chemical Engineering, 2023;179:108440.
  • [44] Surana A, Williams MO, Morari M, Banaszuk A. Koopman operator framework for constrained state estimation. In Conference on Decision and Control (CDC), 2017;94–101, Melbourne, VIC, Australia.
  • [45] Arbabi H, Korda M, Mezić I. A data-driven Koopman model predictive control framework for nonlinear partial differential equations. In Conference on Decision and Control (CDC), 2018;6409–6414, Miami Beach, FL, USA.
  • [46] Budišić M, Mohr R, Mezić I. Applied Koopmanism. Chaos: An Interdisciplinary Journal of Nonlinear Science, 2012;22(4):047510.
  • [47] Rao CV, Rawlings JB, Lee JH. Constrained linear state estimation–a moving horizon approach. Automatica, 2001;37(10):1619–1628.
  • [48] Li X, Law AWK, Yin X. Partition-based distributed extended Kalman filter for large-scale nonlinear processes with application to chemical and wastewater treatment processes. AIChE Journal, 2023;e18229.
  • [49] Li X, Bo S, Qin Y, Yin X. Iterative distributed moving horizon estimation of linear systems with penalties on both system disturbances and noise. Chemical Engineering Research and Design, 2023;194:878-893.
  • [50] Rashedi M, Liu J, Huang B. Triggered communication in distributed adaptive high-gain EKF. IEEE Transactions on Industrial Informatics, 2017;14(1):58–68.
  • [51] Richards LA. Capillary conduction of liquids through porous mediums. Physics, 1931;1(5):318–333.
  • [52] Mualem Y. A new model for predicting the hydraulic conductivity of unsaturated porous media. Water Resources Research, 1976;12(3):513–522.
  • [53] Van Genuchten MT. A closed-form equation for predicting the hydraulic conductivity of unsaturated soils. Soil Science Society of America Journal, 1980;44(5):892–898.
  • [54] Bo S, Sahoo SR, Yin X, Liu J, Shah SL. Parameter and state estimation of one-dimensional infiltration processes: A simultaneous approach. Mathematics, 2020;8(1):134.
  • [55] Carsel RF, Parrish RS. Developing joint probability distributions of soil water retention characteristics. Water resources research, 1988;24(5):755–769.
  • [56] Yeung E, Kundu S, Hodas N. Learning deep neural network representations for Koopman operators of nonlinear dynamical systems. In American Control Conference, 2019;4832–4839, Philadelphia, PA, USA.
  • [57] Han Y, Hao W, Vaidya U. Deep learning of Koopman representation for control. In Conference on Decision and Control (CDC), 2020;1890–1895, Jeju Island, South Korea.
  • [58] Yin X, Liu J. Subsystem decomposition of process networks for simultaneous distributed state estimation and control. AIChE Journal, 2019;65(3):904-914.
  • [59] Tang W, Allman A, Pourkargar DB, Daoutidis P. Optimal decomposition for distributed optimization in nonlinear model predictive control through community detection. Computers & Chemical Engineering, 2018;111:43-54.
  • [60] Yin X, Liu J. Distributed moving horizon state estimation of two-time-scale nonlinear systems. Automatica, 2017;79:152-161.
  • [61] Wu L, Yin X, Pan L, Liu J. Economic model predictive control of integrated energy systems: A multi-time-scale framework. Applied Energy, 2022;328:120187.
  • [62] Kang L, Tang T, Liu Y, Daoutidis P. Control configuration synthesis using agglomerative hierarchical clustering: A graph-theoretic approach. Journal of Process Control, 2016;46:43-54.
  • [63] Heo S, Daoutidis P. Control-relevant decomposition of process networks via optimization-based hierarchical clustering. AIChE Journal, 2016;62(9):3177-3188.

Refer to caption

Figure 1: An illustrative diagram depicting the connection between the subsystem modeling and local estimator design of the proposed distributed estimation framework.
Refer to caption
Figure 2: A schematic of the distributed state estimation scheme developed using Koopman subsystem models.
Refer to caption
Figure 3: A schematic of the four-CSTR process.

Refer to caption

Figure 4: Cross-validation of the Koopman subsystem models for the four-CSTR process.

Refer to caption

Figure 5: Trajectories of the actual system states and state estimates for four vessel CSTRs.
Refer to caption
Figure 6: An illustrative schematic of the agro-hydrological process.

Refer to caption

Figure 7: The compartments of the soil profile after vertical space discretization, and the locations of the sensors (black dots).

Refer to caption

Figure 8: Cross-validation of the Koopman subsystem models (1st subsystem to 4th subsystem) for the argo-hydrological process.

Refer to caption

Figure 9: Cross-validation of the Koopman subsystem models (5th subsystem to 8th subsystem) for the argo-hydrological process.

Refer to caption

Figure 10: Trajectories of the actual system states and state estimates for the argo-hydrological process (1st subsystem to 4th subsystem).

Refer to caption

Figure 11: Trajectories of the actual system states and state estimates for the argo-hydrological process (5th subsystem to 8th subsystem).

Refer to caption

Figure 12: Trajectories of the estimation error norm for the argo-hydrological process.
Table 1: The initial state x0subscript𝑥0x_{0}italic_x start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT and the initial guess x¯0subscript¯𝑥0\bar{x}_{0}over¯ start_ARG italic_x end_ARG start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT for the process
T1subscript𝑇1T_{1}italic_T start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT CA1subscript𝐶𝐴1C_{A1}italic_C start_POSTSUBSCRIPT italic_A 1 end_POSTSUBSCRIPT T1subscript𝑇1T_{1}italic_T start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT CA2subscript𝐶𝐴2C_{A2}italic_C start_POSTSUBSCRIPT italic_A 2 end_POSTSUBSCRIPT T3subscript𝑇3T_{3}italic_T start_POSTSUBSCRIPT 3 end_POSTSUBSCRIPT CA3subscript𝐶𝐴3C_{A3}italic_C start_POSTSUBSCRIPT italic_A 3 end_POSTSUBSCRIPT T4subscript𝑇4T_{4}italic_T start_POSTSUBSCRIPT 4 end_POSTSUBSCRIPT CA4subscript𝐶𝐴4C_{A4}italic_C start_POSTSUBSCRIPT italic_A 4 end_POSTSUBSCRIPT
x0subscript𝑥0x_{0}italic_x start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT 310.9387 3.0317 310.8962 2.8002 312.5167 2.8441 311.2047 3.0141
x¯0subscript¯𝑥0\bar{x}_{0}over¯ start_ARG italic_x end_ARG start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT 311.0766 3.0318 311.1287 2.8003 312.7482 2.8440 311.5002 3.0139
Table 2: Comparison of computation time and estimation norm
Estimation method Computation time (s) RMSE
Distributed MHE based on Koopman subsystem models 0.1858 0.0135
Distributed MHE based on linearized subsystem models 0.0653 1.5204