Journal of Modern Power Systems and Clean Energy

ISSN 2196-5625 CN 32-1884/TK

网刊加载中。。。

使用Chrome浏览器效果最佳,继续浏览,你可能不会看到最佳的展示效果,

确定继续浏览么?

复制成功,请在其他浏览器进行阅读

Distributed Real-time State Estimation for Combined Heat and Power Systems  PDF

  • Tingting Zhang
  • Wen Zhang
  • Qi Zhao
  • Yaxin Du
  • Jian Chen
  • Junbo Zhao
Key Laboratory of Power System Intelligent Dispatch and Control, Ministry of Education, Shandong University, Jinan 250061, China; Department of Electrical and Computer Engineering, Mississippi State University, Starkville, MS 39762, USA

Updated:2021-03-16

DOI:10.35833/MPCE.2020.000052

  • Full Text
  • Figs & Tabs
  • References
  • Authors
  • About
OUTLINE

Abstract

This paper proposes a distributed real-time state estimation (RTSE) method for the combined heat and power systems (CHPSs). First, a difference-based model for the heat system is established considering the dynamics of heat systems. This heat system model is further used along with the power system steady-state model for holistic CHPS state estimation. A cubature Kalman filter (CKF)-based RTSE is developed to deal with the system nonlinearity while integrating both the historical and present measurement information. Finally, a multi-time-scale asynchronous distributed computation scheme is designed to enhance the scalability of the proposed method for large-scale systems. This distributed implementation requires only a small amount of information exchange and thus protects the privacy of different energy systems. Simulations carried out on two CHPSs show that the proposed method can significantly improve the estimation efficiency of CHPS without loss of accuracy compared with other existing models and methods.

I. Introduction

CONSIDERING the growing energy crisis and the environmental pollution, the integrated energy system (IES), including electricity, heat, cold, gas, and other complementary energy forms, has gained more and more attention [

1]. IES can make full use of the flexibility of various energy systems to promote energy efficiency and reduce carbon emissions [2]. The combined heat and power system (CHPS) is the most widely-used IES, which consists of the power system, the heat system, and energy conversion units such as the combined heat and power (CHP) units. CHPS can not only generate both power and heat, but also collect the waste heat to achieve higher energy efficiency [3].

The online control and management of CHPS necessitate the need to gain access to reliable and real-time system operation states. However, due to the lack of real-time measurements and the information barrier between different energy systems, the measurement information collected from the metering devices is incomplete and also suffers from unavoidable noises. Therefore, there is an essential need to develop state estimation for the CHPS for real-time monitoring. Although the state estimation for power system has been extensively investigated, there are few works on CHPS. In [

4], a state estimation model for CHPS considering the operation characteristics of coupling components is first established. A distributed bilinear state estimation method for the IES is proposed in [5].

In practical heat systems, the temperature fluctuation at the inlet of a pipeline slowly spreads to the outlet, causing a time delay in heat transfer [

6]. Such heat dynamics couple the current nodal temperature with the historical ones and the mass flow rates. In this condition, the static model, which treats the heat transfer as a quasi-steady model, will lead to poor estimation accuracy. Hence, the dynamic model of the heat system should be considered in CHPS state estimation. In [7], a node method is proposed for dynamic modeling of the heat system to incorporate the discrete spatio-temporal constraints. However, this backward induction method is quite complicated. To address this issue, a novel difference-based model is developed in this paper to capture the dynamics of heat systems, where only the last historical state snapshot is needed.

Static state estimators may obtain inaccurate state awareness if the heat dynamics are ignored. However, few studies focus on the dynamic state estimation (DSE) for CHPS. In [

8], a two-stage state estimation approach for CHPS is employed considering the time delay of heat transfer. In [9], [10], an alternating estimation strategy for CHPS is proposed to handle the complicated heat dynamic constraints of temperature. These methods incorporate the discrete spatio-temporal constraints into the weighted least squares (WLS) [11], yielding the so-called quasi-dynamic state estimations (QDSEs). For power systems, a static state estimation is performed without state transition matrix. Unlike the quasi-dynamic ones, the nonlinear estimators within the Kalman filter (KF) framework can utilize both the previous states and the current measurements to recursively estimate the system states [12]. The extended KF (EKF) is a widely-used method [13]. Nevertheless, it may have large linearization errors for highly nonlinear systems. The unscented KF (UKF) employs the unscented transformation to approximate the Gaussian distributions of the nonlinear functions without linearization [14]. However, it may suffer from numerical stability issues when applied to high-dimension systems due to the nonpositive definite covariance matrix in case of inappropriate parameters [15]. To address that, the cubature KF (CKF) [16] is adopted in this paper.

Compared with the individual estimations for different energy systems, combined state estimation for CHPS can obtain coordinated estimates, avoiding data mismatches on coupling units. However, due to the independent operation of different energy systems and the requirement for privacy protection, only limited data can be exchanged [

17]. Besides, since the CHPS is typical of large scale, a centralized combined state estimation (CCSE) is too computationally demanding for practical application. Furthermore, the difference of time resolution between the two energy systems also challenges the implementation of the CCSE. Therefore, there is an urgent need to design a distributed scheme for CHPS state estimation.

This paper proposes a distributed CKF-based real-time state estimation (CKF-RTSE) method for CHPS. Firstly, a difference model considering the dynamic characteristics of the heat system is established based on heat conservation. Then, a CKF-RTSE method is proposed for the CHPS. To deal with the time resolution issue, a multi-time-scale asynchronous distributed scheme is proposed. This allows performing real-time state estimations (RTSEs) for the power and heat systems at different frequencies during the system operation. These RTSEs are implemented in an asynchronous distributed manner when they need to be executed simultaneously. The main contributions of this paper are summarized as follows:

1) A difference-based model for the heat system in CHPS is developed considering both the steady and dynamic characteristics of the heat system, where only the last historical state snapshot is required to predict the current temperatures.

2) A CKF-RTSE method for CHPS is proposed by incorporating the proposed CHPS models into a CKF algorithm. It allows tracking CHPS states and its heat dynamic process, which facilitates real-time control and management.

3) A multi-time-scale asynchronous distributed scheme is developed for RTSE implementation, which enhances the computation efficiency and protects the privacy of different energy systems while obtaining the coordinated operation states.

The rest of the paper is organized as follows. The CHPS model considering the thermal dynamics of the heat system is established in Section II. The CKF-RTSE for CHPS is proposed in Section III. The multi-time-scale asynchronous distributed scheme of the RTSE is proposed in Section IV. Simulation results are shown in Section V. Finally, Section VI concludes the paper.

II. CHPS Models Considering Thermal Dynamics

A CHPS consists of the power system, the heat system, and their coupling units. Since the time constant of the heat system is much larger than that of the power system, the heat system has much slower dynamics. Considering the requirement of the online control and management of CHPS, although the measurements can be collected at a high frequency, the state estimation for CHPS is preferred for the implementation at a minute level. Therefore, the quasi-steady model is assumed for power system to capture bus voltage magnitude and angle variations while the dynamic model for the heat system is developed.

A. Models for Heat Systems

In this paper, the heat systems are described by a set of differential and algebraic equations. In particular, the continuity of water flow, loop pressure, heat power, temperature mixture, and heat loss [

8], [18], [19] are described by algebraic equations; while the thermal dynamics are described by a set of differential equations.

1) Algebraic Equations for Heat System

1) Continuity of water flow

The mass flow of the water is governed by the rule that the mass flow of water entering into a bus equals the sum of the mass flow of water flowing out of the bus and the water consumption of the bus, which can be expressed as:

Am=mq (1)

where A is the network incidence matrix that indicates the incidence relationships among the buses and branches in the heat system; m is the vector of mass flow inside each pipeline; and mq is the vector of mass flow injected to the heat load at each bus.

2) Calculation of loop pressure

Due to the friction inside the pipelines, the pressure head of water drops as it flows, and the difference in pressure head between buses will drive the flow of water in turn. Similar to the Kirchhoff’s voltage law, the head loss around a closed hydraulic loop must be zero, which can be described by:

BKfricm*m=0 (2)

where B is the loop incidence matrix that shows the relationships between the closed circles and each pipeline; Kfric is the resistance coefficient of each pipeline; m=mi is used to denote the absolute value of the matrix; and * is the Hadamard product.

3) Calculation of heat power

The heat power generated by a heat source or consumed by a heat load can be calculated via:

Φ=CpmqTs-Tr (3)

where Φ is the heat power at each bus; Cp is the specific heat of water; mq is a single variable of mq; and Ts and Tr are the nodal supply and return temperatures, respectively.

4) Calculation of heat loss

Considering the heat loss, the relationship between the temperature at the beginning and the end of the pipeline can be expressed as:

Tend=Tstart-Tae-λLCpm+Ta (4)

where Tend is the temperature at the end of the pipeline before mixing; Tstart is the temperature at the beginning of the pipeline; Ta is the ambient temperature; m is a single variable of m; λ is the heat transfer coefficient per unit length; and L is the length of the pipeline.

5) Mixture of nodal temperature

The temperature of water flows out of a mixing bus can be calculated through:

moutTout=minTin (5)

where mout and min are the mass flow rates leaving and entering the bus, respectively; Tout is the temperature of the bus after the mixture; and Tin is the water temperature at the end of an incoming pipeline.

2) Differential Equations for Thermal Dynamics Model

The heat energy is transmitted by hot water inside pipelines slowly. The change of the temperature at the inlet of a pipeline will affect the outlet temperature after the time delay of water transport. Considering that most heat systems operate at the constant flow and variable temperature (CF-VT) in practice, the hydraulic dynamics can be ignored and the thermal dynamics will be taken into consideration.

The thermal dynamics of the heat system are mainly reflected in the heat transport delay and the heat loss in the transition process. Based on the heat conservation law, the partial differential equation can be expressed as [

19]:

Tt+mρSTx+λCpρST-Ta=0 (6)

where T is the temperature at any point in a pipeline; S is the cross-sectional area of a pipeline; ρ is the density of water; and t and x are the time and the length from the pipeline inlet, respectively.

By discretization based on the Lax-Wendroff method [

20], (6) can be rewritten as:

Ti,k+1+Tj,k+1-Ti,k-Tj,k+mijΔtρSijLijTj,k+1+Tj,k-Ti,k+1-Ti,k+λΔt2CpρSijTi,k+Tj,k+Ti,k+1+Tj,k+1-4Ta=0 (7)

where k is the kth time instant; and Ti,k and Tj,k are the temperatures at the start and end of the pipeline ij at the kth time instant, respectively.

Compared with the node-based model in [

7], the proposed difference-based model only needs the last historical state snapshot of nodal supply and returns the information of temperatures. Thus, only a small amount of historical information is needed, simplifying the state prediction with dynamic characteristics of the heat system considered.

B. Models for Power System

The AC power flow model is utilized, which can be expressed as:

Pi=VijVjgijcosθij+bijsinθijQi=VijVjgijsinθij-bijcosθij (8)

where Pi and Qi are the active and reactive power at bus i, respectively; Vi is the voltage magnitude of bus i; gij and bij are the conductance and the susceptance of the branch between bus i and bus j, respectively; and θij is the difference between the phase angles of bus i and bus j.

C. Models for CHP Units

The CHP units can generate power and heat simultaneously. According to the internal mechanism, the CHP units can be divided into gas turbines, internal combustion reciprocating engines, and steam turbines. In this paper, the gas turbines and steam turbines are taken as typical CHP units.

For gas turbines, the relationship between their power and heat output is shown as:

cm=ΦCHPPCHP (9)

where ΦCHP is the heat output of the CHP unit; PCHP is the power output; and cm is the heat-to-power ratio.

For steam turbines, if the fuel consumption is constant, the increase of the power output will result in a decrease of heat output. Under this condition, the relationship between the power and heat output can be expressed as:

ηs=ΔΦΔP=ΦCHPPcap-PCHP (10)
Pcap=ηeGin (11)

where ηs is the ratio that describes the heat and power output of the CHP unit; Pcap is the maximum power output with zero heat output; ηe is the efficiency of the CHP unit; and Gin is the fuel input.

III. CKF-RTSE for CHPS

In this section, the CKF-RTSE method is proposed for CHPS considering the heat dynamics. First, the RTSE problem for CHPS is formulated. Then, the CKF-RTSE method is proposed.

A. Problem Formulation of RTSE for CHPS

Although an increasing number of real-time measurements are available such as phasor measurement units (PMUs), there is still a lack of real-time measurements in power systems. In this condition, the power system state estimation is still nonlinear due to the utilization of supervisory control and data acquisition (SCADA) and pseudo-measurements. Hence, a nonlinear filtering method is needed to provide accurate state awareness for the CHPS.

The conventional state estimation based on WLS extracts the current states from a single measurement snapshot. It may get trapped in local minimum if unsuitable initial states are utilized when severe power fluctuations occur. In this case, the WLS-based state estimation will suffer from low accuracy and efficiency, thus being inappropriate for CHPS real-time monitoring. To address this issue, the nonlinear estimators developed within the KF framework are preferred, since they can leverage both previous states and current measurements to recursively estimate states, yielding higher computation efficiency without loss of accuracy.

Though the measurements can be collected at a high frequency, the state estimator of CHPS will suffer a great deal of unnecessary communication and computation overhead if executed at such a frequency. In this paper, the RTSE for the CHPS aims to provide reliable and real-time state awareness in normal cases for advanced applications such as contingency analysis. Under this condition, the intervals of state estimation for both the power and heat systems should be consistent with those of system operation (5-15 min), and the measurements can be updated only when the state estimations of the power and heat systems are carried out. Therefore, the RTSE for CHPS is preferred to be executed every few minutes. In this way, the operation state of CHPS can be tracked well while avoiding the waste of computation resources.

Note that due to the distinction in time scale between the different energy systems, the transient process of the power system is neglected, but the slow dynamics of the heat system will be tracked. In other words, the dynamics of the power system is regarded as a quasi-steady behavior while the slower thermal dynamics of the heat system is considered. Therefore, a forecasting-aided state estimation (FASE) will be performed in the power system and a DSE will be developed for the heat system.

Generally, the state-space model for a CHPS, both for the FASE and DSE, can be formulated by a set of discrete-time nonlinear functions as [

12]:

x˜k=f(xk-1,uk,ωk)zk=h(x˜k,uk,ξk) (12)

where xk=[xe,k,xh,k]T and uk=[ue,k,uh,k]T are the state and the input vectors at time instant k, respectively, and the subscripts e and h denote the power and the heat systems, respectively; ωk and ξk are the white noise vectors of system process and measurements, which are assumed to be white Gaussian in this paper, and their covariance matrices are denoted by Qk and Rk, respectively; f() is a set of nonlinear functions relevant to the predicted states x˜k; and h() is a set of vector-valued nonlinear functions relating xk to the measurement vector zk=[ze,k,zh,k]T.

1) FASE for Power System

For the power system, its state vector xe includes the nodal voltage magnitudes V and phase angles θ, i.e., xe=[V,θ]T. The input vector ue is composed of the nodal active and reactive power injections. ze denotes the measurement vector of the power system, including the real-time measurements provided by SCADAs and PMUs, and pseudo-measurements.

As mentioned above, the fast dynamics of the power system are neglected. Thus, the state transition in power systems is treated as a quasi-steady behavior. In this paper, the state transition of the power system is formulated by a load-flow-based state prediction, which is expressed as:

xe,k+1=xe,k-Je-1Δue (13)

where Je is the Jacobian matrix of the power system; and Δue is the change of nodal loads obtained by load forecasting.

The real-time measurements consist of the nodal voltage magnitude measurements, the branch current measurements, and the active as well as reactive power flows on the branches from SCADA. The voltage magnitudes and phase angles of the buses equipped with PMUs are also used as real-time measurements. All the nodal active and reactive power injections are leveraged as pseudo-measurements. Hence, the measurement functions can be expressed as follows.

1) Real-time measurements

Pij,meas=ViVjgijcosθij+bijsinθij-Vi2gijQij,meas=ViVjgijsinθij-bijcosθij+Vi2bijVi,meas=Viθi,meas=θiIij,meas=Vi-VjZij (14)

2) Pseudo-measurements

Pi,meas=VijVjgijcosθij+bijsinθijQi,meas=VijVjgijsinθij-bijcosθij (15)

where the subscript meas refers to the measurements employed in the state estimation; θi is the phase angle of bus i with Vi as the corresponding complex voltage phasor; Pij and Qij are the active and reactive power flows of branch ij, respectively; Iij is the magnitude of current of branch ij; and Zij is the impedance of the branch between bus i and bus j.

2) DSE for Heat System

As for the heat system, the state variables include the supply and return temperatures at each bus denoted by xh=[Ts,Tr]T. The measurements are the heat energy consumed at each bus as well as the supply and return temperatures, i.e., zh=[Φ,Ts,Tr]T. Both the state transition and the measurement equations can be derived based on the model of the heat system shown in Section II-A.

B. CKF-RTSE

The quality of the approximation of the nonlinear functions has a great impact on the performance of nonlinear estimators. A poor approximation will lead to inaccurate estimation. Due to the nonlinearity of the CHPS, especially that of the power systems, an effective nonlinear estimator is needed.

The CKF is a useful recursive estimation method to address this issue. Similar to the UKF, the CKF also employs the sampling technique to approximate the Gaussian distributions of the nonlinear functions [

15], [16], [21]. In this way, the higher order information can be used, which enhances the estimation precision.

The CKF involves three main steps: cubature point calculation, time update, and measurement update, which are described as follows [

21].

1) Cubature Point Calculation

The cubature point calculation is performed based on the spherical-radical rule. A set of 2n cubature points are generated with corresponding weights wi to capture the statistics, i.e., the mean and covariance, of the previous state estimate xk according to:

Pk=SkSkT (16)
χkk(i)=Skξi+xk    i=1,2,,2n (17)
ξi=nei       i=1,2,,n-nei    i=n+1,n+2,,2n (18)
wic=wim=12n (19)

where Sk is the square-root matrix of the estimation error covariance Pk at time instant k; χkk(i) is the ith cubature point; ξi is the basic data point set; ei is the ith column of a unity matrix; and wm and wc are the weights of each cubature point for the mean and the covariance, respectively. Since the weight wc of each cubature point keeps positive, the positive definitiveness of the covariance matrix can be guaranteed.

2) Time Update

Based on the state transition model in (12), the set of cubature points can be propagated via the nonlinear transition function, which yields:

χk+1k(i)=f(χkk(i),k)    i=1,2,,2n (20)

Then the predicted state x˜k+1k, i.e., x˜k in (12), and its covariance Pk+1k can be calculated by:

x˜k+1k=i=12nwimχk+1k(i) (21)
Pk+1k=i=12nwicχk+1k(i)(χk+1k(i))T-x˜k+1kx˜k+1kT+Qk-1 (22)

3) Measurement Update

With the predicted state x˜k+1k and its covariance Pk+1k, another set of χk+1k_(i) can be generated similarly:

Pk+1k=Sk+1kSk+1kT (23)
χk+1k_(i)=Sk+1kξi+x˜k+1k    i=1,2,,2n (24)

These cubature points are then propagated through the measurement model in (12), which leads to:

zk+1k(i)=h(χk+1k_(i))    i=1,2,,2n (25)

The mean of the propagated cubature points, the measurement covariance, and the cross covariance of the state and measurement can be expressed as:

z¯k+1k=i=12nwimzk+1k(i) (26)
Pk+1kzz=i=12nwiczk+1k(i)(zk+1k(i))T-z¯k+1k(z¯k+1k(i))T+Rk+1 (27)
Pk+1kxz=i=12nwicχk+1k_(i)(zk+1k(i))T-x˜k+1kz¯k+1kT (28)

where the superscript zz indicates the measurement covariance; and the superscript xz indicates the cross covariance of the state and measurement.

Thereafter, the Kalman gain can be calculated as:

Kk+1=Pk+1kxz(Pk+1kzz)-1 (29)

Finally, the estimated state x̂k+1k+1 and its covariance matrix Pk+1k+1 can be calculated as:

x̂k+1k+1=x˜k+1k+Kk+1(zk-z¯k+1k) (30)
Pk+1k+1=Pk+1k-Kk+1Pk+1kzzKk+1T (31)

IV. Distributed RTSE for Multi-time-scale CHPS

In this section, a multi-time-scale asynchronous distributed scheme is designed for the CKF-based RTSE implementation to enhance the computation efficiency. First, the background of the distributed RTSE for CHPS is introduced. Second, a multi-frequency state estimation framework is established considering the different time scales of the diverse energy systems. Third, an asynchronous parallel computation scheme is designed for the time instant when the RTSEs of the power and the heat systems are performed simultaneously. Finally, the overall distributed scheme is presented.

A. Background of Distributed RTSE for CHPS

To monitor the CHPS accurately, a CCSE should be performed to obtain the coordinated operation states. However, it confronts three main challenges.

1) Since CHPS is typical of large scale, a CCSE is very time-consuming and inapplicable for real-time monitoring.

2) As the time scales of the diverse energy systems are quite different, a CCSE may suffer from unnecessary computation overhead or failure in tracking the CHPS dynamics in case of inappropriate time interval.

3) The privacy concern of the power and the heat systems may challenge the implementation of a CCSE for CHPS. Since the two systems are operated independently by different utilities, it is quite difficult to obtain the measurement information from the other system.

To overcome these challenges, a multi-frequency asynchronous parallel distributed scheme is proposed for the RTSE of the CHPS based on the proposed decomposition and coordination method. It should be noted that the multi-frequency estimation scheme is carried out all the time for the operation of the CHPS; while the asynchronous parallel scheme is employed when the state estimations for both the power and the heat systems are performed simultaneously at one time instant. The detailed multi-frequency estimation scheme and the asynchronous parallel scheme are introduced in the following two subsections.

The system partition is a prerequisite for a distributed state estimation. In this paper, the whole CHPS is separated into a power system and a heat system, where the two systems are overlapped by their coupling CHPs, as shown in Fig. 1. However, a deeper network partition of the power or the heat system is not considered in this paper. When the scales of the two systems increase, they can be further decomposed into smaller regions based on the network partition method proposed in [

22].

Fig. 1 System partition scheme for CHPS.

B. Multi-frequency State Estimation Framework

As mentioned above, the dynamics of the heat system are much slower than those of the power system. When executing the RTSE of CHPS at a unified frequency, it may fail to track the states accurately if the execution frequency is too low, or may suffer from a communication bottleneck with high execution frequency. Hence, the RTSE of these two different systems should be performed at different frequencies.

Since the time scale of the power system is much shorter than that of the heat system, the FASE of the power system is required to be carried out more frequently than the heat system. As shown in Fig. 2, the DSE of the heat system is performed once every few implementations of the FASE of the power system. In this paper, considering the real-time monitoring of the CHPS in normal cases, the FASE for the power system will be performed every 5 min, which is consistent with the time interval of the online dispatch and control of the power system; while the DSE for the heat system is implemented every 15 min. When it is time to perform the state estimations for the two systems simultaneously (t=3KΔT, KN+), an asynchronous parallel scheme is adopted, as described in Section IV-C, during which the states of the heat system are estimated and updated. Otherwise, the FASE for the power system will be carried out individually, where the states of the heat system are kept consistent with their last updated values.

Fig. 2 Multi-frequency state estimation framework.

The implementation of this multi-frequency state estimation framework requires that the state estimations for the power and the heat systems are synchronized every 3ΔT. To achieve this, a buffer is employed at each local estimator for the two systems. Besides, a short latency time is set to wait for the measurements due to their transmission delay. For the measurements with time stamp such as the PMU measurements, they can be easily aligned according to their time stamps. For the measurements without time stamp such as the SCADA measurements received during the latency time, they will be regarded as the ones updated at the same time instant, which will then be used for the estimation together with those with time stamp. For future work, the event-trigger strategy [

23] will be investigated to achieve the synchronous aligning, where only the measurements containing innovational information will be updated to the estimators.

C. Asynchronous Parallel Computation Scheme

When performing the state estimations of the power system and the heat system simultaneously, as shown in Fig. 2, the constraints of CHPs’ power and heat generation (9) or (10) should be satisfied, calling for data communication between local estimators to avoid data mismatch. However, since the power system adopts a nonlinear model while the proposed model for the heat system is a linear one, there may be a difference between the execution time for their individual estimations. If the estimations are performed synchronously, a large amount of time will be wasted to wait for the updates from each other. Therefore, an asynchronous parallel computation scheme is proposed in this paper, which allows to reduce the communication burden significantly.

When performing the local state estimations simultaneously, the local CKF-RTSEs are performed individually by using the measurements within their systems. When the CKF-RTSE in one system has converged, the power or heat outputs of the CHPs, P̂CHP,e or Φ̂CHP,h, will be calculated using their state estimates. Thereafter, the outputs will be transformed into the equivalent ones in the other energy form, Φ̂CHP,e or P̂CHP,h, and then transferred to the estimators in the other system as augmented measurements. This can be used for another iteration of local estimations. The global convergence criterion for the two systems can be expressed as:

maxP̂CHP,e-P̂CHP,h<εe (32)
maxΦ̂CHP,h-Φ̂CHP,e<εh (33)

where εe and εh are the thresholds of the difference of power and heat generation for the CHPs attained from different local estimators, respectively.

The proposed asynchronous parallel computation scheme for the RTSE of CHPS is illustrated in Fig. 3.

Fig. 3 Asynchronous parallel computation scheme for RTSE of CHPS.

D. Overall Distributed Scheme for RTSE of CHPS

By combining the aforementioned techniques, a multi-time-scale asynchronous distributed scheme is developed for the RTSE of CHPS, which is presented in Algorithm 1.

Algorithm 1  : multi-time-scale asynchronous distributed RTSE of CHPS

1:

Initialize k=1. Obtain the initial states and perform network partition.

2:

for k=1 to kmax do

Determine whether it is time to perform RTSEs for the power and the heat systems simultaneously.

3:

if k=3K, KN+ then

4:

Perform RTSEs for the power and the heat systems simultaneously in an asynchronous parallel manner as described in Section IV-C.

5:

else

6:

Perform RTSEs for the power system individually.

7:

end if

8:

end for

The implementation of the proposed method in reality calls for the availability of the measurements and the information interchange between these two systems. As for the availability of the measurements, with the development of technology, more real-time measurement devices are installed in the CHPS [

24]. With additional pseudo-measurements, the system observability for the CHPS can be guaranteed. As for the information interchange, by adopting the proposed distributed scheme, only a small volume of information, i.e., the power and the heat output of the coupling units, needs to be exchanged between the two systems. Therefore, it is feasible to implement the state estimation in the CHPS.

V. Case Study

The proposed distributed CKF-RTSE method is tested on a 26-bus CHPS. A 65-bus CHPS is also used to demonstrate its scalability. The tests are run on a 4-core 3.3 GHz laptop based on MATLAB R2017b.

To evaluate the effectiveness of the proposed method, the following indices are used:

RMSE=1Mi=1M1nj=1nx̂j,est-xj,real2 (34)
SHSM=1Mi=1M1mj=1mh(x̂j,est)-h(xj,real)σj21Mi=1M1mj=1mzj,meas-h(xj,real)σj2 (35)
ηeffic=Tc-TdTc×100% (36)

where RMSE is the root-mean-square error for state estimation; M is the number of Monte Carlo simulations, which is set to be 200 for each test; n is the number of state variables; x̂j,est and xj,real are the estimated and true values of the state xj, respectively; SM and SH are the average value of the weighted residuals of the RTSE data and the measurement errors, respectively; SH/SM is the magnitude of the residuals of the estimates relative to the true values [

10], and a smaller value corresponds to lower residuals of the estimation results and more accurate state estimation results; σi is the standard deviation of measurement noise for the ith measurement; m is the number of measurements; ηeffic is the efficiency improvement when distributed RTSE is applied instead of a centralized one; and Tc and Td are the execution time for centralized and distributed RTSE methods, respectively.

The 26-bus CHPS consists of a 13-bus power system, two CHP units, and a 13-bus heat system, as shown in Fig. 4, whose parameters can be found in Appendix A.

Fig. 4 Schematic diagram of 26-bus CHPS.

Two CHPs are installed in the CHPS. CHP1, which is a gas turbine, is placed on bus 3 in the power system and bus 12 in heat system, respectively. CHP2, which is a steam turbine, is placed on bus 2 in the power system and bus 13 in heat system, respectively. It is assumed that the supply temperature for each heat source is 100 ℃ while the outlet temperature (return temperature before mixing) at each heat load is 50 ℃. The measurement deployment can be found in Table I. The proposed RTSE is implemented every 5 min for the power system and every 15 min for the heat system during a 24 hour simulation period. The time step in the following figures is chosen as 5 min. The load scaling factor of the power system for one day can be found in [

25] while that of the heat system is demonstrated in Fig. 5. The noises for real-time measurement and pseudo-measurement follow Gaussian distributions, whose triple standard deviations are set to be 1% and 20% of their true values, respectively. The thresholds of the individual estimations for both the power system and heat system are set to be 1×10-4; while that of global consistency is set to be 1×10-3.

TABLE I Measurement Deployment for 26-bus CHPS
SystemMeasurement typeVariableMeasurement location
Power system PMU Vbus, θbus Bus 2
SCADA Pbranch Lines 4-5 and 9-10
Ibranch Lines 6-7 and 8-12
Vbus Bus 3
Pseudo-measurement Pbus All buses
Heat system Nodal temperature Ts, Tr Buses 1, 7, 9, and 12
Pseudo-measurement Ф All buses

Fig. 5 Load scaling factor for heat system.

A. Comparison Among Different State Estimation Algorithms

To verify the validity of the proposed CKF-RTSE method, tests are carried out on the 26-bus CHPS using different estimation methods, as listed in Table II. The comparison of the state estimates among these methods for one Monte Carlo run is demonstrated in Figs. 6-9.

TABLE II Five Different State Estimation Methods for Comparison
MethodEstimation algorithmModeling of heat system
Static state estimation (SSE) WLS Without heat dynamics
Quasi-dynamic state estimation 1 (QDSE1) [8] WLS Heat dynamic model based on node method [7]
Quasi-dynamic state estimation 2 (QDSE2) WLS Proposed difference-based heat dynamic model
EKF-RTSE EKF Proposed difference-based heat dynamic model
CKF-RTSE CKF Proposed difference-based heat dynamic model

Fig. 6 Estimated and true values of voltage magnitude of bus 7.

Fig. 7 Estimated and true values of voltage angle of bus 7.

Fig. 8 Estimated and true values of supply temperature of bus 7.

Fig. 9 Estimated and true values of return temperature of bus 7.

The corresponding statistical results in terms of RMSE as well as execution times over 200 Monte Carlo runs are presented in Table III. The comparison of SH/SM is shown in Fig. 10.

TABLE III Comparison of Performance Based on Different Methods
MethodAverage RMSE (10-3 p.u.)Average execution time (s)
VθTsTr
SSE 1.671 1.632 8.631 4.619 1.211
QDSE1 1.653 1.644 5.485 4.032 1.509
QDSE2 1.659 1.631 4.817 3.164 1.325
EKF-RTSE 0.462 0.663 3.042 2.721 0.256
CKF-RTSE 0.031 0.099 3.011 1.503 0.454

Fig. 10 Average SH/SM values for different methods.

It can be found that compared with other methods, the proposed CKF-RTSE has the best performance in terms of computation accuracy. Leveraging CKF, the states of the heat system can be estimated recursively during the simulation process. Moreover, by employing sampling techniques, the nonlinearity can be handled without linearization. This allows utilizing both the previous and the current measurements to estimate the system operation states. With the developed heat dynamic model and the dynamic estimator, the proposed CKF-RTSE can accurately track the system dynamic behavior caused by load fluctuations, yielding better performance than the others. Note that the estimation results of QDSE2 are close to or even better than those of QDSE1. The latter adopts the widely-used node method to build the thermal dynamic model. Therefore, the validity of the proposed difference-based model for CHPS to describe the heat dynamic characteristics is verified.

Besides, the DSE methods, the EKF-RTSE and CKF-RTSE, have better computation efficiency than the SSE and QDSEs. The reason is that by incorporating state prediction, the EKF-RTSE and CKF-RTSE can leverage better initial states for the filtering stage, which facilitates the convergence of the state estimation. However, CKF-RTSE is less efficient than EKF-RTSE. This is because 2n cubature points need to be propagated, and the implementation of a CKF is more computationally expensive than that of an EKF.

B. Comparison of Numerical Stability

To show the superiority of proposed CKF-RTSE over UKF-based one in numerical stability, simulations are carried out by using UKF estimators [

14] with different scaling parameters. The tuning parameters α for the three UKF estimators are set to be 0.1, 0.5, and 1, respectively; while the other parameters β and κ are set to be 0. The comparison of the estimation accuracy in terms of RMSE for one Monte Carlo run among the three UKF estimators and the proposed CKF-RTSE is shown in Figs. 11-14.

Fig. 11 RMSE for voltage magnitude.

Fig. 12 RMSE for phase angle.

Fig. 13 RMSE for supply temperature.

Fig. 14 RMSE for return temperature.

It can be found that the UKF estimator with α=0.1 (UKF1) halts its operation at 14:25 (the 173rd time step) because of the unavailability of the square-root covariance, while the UKF with α=0.5 (UKF2) also has reduced accuracy. This indicates that the parameters of UKF have a great impact on its performance, especially its numerical stability. If appropriate parameters are utilized, the UKF can obtain state estimates with good accuracy, as shown by the performance of the UKF with α=1 (UKF3). Otherwise, the UKF may have poor accuracy or even encounter numerical instability since its estimation error covariance matrix Pk may fail to keep positive definite. Besides, the appropriate parameters are not generally applicable and have to be tailored in each case, which brings difficulties to the application of the UKF. On the contrary, the CKF can ensure the positive definitiveness of Pk, and thus has better numerical stability than the UKF. In a sense, the CKF can be regarded as a particular case of the UKF, whose scaling parameters satisfy α=1, β=0, and κ=0 [

26].

C. Sensitivity Analysis of Estimation Accuracy to Measurement Noises

Different levels of noises are added to the measurements to analyze the sensitivity of the proposed method to measurement noises. The test results are shown in Table IV, where the RMSE is used as the evaluation index. It can be observed that the proposed CKF-RTSE method has good robustness to measurement noises. The reason is that the use of developed differential equations allows to predict the system states one-step ahead and the prediction information increases the measurement redundancy. The latter is one of the most important factors to enhance the robustness of a method against measurement noises.

TABLE IV Comparison with Different Levels of Noises
Noise level (%)Average RMSE (10-3 p.u.)
Real-time measurement Pseudo-measurement V θ Ts Tr
1 20 0.031 0.099 3.011 1.503
50 0.033 0.112 3.223 1.584
5 20 0.037 0.123 7.562 5.405
50 0.040 0.126 9.051 7.712

D. Scalability of Proposed RTSE for CHPS

To verify the scalability of the proposed method, the tests on different CHPSs with diverse scales are carried out, namely the 26-bus CHPS and a 65-bus CHPS comprising a modified IEEE 33-bus distribution system [

27] and the 32-bus Barry Island district heating network [28]. Three CHPs are installed in the 65-bus CHPS. CHP1, which is a gas turbine, is placed on bus 7 in the power system and bus 1 in the heat system. CHP2, a steam turbine, is placed on bus 24 in the power system and bus 31 in the heat system. CHP3, which is also a gas turbine, is placed on bus 32 in the power system and bus 32 in the heat system. The measurement configuration can be found in Table V.

TABLE V Measurement Configuration for 65-bus CHPS
SystemMeasurement typeVariableMeasurement location
Power system PMU Vbus, θbus Buses 2 and 12
SCADA Pbranch Lines 6-7, 20-21, and 24-25
Ibranch Lines 16-17 and 29-30
Vbus Buses 17 and 30
Pseudo-measurement Pbus All buses
Heat system Nodal temperature Ts, Tr Buses 2, 7, 15, 19, and 28
Pseudo-measurement Ф All buses

The comparison results between the proposed distributed CKF-RTSE (D-RTSE) and its centralized counterpart (C-RTSE) are shown in Table VI. It is observed that the proposed method achieves much higher computation efficiency, i.e., 39.95% less computing time than the centralized one in the 26-bus system. This is expected as the distributed implementation allows to simplify the scale of the whole problem. That benefit can be further highlighted for larger-scale systems, i.e., over 50% less computation time compared with the centralized one in the 65-bus CHPS. A side effect is the slight reduction of estimation accuracy of the proposed method since only measurements within their own systems are leveraged. Nevertheless, the proposed method can gain high computation efficiency without loss of accuracy.

TABLE VI Performance of CKF-RTSE in Different Computation Schemes
Test systemMethodAverage RMSE (10-3 p.u.)Average execution time (s)ηeffic (%)
VθTsTr
26-bus D-RTSE 0.031 0.099 3.011 1.503 0.454 39.95
C-RTSE 0.028 0.076 2.713 1.281 0.756
65-bus D-RTSE 0.140 0.476 2.173 3.938 2.503 53.20
C-RTSE 0.123 0.420 1.875 3.407 5.348

VI. Conclusion

This paper proposes a distributed CKF-RTSE method for the real-time monitoring of CHPS. The proposed method can significantly improve computation efficiency without loss of estimation accuracy. Important conclusions are summarized as follows:

1) A difference-based model is derived for the heat system to describe its heat dynamic characteristics, making it feasible to predict the system states based on the last state snapshot.

2) A CKF-RTSE method is proposed for CHPS considering heat dynamics. Compared with the static and the quasi-dynamic state estimations, it can significantly improve the accuracy and efficiency of CHPS state estimation.

3) The computation burden is relieved by the proposed asynchronous multi-time distributed implementation scheme. Meanwhile, the privacy of different energy systems can be protected.

In our future works, the measurement placement problem of the CHPS will be investigated. The event-trigger strategy will also be considered to achieve the synchronous alignment of measurements.

Appendix

Appendix A

The parameters for the 26-bus CHPS can be found in [

29], which are listed in detail in Tables AI-AIV.

TABLE AI Line Segment Data of Power System
Line No.From busTo busImpedance (p.u.)
1 13 1 0.020000+j0.016000
2 1 2 0.008205+j0.019207
3 1 3 0.008205+j0.019207
4 3 4 0.008205+j0.019207
5 1 5 0.008205+j0.019207
6 5 6 0.008205+j0.019207
7 6 7 0.008205+j0.019207
8 6 8 0.008205+j0.019207
9 5 9 0.008205+j0.019207
10 9 10 0.008205+j0.019207
11 5 11 0.008205+j0.019207
12 2 12 0.008205+j0.019207
TABLE AII Electrical Load Values of Power System
Bus No.Load value
Active power (MW)Reactive power (Mvar)
1 0.200 0.116
2 0.500 0.125
3 0.800 0.400
4 0.800 0.290
5 1.155 0.660
6 0.800 0.400
7 0.170 0.080
8 0.128 0.086
9 0.170 0.151
10 0.582 0.462
11 0.100 0.050
12 0.230 0.132
13 0.000 0.000
TABLE AIII Pipe Parameter of Heat System
Pipe No.From busTo busLength (m)Diameter (mm)
1 13 1 500 200
2 1 2 400 200
3 2 3 600 200
4 4 3 500 200
5 12 4 600 200
6 1 5 200 200
7 1 6 150 200
8 2 7 180 200
9 2 8 150 200
10 3 9 100 200
11 3 10 110 200
12 4 11 90 200
TABLE AIV Heat Load Value of Power System
Bus No.Heat load (MW)
1 0.2
2 0.2
3 0.2
4 0.2
5 0.2
6 0.2
7 0.1
8 0.1
9 0.3
10 0.2
11 0.2

REFERENCES

1

J. Wu, J. Yan, H. Jia et al., “Integrated energy systems,” Applied Energy, vol. 167, pp. 155-157, Apr. 2016. [百度学术

2

A. Shabanpour-Haghighi and A. R. Seifi, “An integrated steady-state operation assessment of electrical, natural gas, and district heating networks,” IEEE Transactions on Power Systems, vol. 31, no. 5, pp. 3636-3647, Sept. 2016. [百度学术

3

J. Gustafsson, J. Delsing, and J. V. Deventer, “Improved district heating substation efficiency with a new control strategy,” Applied Energy, vol. 87, no. 6, pp. 1996-2004, Jun. 2010. [百度学术

4

J. Dong, H. Sun, Q. Guo et al., “State estimation for combined electricity and heat networks,” Power System Technology, vol. 40, no. 6, pp. 1635-1641, Jun. 2016. [百度学术

5

Y. Du, W. Zhang, and T. Zhang, “ADMM based distributed state estimation for integrated energy system,” CSEE Journal of Power and Energy Systems, vol. 5, no. 2, pp. 275-283, Jun. 2019. [百度学术

6

I. Gabrielaitiene, B. Bøhm, and B. Sunden, “Modelling temperature dynamics of a district heating system in Naestved, Denmark–a case study,” Energy Conversion and Management, vol. 48, no. 1, pp. 78-86, May 2006. [百度学术

7

Z. Li, W. Wu, J. Wang et al., “Transmission-constrained unit commitment considering combined electricity and district heating networks,” IEEE Transactions on Sustainable Energy, vol. 7, no. 2, pp. 480-492, Apr. 2016. [百度学术

8

T. Sheng, Q. Guo, H. Sun et al., “Two-stage state estimation approach for combined heat and electric networks considering the dynamic property of pipelines,” in Proceedings of 2017 International Conference on Applied Energy, Cardiff, UK, Dec. 2017, pp. 3014-3019. [百度学术

9

T. Zhang, Z. Li, Q. Wu et al., “Dynamic state estimation of combined heat and power system considering quasi-dynamics of temperature in pipelines,” in Proceedings of 2018 International Conference on Power System Technology, Guangzhou, China, Nov. 2018, pp. 232-237. [百度学术

10

T. Zhang, Z. Li, Q. Wu et al., “Decentralized state estimation of combined heat and power systems using the asynchronous alternating direction method of multipliers,” Applied Energy, vol. 248, pp. 600-613, May 2019. [百度学术

11

F. C. Schweppe, “Power system static-state estimation, Part III: implementation,” IEEE Transactions on Power Apparatus and Systems, vol. PAS-89, no. 1, pp. 130-135, Jan. 1970. [百度学术

12

J. Zhao, A. Gomez-Exposito, M. Netto et al., “Power system dynamic state estimation: motivations, definitions, methodologies and future work,” IEEE Transactions on Power Systems, vol. 34, no. 4, pp. 3188-3198, Jul. 2019. [百度学术

13

P. Du, Z. Huang, Y. Sunet et al., “Distributed dynamic state estimation with extended Kalman filter,” in Proceedings of 2011 North American Power Symposium, Boston, USA, Aug. 2011, pp. 1-6. [百度学术

14

G. Valverde and V. Terzija, “Unscented Kalman filter for power system dynamic state estimation,” IET Generation, Transmission & Distribution, vol. 5, no. 1, pp. 29-37, Jan. 2011. [百度学术

15

Y. Zhao, “Performance evaluation of cubature Kalman filter in a GPS/IMU tightly-coupled navigation system,” Signal Processing, vol. 119, pp. 67-79, Feb. 2016. [百度学术

16

S. Li, Y. Hu, L. Zheng et al., “Stochastic event-triggered cubature Kalman filter for power System dynamic state estimation,” IEEE Transactions on Circuits and Systems II: Express Briefs, vol. 66, no. 9, pp. 1552-1556, Sept. 2019. [百度学术

17

Z. Pan, Q. Guo, and H. Sun, “Feasible region method based integrated heat and electricity dispatch considering building thermal inertia,” Applied Energy, vol. 192, no. 15, pp. 395-407, Apr. 2017. [百度学术

18

X. Liu, N. Jenkin, J. Wu et al., “Combined analysis of electricity and heat networks,” Energy Procedia, vol. 64, pp. 155-159, Mar. 2015. [百度学术

19

J. Zheng, Z. Zhou, J. Zhao et al., “Function method for dynamic temperature simulation of district heating network,” Applied Thermal Engineering, vol. 123, no. 1, pp. 682-688, Aug. 2017. [百度学术

20

J. Fang, Q. Zeng, X. Ai et al., “Dynamic optimal energy flow in the integrated natural gas and electrical power systems,” IEEE Transactions on Sustainable Energy, vol. 9, no. 1, pp. 188-198, Jun. 2018. [百度学术

21

I. Arasaratnam and S. Haykin, “Cubature Kalman filters,” IEEE Transactions on Automatic Control, vol. 54, no. 6, pp. 1254-1269, Jun. 2009. [百度学术

22

T. Zhang, P. Yuan, Y. Du et al., “Robust distributed state estimation of active distribution networks considering communication failures,” International Journal of Electrical Power and Energy Systems, vol. 118, pp. 1-11, Jun. 2020. [百度学术

23

S. Li, Y. Li, Z. Li et al., “Event-trigger heterogeneous nonlinear filter for wide-area measurement systems in power grid,” IEEE Transactions on Smart Grid, vol. 10, no. 3, pp. 2752-2764, May 2019. [百度学术

24

T. Fang and R. Lahdelma, “State estimation of district heating network based on customer measurements,” Applied Thermal Engineering, vol. 73, no. 1, pp. 1211-1221, Dec. 2014. [百度学术

25

Southern California Edison. (2011, Dec.). Regulatory information-SCE load profiles. [Online]. Available: https://www.sce.com/wps/wcm/connect/sce_content_zh/content/regulatory/sce%20load%20profiles [百度学术

26

R. Turner and C. E. Rasmussen, “Model based learning of sigma points in unscented Kalman filtering,” Neurocomputing, vol. 80, pp. 47-53, Mar. 2012. [百度学术

27

IEEE. (2012, Dec.). 2012 IEEE test feeder specifications. [Online]. Available: http://ewh.ieee.org/soc/pes/dsacom/testfeeders/index.html [百度学术

28

X. Liu, “Combined analysis of electricity and heat networks,” Ph.D. dissertation, Institute of Energy, Cardiff University, Cardiff, UK, 2013. [百度学术

29

Y. Wang, B. Zeng, J. Guo et al., “Multi-energy flow calculation method for integrated energy system containing electricity, heat and gas,” Power System Technology, vol. 40, no. 10, pp. 2942-2950, Oct. 2016. [百度学术