Low-carbon optimization scheduling of multi-park integrated energy systems considering hybrid game theory mechanisms

1 Introduction

Under the “dual carbon” goals, expediting the development of a new power system is critical. The Park Integrated Energy System (PIES) imports energy from power and natural gas pipeline networks. Relying on flexible conversion equipment inside the park, the system provides energy to park users. It realizes multi-energy complementarity of flexible cross-regional resources, supports low-carbon cost-effective operation, and fully meets daily energy consumption requirements of users (Mancarella, 2014; Zhihong et al., 2018; Karki et al., 2007; Zhang et al., 2024).

Most existing studies focus on the optimal operation strategies inside PIES, yet overlook the coupling relationships between various park subsystems as well as their connections with the upstream power grid and natural gas network. By investigating distributed collaborative optimization strategies between upstream energy markets and PIES, game theory can maximize the benefits of all participants while accounting for the overall system revenue. Shui et al. (2024) puts forward an optimization strategy based on real-time dynamic games between upstream energy markets and park operators. Ma et al. (2024) constructs a Stackelberg game framework between energy suppliers and multiple PIES. Stackelberg games have inherent limitations: lower-level followers can only make responsive decisions according to the plans formulated by upper-level leaders. Zheng et al. (2025) converts the multi-agent power trading problem into a cooperative game model, realizing rational benefit allocation between collective and individual participants. Jin et al. (2025); Xingyu et al. (2023); Wang et al. (2026) build multi-stakeholder cooperative game models and distribute revenues among all participants via the Shapley value method. Nevertheless, the Shapley value only partially reflects participants’ contributions and cannot fully quantify each agent’s input throughout the entire game process. Yang et al. (2026) develops a cooperative scheduling model with power sharing among multiple micro grids, and achieves equal benefit distribution across all stakeholders via Nash bargaining. The aforementioned literature only discusses energy interactions between upstream energy markets and individual parks, ignoring the exchange of flexible resources among different parks. Accordingly, simply adopting game theory alone fails to maximize collective benefits for the energy market and all participating parks. Yang et al. (2026) also proposes a game-theory-based peer-to-peer distributed trading strategy. Dynamic electricity and heat prices are delivered to park operators at the lower tier via the dynamic pricing mechanism derived from upper-level games. On this basis, a bi-level optimization model is constructed, which integrates the Stackelberg dynamic pricing mechanism and cooperative games for each park integrated energy system. Gou et al. (2024) establishes a bi-level game-based optimal operation strategy. At the upper layer, a cooperative game is formulated to minimize the total operating cost of multi-region integrated energy systems. At the lower layer, each park implements a non-cooperative game to minimize its own energy consumption cost. Simulation results verify that the proposed strategy can realize win-win outcomes for all stakeholders in the system.

In addition, most current research on the game interaction of multi-park integrated energy systems focuses on economic aspects and lacks consideration for the low-carbon economic operation of the system. Therefore, how to achieve the low-carbon economic operation of the system is the focus of future research. The introduction of carbon emission rights and green certificate trading mechanisms can realize the absorption of renewable energy by establishing a carbon trading market, thereby achieving the goal of reducing carbon emissions. Zhou and Wei, 2025; Xu et al. (2024); Li et al. (2026a); Liu et al. (2025) introduce the tiered carbon trading mechanism into the integrated energy system and consider the impact of tiered carbon trading parameters on the system. Sun and Qin (2023) applies the green trading mechanism to the block chain, realizing the efficiency and transparency of the green certificate trading mechanism. Zheng et al. (2024); Li et al. (2026b); Yan et al. (2026) introduce both carbon trading mechanism (CET) and carbon trading mechanism (GCT) models into the low-carbon optimization model of the integrated energy system, promoting the low-carbon economic transformation of the park.

In summary, this paper proposes a master-follower game model with integrated energy operators in industrial parks as leaders and users as followers, while simultaneously implementing a two-layer hybrid game optimization model for electricity and heat sharing among PIES. First, the operational framework of a multi-park integrated energy system is introduced, and the electricity and heat transactions between integrated energy operators and various park integrated energy systems, as well as the behavior of stakeholders in the PIES, are described in detail, addressing the issue of conflicting interests. Second, in the context of the current GCT-CET market, a GCT-CET collaborative trading mechanism is proposed, maximizing system benefits through the conversion of green certificates into carbon allowances, while exploring the positive impact of the carbon trading mechanism on the integrated energy system. Finally, based on Nash negotiation theory, the problem of profit distribution among PIES is transformed into two easily solvable sub-problems: minimizing cooperation costs and maximizing cooperation profit distribution. Simulation examples verify that the model can effectively reduce overall carbon emissions and balance the interests of all stakeholders.

2 Multi-park integrated energy system operation architecture

The operational framework of the multi-park integrated energy system proposed in this paper is shown in Figure 1. The system mainly includes the upper-level power grid, the natural gas network, and multiple park integrated energy systems. The integrated energy system of a single park includes a combined heat and power unit (CHP), gas boiler (GB), electric energy storage, thermal energy storage, gas energy storage, photovoltaic, wind power and a refined hydrogen production system.

During system operation, the energy within a PIES is not perfectly balanced. Multiple PIES can also communicate with each other using their own resources and with external parks, resulting in differences in supply and demand at different times. PIES interact through electrical and thermal power lines, forming cooperative alliances through energy sharing. PIES can choose trading partners based on their internal needs in energy transactions. By introducing a DSO (Integrated Energy System Operator) as an intermediary energy supplier, PIES transmit the energy information they need to trade to the DSO for transaction. The DSO conducts transactions with the upstream power grid based on the electricity purchase and sales volumes of the PIES alliance; it settles accounts using the difference between grid electricity prices and feed-in tariffs to generate profit. Regarding electricity trading, when internal load demand within the park exceeds local supply, priority is given to exchanging power with other parks; if requirements still cannot be met, electricity is purchased from the upstream grid. The hybrid game model of the DSO and PIES alliance constructed in this paper is shown in Figure 2.

The model proposed in this paper can be divided into two stages: Stage 1 is a master-slave game model between the DSO and the PIES cooperative alliance; the DSO, as the upper-level leader, can determine the optimal electricity price by considering the energy price constraints of the upper-level distribution network and the electricity trading needs among the lower-level PIES alliances, thereby maximizing the system’s revenue. Stage 2 is a cooperative game model among the PIES; based on the electricity trading price information provided by the DSO, the electricity and heat energy trading volume and price among the PIES can be calculated, and the electricity purchase and sale plan can be transmitted to the DSO while satisfying the needs of the PIES alliance and maximizing the cooperative benefits.

3 Carbon-green certificate joint trading and integrated demand response model

3.1 Green certificate trading model (GCT)

Green certificates, also known as green utility certificates, are positively correlated with the amount of renewable energy generated. They represent national certification of renewable energy power generation and serve as proof of demand-side consumption of green electricity. Similar to carbon trading mechanisms, green certificate trading mechanisms leverage the market’s role in optimizing resource allocation through trading. When the number of green certificates in an energy system exceeds its free green certificate quota, the surplus can be sold on the green certificate trading market to generate revenue. Conversely, insufficient green certificates must be purchased from the market to meet the quota. Formula for calculating the transaction cost of green certificates as shown in Formulas 1–4:

N q u = δ q u ∑ t = 1 T P load t (1)

N r e = δ re ∑ t = 1 T P r e t (2)

C i g r e e n = η green N re − N q u (3)

I g r e e n = γ N r e − N q u (4)

In the formula: N q u represents the quota quantity of green certificates allocated by the government based on system load; δ q u represents the quota coefficient of green certificate quantity; P load t represents system load; N re represents the quantity of green certificates obtained by the system through its own consumption of new energy; δ re represents the conversion coefficient for converting renewable energy into green certificates, with 1 green certificate corresponding to 1 MWh of renewable energy settlement; P r e t represents the amount of renewable energy consumed by the system; C i g r e e n represents the transaction cost of system green certificates; η green represents the transaction unit price of green certificates, set at 8 CNY/per certificate; I g r e e n represents the carbon emission reduction achieved through the green certificate trading mechanism; γ represents the quantity of green certificates that can offset 1 ton of carbon dioxide emissions; T represents the system operation cycle.

3.2 Carbon emission trading model (CET)

This paper adopts a carbon emission market based on the benchmark method for the free allocation of initial carbon emission allowances. The heat load in the PIES system is mainly provided by gas-fired boilers and cogeneration units, while the electrical load comes only from the upstream power grid and cogeneration units. In this paper, upstream power generation is uniformly considered as thermal power generation. Since the carbon allowance coefficients between CHP and GB are relatively similar, a unified carbon allowance coefficient is used for calculation. Therefore, the initial carbon emission quota model is shown in Formula 5:

I I E S q = α h ∑ t = 1 T H i , t G B , H + P i , t C H P + λ c h H i , t C H P + α e ∑ t = 1 T P i , t g r i d (5)

In the formula: I I E S q is the free carbon emission allowance of PIES; α e and α h are the carbon emission allowance coefficients per unit of purchased electricity and per unit of purchased heat, respectively; λ c h is the converted power of electricity and heat; P i , t C H P is the electric power output by the CHP in PIES i at time t; H i , t G B , H is the thermal power output by the GB in PIES i at time t; H i , t C H P is the thermal power output by the CHP in PIES i at time t; P i , t g r i d represents the amount of electricity purchased from the upper-level power grid at time t (PIES i).

The actual carbon emission of the system is calculated as shown in Formula 6:

I I E S e m = β g ∑ t = 1 T H i , t G B + P i , t C H P + λ c h H i , t C H P + β e ∑ t = 1 T P i , t g r i d (6)

In the formula: I I E S e m is the actual carbon emissions of the system; β e is the unit carbon emission coefficient of the power purchased from the upstream power grid; β g is the actual carbon emission coefficient of the system per unit heat power.

Therefore, the actual carbon emission trading share of the system is as shown in Formulas 7, 8:

I I E S = I I E S e m − I I E S q − φ I g r e e n (7)

κ = x C E T φ (8)

In the formula: I I E S represents the actual share of carbon emissions traded by PIES; x C E T represents the base price for carbon trading; φ represents the unit price of green certificate trading; κ is the conversion coefficient.

3.3 Joint mechanism for green certificates and carbon trading

In traditional carbon trading mechanisms, the base price is fixed; in contrast, tiered carbon trading mechanisms divide the volume of traded carbon emissions into several purchasing brackets, meaning that carbon trading costs rise as the system’s trading volume increases.

Since the carbon emission reductions from new energy power generation can be quantified, a joint green certificate and carbon trading mechanism can be implemented through tiered carbon trading and green certificate trading systems.

The tiered carbon trading cost is shown in Formula 9:

C i c o 2 = λ I I E S , I I E S ≤ l λ 1 + α I I E S − l + λ l , l < I I E S ≤ 2 l λ 1 + 2 α I I E S − 2 l + 2 + α λ l , 2 l < I I E S ≤ 3 l λ 1 + 3 α I I E S − 3 l + 3 + 3 α λ l , 3 l < I I E S ≤ 4 l λ 1 + 4 α I I E S − 4 l + 4 + 6 α λ l , 4 l < I I E S (9)

In the formula: C i c o 2 represents the carbon trading cost at time t; λ is the carbon trading benchmark value; l is the length of the carbon emission interval; and α is the price growth rate of carbon trading.

3.4 System internal equipment model

• The CHP model is shown in Formulas 10–13

P i , t C H P = η P C H P P i , t G T (10)

P i , t C H P , H = η H C H P P i , t G T (11)

P i , ⁡ min G T ⩽ P i , t G T ⩽ P i , ⁡ max G T (12)

− Δ P i , ⁡ max G T ⩽ P i , t G T − P i , t − 1 G T ⩽ Δ P i , ⁡ max G T (13)

In the formula: P i , t C H P and H i , t C H P are the electrical power and thermal power output of CHP in PIES i at time t; η P C H P and η H C H P are the power generation efficiency and thermal generation efficiency of the gas turbine; P i , t − 1 G T and P i , t G T are the amount of natural gas input at times t and t-1; P i , ⁡ max G T and P i , ⁡ min G T are the upper and lower limits of the input natural gas, respectively; Δ P i , ⁡ max G T is the maximum ramp power of the gas turbine input natural gas.

• 2. The GB model is shown in Formulas 14–16

H i , t G B = η G B G i , t G B (14)

H i , ⁡ min G B ⩽ H i , t G B ⩽ H i , ⁡ max G B (15)

− Δ H i , ⁡ max G B ⩽ H i , t G B − H i , t − 1 G B ⩽ Δ H i , ⁡ max G B (16)

In the formula: H i , ⁡ max G B and H i , ⁡ min G B are the upper and lower limits of the output power of the gas boiler; η G B is the heat conversion efficiency of the gas boiler; G i , t G B is the power of natural gas consumed by the gas boiler; Δ H i , ⁡ max G B is the maximum ramp-up power of the gas boiler; H i , t G B and H i , t − 1 G B are the GB heat output at time t and t-1.

• 3. Energy storage devices

The energy storage model is shown in Formulas 17–21:

0 ≤ P i , t B T ⩽ P i , ⁡ max B T (17)

0 ⩽ P i , t c h r ⩽ u i , t c h r P i , c h c h r , ⁡ max (18)

0 ⩽ P i , t d i s ⩽ u i , t d i s P i , c h d i s , ⁡ max (19)

u i , t c h r + u i , t d i s ≤ 1 (20)

P i , t B T = P i , t − 1 B T 1 − λ α + η c h r E P i , t c h r − P i , t d i s η d i s E (21)

In the formula: P i , t B T and P i , t − 1 B T represent the electrical energy storage capacity at times t and t-1, respectively; P i , ⁡ max B T represents the upper limit of the energy storage capacity; P i , c h c h r , ⁡ max and P i , c h d i s , ⁡ max represent the electrical energy storage charging and discharging power; u i , t c h r and u i , t d i s are 0–1 variables used to limit simultaneous charging and discharging; P i , c h c h r , ⁡ max and P i , c h d i s , ⁡ max represent the maximum charging and discharging efficiency of the electrical energy storage; η c h r E and η d i s E represent the charging and discharging efficiency of the electrical energy storage; λ α is the energy consumption coefficient.

The energy storage devices described in this paper include electrical energy storage, thermal energy storage, and hydrogen energy storage. Since they operate on the same principle, the same model can be used.

• 4. The EL model is shown in Formulas 22, 23:

H y i , t EL = η EL P i , t E L (22)

P min E L ≤ P i , t E L ≤ P max E L (23)

In the formula: H y i , t EL is the hydrogen production power of EL; η EL is the hydrogen production efficiency of EL; P i , t E L is the electrical power of EL; P min E L and P max E L are the minimum and maximum electrical power of EL.

• 5. The HFC model is shown in Formulas 24–26:

P i , t H F C = η HFC e H y i , t H F C (24)

H i , t H F C = η HFC h H y i , t H F C (25)

H y min H F C ≤ H y i , t H F C ≤ H y max H F C (26)

In the formula: H y i , t H F C is the hydrogen power input to the HFC; η HFC e and η HFC h are the electrical and thermal power conversion efficiencies of the HFC, respectively; P i , t H F C and H i , t H F C are the electricity and heat production power of the HFC, respectively; H y max H F C and H y min H F C are the maximum and minimum values ​​of the hydrogen power input to the HFC, respectively.

• 6. The MR model is shown in formulas Formulas 27, 28:

Q i , t M R = η M R H y i , t M R (27)

H y min M R ≤ H y i , t M R ≤ H y max M R (28)

In the formula: η M R is the methanation efficiency of MR; η M R is the electro-conversion efficiency of MR; Q i , t M R and H y i , t M R are the hydrogen power entering and exiting MR and the gas power output, respectively; H y max M R and H y min M R are the upper and lower limits of the hydrogen power input to MR, respectively.

4 Two-level game model

4.1 Master-slave game between DSO and PIES alliance

• Leader Effectiveness Model

As the upper-level leader in the game theory model, the DSO aims to maximize its own benefits, which includes the costs and benefits of the upper-level distribution network and the PIES alliance in purchasing and selling electricity.

max ⁡ C D S O = ∑ t = 1 24 u t I E , b − u t P E , b ∑ i = 1 N P i , t P E , b − ∑ t = 1 24 u t I E , s − u t P E , s ∑ i = 1 N P i , t P E , s (29)

In the formula: C D S O represents the transaction revenue of the DSO; N represents the number of PIES participating in the game; u t P E , s represents the electricity price sold by the DSO to the PIES alliance at time t; u t P E , b represents the electricity purchase price of the DSO to the PIES alliance at time t; u t I E , b represents the electricity purchase price of the upstream distribution network to the DSO at time t; u t I E , s represents the electricity price sold by the upstream distribution network to the DSO at time t; P i , t P E , s represents the electricity sold by the DSO to the i PIES at time t; and P i , t P E , b represents the electricity purchased by the DSO from the i PIES at time t.

To ensure the benefits for all parties in the system and to prevent the PIES alliance from setting excessively high energy prices due to the DSO maximizing its benefits, and to minimize the energy prices for the PIES alliance, a price constraint is imposed, namely,:

u min P E , b ≤ u t P E , b ≤ u max P E , b , ∑ t = 1 24 u t P E , b / 24 ≤ u ave P E , b u min P E , s ≤ u t P E , s ≤ u max P E , s , ∑ t = 1 24 u t P E , s / 24 ≤ u ave P E , s (30)

In the formula: u max P E , b and u min P E , b are the upper and lower limits of the DSO’s electricity purchase price, respectively; u max P E , s and u min P E , s are the upper and lower limits of the DSO electricity sales price, respectively; u ave P E , b and u ave P E , s are the average values ​​of the electricity purchase and sales prices set by the DSO within a period.

• 2. Follower Benefit Model

As a lower-level follower, the PIES Alliance aims to maximize its own profits, including energy purchase and sale revenue between the PIES Alliance and DSO, PIES inter-electricity and heat trading models, natural gas costs, equipment operation and maintenance costs, and green certificate carbon trading costs.

max ⁡ C i PIES = C i p i + C i trade − C i g r e e n − C i n g − C i o m − C i c o 2 (31)

C i p i = u t P E , b ∑ t = 1 24 P i , t P E , b − u t P E , s ∑ t = 1 24 P i , t P E , s (32)

C i trade = ∑ t = 1 24 ∑ j = 1 , j = i N u i − j , t E , tra P i − j , t E , t ra + ∑ t = 1 24 ∑ j = 1 , j + i N u i − j , t H , tra H i − j , t H , t r a (33)

C i n g = ∑ t = 1 24 u t n g P i , t C H P + λ c h H i , t C H P η chpe L g + H i , t G B η g b L g (34)

C i o m = ∑ t = 1 24 c i n P i , t n (35)

In the formula: C i PIES is the operating revenue of PIES; C i p i is the revenue from the purchase and sale of electricity and heat between PIES and DSO; P i − j , t E , t r a is the electricity trading volume between PIES i and PIES j at time t. H i − j , t H , t r a is the thermal energy trading volume between PIES i and PIES j at time t. C i trade is the revenue from electricity and heat transactions between PIES; C i n g is the cost of natural gas; η chpe represents CHP electrical energy conversion efficiency, L g represents natural gas calorific value, and η g b represents GB human thermal energy conversion efficiency; C i o m is the equipment operation and maintenance cost; u i − j , t E , t r a and u i − j , t H , t r a are the electricity and heat transaction volume of PIES i at time t; u t n g is the price of natural gas; is the lower heating value of natural gas; c i n is the maintenance cost coefficient of equipment i; P i , t n is the output power of equipment i at time t.

Formulas 29–35 define the electro thermal power balance constraints and energy trading constraints of PIES, where the specific electro thermal power balance constraints are, as shown in Formulas 36–41:

P i , t C H P + P i , t P V + P i , t W T + P i , t d i s + P i , t P E , s + P i , t H F C = P i , t P E , b + P i , t E , l o a d + P i , t c h r + P i , t E L + ∑ j = 1 , j ≠ i N P i − j , t E , t r a (36)

H i , t G B + H i , t C H P + H i , t d i s + H i , t H F C = P i , t H , l o a d + H i , t c h r + ∑ j = 1 , j ≠ i N H i − j , t H , t r a (37)

Q i , t M R + Q i , t d i s + Q i , t b u y = P i , t G , l o a d + Q i , t c h r + Q i . t G B + Q i , t C H P (38)

H y i , t M R + H y i , t H E S + H y i , t H F C = H y i , t EL (39)

P i − j , t E , t r a ≤ P i − j , ⁡ max E , t r a P i − j , t E , t r a + P j − i , t E , t r a = 0 (40)

H i − j , t H , t r a ≤ H i − j , ⁡ max H , t r a P i − j , t H , t r a + H j − i , t H , t r a = 0 (41)

In the formula: P i , t P V is the PV electric power of PIES i at time t; P i , t W T is the WT electric power of PIES i at time t; P i , t d i s is the electric power output by the BT in PIES i at time t; P i , t E , l o a d is the electric load of PIES i at time t; P i , t c h r is the electric power absorbed by the battery in PIES i at time t; P i , t E L is the electric power consumed by the EL in PIES i at time t; H i , t d i s is the thermal power output by the HST in PIES i at time t; P i , t H , l o a d is the thermal load of PIES i at time t; H i , t c h r is the thermal power absorbed by the HST in PIES i at time t; Q i , t M R is the natural gas power output by the MR in PIES i at time t; G i , t d i s is the natural gas power output by the GST in PIES i at time t; Q i , t b u y is the natural gas power purchased from the superior grid by PIES i at time t; P i , t G , l o a d is the natural gas load of PIES i at time t; Q i , t c h r is the natural gas power absorbed by the GST in PIES i at time t; Q i . t G B is the natural gas power consumed by the GB in PIES i at time t; Q i , t C H P is the natural gas power consumed by the CHP in PIES i at time t. H y i , t M R represents the hydrogen consumed by MR in PIES i during time t; H y i , t H E S represents the hydrogen consumed by MR in PIES i during time t; H y i , t H F C represents the hydrogen consumed by MR in PIES i during time t; H y i , t EL represents the hydrogen gas generated by EL in PIES i at time t; P i − j , ⁡ max E , t r a and H i − j , ⁡ max H , t r a are the maximum limits for electric and thermal energy trading power between PIES, respectively.

4.2 PIES coalition nash bargaining model

The coalition Nash bargaining model constructed in this paper is shown in Formula 42:

max ∏ i = 1 N C i − C i 0 s . t . C i − C i 0 ≥ 0 (42)

In the formula: N represents the total number of PIES participating in the game; C i 0 represents the negotiation breakdown point of PIES i, the benefit obtained before participating in the cooperative game; C i represents the gain obtained by PIES i after participating in the cooperative game. Since the constructed model is a non-convex nonlinear optimization problem, it is difficult to solve. Therefore, the model is transformed into two easier-to-solve sub problems, which are solved sequentially.

• The problem of maximizing the benefits of the PIES cooperation alliance is shown in Formula 43:

max ⁡ C PIES = ∑ i = 1 N C i p i − C i g r e e n − C i n g − C i a m − C i d r − C i c o 2 s . t . formula 16 − 20 (43)

• 2. Sub-problem of revenue distribution considering overall contribution rate

In a PIES alliance, each PIES maximizes its overall benefit by cooperating and choosing to share energy with other PIES. The standard Nash equilibrium solution ensures that all participants in the alliance have equal bargaining power, but it does not consider the differences in their contributions during cooperation. Therefore, by quantifying the overall energy contribution of each PIES during cooperation, and differentiating their bargaining power in the distribution of cooperative benefits, a fairer distribution of cooperative benefits can be achieved.

First, calculate the electrical energy E i sup and heat energy H i sup provided by each PIES during its participation in energy sharing, as well as the electrical energy E i r e c and heat energy H i rec obtained, as shown in Formulas 44–47:

E i sup = ∑ t = 1 T ⁡ max 0 , P i − j , t E , t r a (44)

E i r e c = − ∑ t = 1 T ⁡ min 0 , P i − j , t E , t r a (45)

H i sup = ∑ i = 1 T ⁡ max 0 , H i − j , t H , tra (46)

H i r e c = − ∑ t = 1 T ⁡ min 0 , H i − j , t H , t r a (47)

Secondly, based on the nonlinear energy sharing mapping method, the electrical energy contribution factor d i , e and thermal energy contribution factor d i , h of each PIES are quantified according to their contribution to the sharing of electrical and thermal energy, as shown in Formulas 48, 49:

d i , e = e E i m e p / E max m p − e − E i r e / E max r e (48)

d i , h = e H i a p / H max a p − e − H i r c e / H max m e (49)

Therefore, the comprehensive contribution factor d i for each park can be obtained as shown in Formula 50:

d i = ζ e d i , e + ζ h d i , h (50)

In the formula: ζ e and ζ h represent the weights of electrical energy contribution and thermal energy contribution, respectively, satisfying condition ζ e + ζ h = 1. Therefore, the cooperative benefit distribution model is obtained as shown in Formulas 51–54:

max ∑ i = 1 N d i ⁡ ln I i PIES + ∑ t = 1 24 ∑ j = 1 , j ≠ i N u i − j , t E , tra P i − j , t E , t r a * + ∑ t = 1 24 ∑ j = 1 , j ≠ i N u i − j , t H , t r a H i − j , t H , t r a * − I i 0 (51)

s . t . C i * + ∑ t = 1 24 ∑ j = 1 , j + i N u i − j , t E , t r a P i − j , t E , t r a * + u i − j , t H , t r a H i − j , t H , t r a * ≥ C i 0 (52)

u i , t P E , s * ≤ u i − j , t E , t r a ≤ u i , t P E , b * (53)

C i PIES = C i p i * − C i g r e e n * − C i n g * − C i a m * − C i d r * − C i c o 2 * (54)

In the formula: the superscript “ *” represents the optimal solution to sub problem one.

4.3 Model solution

As the game-theoretic framework proposed in this paper constitutes a bilevel optimization model, a two-stage solution approach is adopted. The upper-level model employs the bisection method to solve the leader-follower game, while the lower-level cooperative model utilizes the ADMM distributed algorithm to preserve the decision-making independence and information privacy of the individual subsystems (Qu et al., 2025; Li et al., 2026c).

• Phase 1 Model Solution

The bisection method is employed for the solution; let X k , t represent the energy price at time t during the k iteration, and assume X k , t is the upper bound for this point.

(1) When X k , t ≠ X k − 2 , t, add the constraint, as shown in Formula 55:

X k − 1 , t ≤ X k + 1 , t ≤ X k , t (55)

(2) When it is X k , t = X k − 2 , t, add the constraint, as shown in Formula 56:

x k + 1 , t ≡ x k − 1 , t + x k , t / 2 (56)

By means of constraint Equations 36, 37, the computational boundaries can be rapidly reduced during the respective iterative processes, thereby enhancing the computational efficiency of the algorithm. As shown in Formulas 57, 58:

u k + 1 , t P E , b − u k , t P E , b / u k , t P E , b ≤ δ (57)

u k + 1 , t P E , s − u k , t P E , s / u k , t P E , s ≤ δ (58)

During the respective iterative processes, convergence is determined by checking whether the convergence conditions defined by Equations 36, 37 are met. The energy trading prices determined by the IEO and the energy trading volumes between PIES units are used as initial parameters for Stage 2 to calculate the electricity and heat trading prices between the PIES units.

• 2. Phase 2 Model Solution

To protect the privacy of each stakeholder and to improve computational speed and convergence, this paper employs the ADMM to solve sub problem 2 in a distributed manner. As shown in Formulas 59, 60:

u i − j , t E , t r a = u j − i , t E , t r a = z i − j , t E , t r a (59)

u i − j , t H , t r a = u j − i , t H , t r a = Z i − j , t H , t r a (60)

In the formula: Z i − j , t E , tra and Z i − j , t H , tra define a shared variable between PIES i and PIES j.

After implementing the decoupling transformation, the solution is derived by considering the bargaining power based on comprehensive contribution during the transaction price negotiation process; the specific steps are as follows:

Step 1: First, calculate the comprehensive bargaining factor for each PIES based on the results obtained in Stage d i.

Step 2: Formulate the Lagrange function for Sub-problem 2, as shown in Formula 61:

min ⁡ L i = − d i ⁡ ln I i 0 − ∑ t = 1 24 ∑ j = 1 , j ≠ i N u i − j , t E , t r a P i − j , t E , t r a * + ∑ t = 1 24 ∑ j = 1 , j ≠ i N u i − j , t H , t r a P i − j , t H , t r a * − I i * + ∑ t = 1 24 ∑ j = 1 , j ≠ i N λ i − j , t E , tra u i − j , t E , tra − z i − j , t E , tra + ρ 2 u i − j , t E , tra − z i − j , t E , tra 2 2 + ∑ t = 1 24 ∑ j = 1 , j ≠ i N λ i − j , t H , tra u i − j , t H , tra − z i − j , t H , tra + ρ 2 u i − j , t H , tra − z i − j , t H , tra 2 2 (61)

Where: λ i − j , t E , tra and λ i − j , t H , tra are the Lagrange multipliers for electricity and heat trading prices, respectively; and ρ is the penalty factor, set to 1.

Step 3: PIES i updates decision variables, as shown in Formulas 62, 63:

u i − j , t E , t r a k + 1 = argmin L i u i − j , t E , t r a , z i − j , t E , t r a , λ i − j , t E , t r a (62)

u i − j , t H , t r a k + 1 = argmin L i u i − j , t H , t r a , z i − j , t H , t r a , λ i − j , t H , t r a (63)

Step 4: PIES j updates its own decision variables after receiving the decision variables from PIES i, as shown in Formulas 64, 65:

u j − i , t E , t r a k + 1 = argmin L i u j − i , t E , t r a , Z j − i , t E , t r a , λ j − i , t E , t r a (64)

u j − i , t H , t r a k + 1 = argmin L i u j − i , t H , t r a , z j − i , t H , t r a , λ j − i , t H , t r a (65)

Step 5: Repeat Steps 3 and 4 until each PIES has updated its trading strategy for this iteration and updated the shared variables, as shown in Formulas 66, 67:

z i − j , t E , tra k + 1 = u i − j , t E , tra k + 1 + u j − i , t E , tra k + 1 / 2 (66)

z i − j , t H , tra k + 1 = u i − j , t H , tra k + 1 + u j − i , t H , tra k + 1 / 2 (67)

Step 6: Update the Lagrange multipliers, as shown in Formulas 68, 69:

λ i − j , t E , tra k + 1 = λ i − j , t E , tra k + ρ u i − j , t E , tra + z i − j , t E , tra / 2 (68)

λ i − j , t E , tra k + 1 = λ i − j , t E , tra k + ρ u i − j , t E , tra + z i − j , t E , tra / 2 (69)

Step 7: Assess the convergence of the algorithm, as shown in Formulas 70, 71:

∑ t = 1 T u i − j , t E , t r a − Z i − j , t E , t r a 2 2 ≤ ω (70)

∑ t = 1 T u i − j , t H , t r a − Z i − j , t H , t r a 2 2 ≤ ω (71)

Where: ω is the convergence coefficient, set to 0.001.

5 Simulation examples

5.1 Basic data

To verify the effectiveness of the proposed model, this paper selects three typical micro grids (PIES1, PIES2, and PIES3) in a region of northern China to validate the effectiveness of the proposed method. As shown in Figure 3, the load of each micro grid is also shown.

The renewable energy load of each integrated energy system is shown in Figure 4.

The parameters of each piece of equipment are shown in Table 1 and the carbon trading benchmark price is 90 CNY/ton.

TABLE 1. Basic parameters.

The time-of-use pricing adopted by each PIES is shown in Figure 5.

5.2 PIES low-carbon operation analysis

To verify that the introduction of the green certificate-carbon trading mechanism and cooperation on electricity and heat energy among the PIES alliance results in better low-carbon performance, five scenarios were compared and analyzed.

Scenario 1: Do not consider the master-slave game between DSO and PIES alliance, and do not consider the electro thermal energy cooperation among alliance PIES.

Scenario 2: Considering the master-slave game between DSO and the PIES alliance, but not considering cooperation on electricity and heat energy among the PIES alliance.

Scenario 3: Considering the master-slave game between DSO and the PIES alliance, and considering cooperation on electricity and heat energy among the PIES alliance.

Scenario 4: Considering the master-slave game between DSO and the PIES alliance, considering cooperation on electricity and heat energy among the PIES alliance, and considering a tiered carbon trading mechanism.

Scenario 5: Considering the master-slave game between DSO and the PIES alliance, considering cooperation on electricity and heat energy among the PIES alliance, and considering the green certificate-carbon trading mechanism.

5.3 Results analysis

Comparing Scenario 5 with Scenarios 1-4, the costs and carbon emissions for each PIES considering various factors—including master-slave game dynamics, cooperation on electricity and heat energy between alliance PIES, master-slave game dynamics between DSO and PIES alliances, cooperation on electricity and heat energy between alliance PIES, a tiered carbon trading mechanism, and a green certificate-carbon trading mechanism—are shown in Table 2.

TABLE 2. Results under different operation schemes.

In Table 1, the total revenues for Scenarios 1-5 are RMB 23812.16, RMB 26800.1, RMB 27621.07, RMB 29340.45, and RMB 32899.29, respectively. Compared to Option 1, Option 2 yields a total revenue increase of 2,987.94 yuan and a total carbon emission reduction of 131.78 kg; Option 3 yields a total revenue increase of 820.97 yuan and a total carbon emission reduction of 903.61 kg; Option 4 yields a total revenue increase of 1719.38 yuan and a total carbon emission reduction of 269.29 kg; and Option 5 yields a total revenue increase of 3,558.84 yuan and a total carbon emission reduction of 412.53 kg. The results indicate that by introducing cooperation in electro thermal energy among PIES alliance members, the total revenue of the PIES alliance after cooperative operation is significantly improved compared to before independent operation. Furthermore, the green certificate-carbon trading mechanism can bring greater revenue to the PIES alliance while reducing the total carbon emissions of the system.

In conclusion, Option 5 offers the highest total revenue for the PIES alliance and the lowest user costs and carbon emissions.

• Analysis of the Master-Slave Game Result between DSO and PIES Alliance

In the hybrid game-based bi-level optimization scheduling framework proposed in this paper, the first stage involves upper-level cooperative game optimization to determine the electricity purchase and sales volumes—along with their corresponding prices at each time step—between the Distribution System Operator (DSO) and the various Prosumer-based Integrated Energy Systems (PIES) within the cooperative alliance. The second stage involves lower-level game optimization to determine the trading prices and volumes of electricity and thermal energy among the integrated energy systems.

In this study, the convergence factor ω is set to 0.001, the maximum number of iterations is set to 30, and the iteration count ρ is initialized to 1.

The equilibrium solution of the game stabilizes after the 5th iteration and converges after the 12th iteration.

The distribution of iteration results is shown in Figure 6.

Figure 7 illustrates the dynamic electricity and heat prices set by the DSO (Distributed Energy Solicitation Office). As a priority for consuming renewable energy within the system, the DSO’s dynamic electricity and heat prices are always within the external energy pricing range, providing users with a more favorable price for energy. Simultaneously, users utilize the DSO’s prices for peak shaving and valley filling, reducing their own energy costs. Furthermore, the higher electricity and heat prices during certain periods in Figure 6 are due to peak user energy consumption. The price cap mechanism ensures that the DSO maximizes its revenue during periods of high demand and high prices. Combined with Table 2, the dynamic electricity pricing strategy effectively improves the overall benefits for both the DSO and users.

• 2. Analysis of Optimization Results of Electricity and Heat Trading among PIES Alliances

Figure 8 shows the optimization results of electricity and heat energy among the PIES alliances. During the periods of 00:00–08:00 and 17:00–24:00, PIES1 and PIES3 have surplus electricity due to low user load and high wind power generation. Therefore, during these periods, the surplus electricity is sold to PIES2, which experiences a nighttime power shortage. From 09:00 to 15:00, due to a sharp increase in daytime electricity demand, PIES1 and PIES3’s own power supply capacity is insufficient to meet their own needs. However, during these periods, PIES2’s photovoltaic power generation is at its peak, producing a large amount of electricity. Through reasonable power dispatch and allocation, PIES2’s surplus electricity can be transmitted to PIES1 and PIES3.

The optimization of thermal energy shows that the system’s internal heat load is mainly supplied by the CHP units. During 00:00–07:00 and 17:00–22:00, PIES1, while meeting its own load demand, can trade surplus heat energy. PIES3 and PIES2 have higher heat loads at night and lower heat loads during the day. Since most of the electricity is supplied by PIES1 during 00:00–08:00 and 17:00–24:00, the output of the CHP units in PIES2 and PIES3 is hindered. Therefore, heat trading with PIES1 is necessary to meet the heat load demand during these periods. During 00:00–07:00 and 17:00–22:00, PIES3 has a higher heat load demand, requiring additional output from the gas-fired boiler to meet the system’s heat energy needs.

The results show that the electricity and heat energy trading mechanism among the alliances effectively balances the electricity supply and demand among the PIES members, improves energy utilization efficiency, and achieves mutual benefit and maximizes the overall system benefits through cooperation.

5.4 Analysis of optimization results

The optimization results for PIES1 are shown in Figure 9. As can be seen from Figure 9, during the periods 1–6 and 23–24, PIES1’s electricity demand is low, and the electricity price set by the DSO is low. PIES1 chooses to purchase electricity from the DSO based on price rather than generating its own. Excess electricity is stored through its own energy storage and sold to PIES2, which has a higher electricity demand. Due to the low electricity price at this time, the heat load is mainly supplied by the heat pump alone. During the periods 7–9, 15–17, and 21–22, PIES1’s electricity and heat load demands are higher than the price set by the DSO. At this time, PIES1’s electricity load is mainly met by gas turbine power generation, which generates heat to meet the heat load demand. Excess electricity generated by the gas turbine is absorbed through thermal energy storage. During the periods of 10:00–14:00 and 18:00–20:00, when electrical and thermal load demands are highest and the DSO sets the highest electricity price, PIES1 still hopes to meet its load demand through its own output. However, due to the limited output of its gas turbines, PIES1 can only meet its own electrical load requirements through inter-PIES1 electricity trading, which increases the trading price accordingly. During this time, PIES1’s energy storage releases heat and discharges. During the most demanding period of 10:00–14:00, PIES1 purchases electricity from PIES2 to meet its load demand. During the less demanding period of 18:00–20:00, PIES1 is willing to sell electricity to PIES2 and PIES3 at a higher price to make a profit. PIES1 improves the overall efficiency of the PIES alliance by purchasing electricity from the DSO during the 02:00–7:00 period to facilitate energy exchange among PIES members, but this results in a loss of efficiency for PIES1 itself. This loss is compensated through cooperative revenue sharing (P2) to ensure that PIES does not suffer losses due to cooperation.

The analysis of PIES2 and PIES3 is similar to that in this stage and will not be repeated. Simulation results show that the hybrid game model proposed in this paper effectively handles the optimal scheduling and revenue sharing of the PIES alliance, achieving cooperative operation of the PIES alliance. The optimization results can be used for daily operation scheduling of the PIES alliance.

5.5 Sensitivity analysis of the combined GCT-CET mechanism

Figure 10 shows the impact of CET price on the revenue of the GCT-CET mechanism.

As shown in Figure 10, the increase in CET price has a positive incentive effect on reducing carbon emissions, especially when the CET price rises from 90 CNY/t to 110 CNY/t, the carbon emission reduction is significant. The total revenue of the system and the revenue of the GCT-CET mechanism also show a steady upward trend with the increase in CET price. Similarly, the GCT price also has a positive incentive effect on the total revenue of the system and the revenue of the GCT-CET mechanism. This is because as the CET price and GCT price continue to rise, the consumption of renewable energy increases, and the curtailment rate of wind and solar power decreases accordingly. Therefore, comprehensively considering the CET price and GCT price can improve the low-carbon economics of the system.

As shown in Figure 11, as the price of green certificates increases, the system’s operating costs gradually decrease and stabilize. At this point, the reward coefficient alters the system’s sensitivity to carbon emissions; a larger reward coefficient increases sensitivity. However, constrained by generator output, the reduction in carbon emissions is slow. Simultaneously, because the system is constantly selling carbon emission rights, a continuously increasing reward coefficient leads to further increases in carbon emissions as the system strives to acquire more carbon emission allowances.

The above analysis shows that the system’s carbon emissions decrease and stabilize as the price of green certificates increases because wind and solar power output has been fully absorbed. Green certificate prices alone cannot reduce carbon emissions. When the carbon trading benchmark price exceeds a certain value, the system’s low-carbon generators are already operating at full capacity, and carbon emissions cannot be reduced further. The compensation coefficient plays a crucial role in green certificate-carbon linkage trading, significantly increasing carbon revenue, greatly motivating enterprises to proactively reduce carbon emissions, and rationally resolving the contradiction between carbon emissions and economic efficiency. Therefore, when formulating the parameters of the green certificate-carbon trading linkage model, the values ​​of green certificate prices, carbon trading benchmark prices, and reward coefficients must be matched to better improve the system’s low-carbon economic efficiency.

6 Conclusion

This paper addresses the diverse demands of multiple stakeholders and, based on master-slave game theory, establishes a dynamic multi-energy pricing model with DSOs as leaders and users as followers. To meet the needs of low-carbon operation of the power system, a CET-GCT collaborative interaction mechanism is proposed. Furthermore, this paper proposes a cooperative and shared model for electricity-heat trading among PIES (Power Systems Integrated Ecosystems). The effectiveness of the model is verified through scenario comparison analysis. The main conclusions are as follows:

• To address the difficulty of simultaneously balancing the interests of PIES and users, a master-slave game-based multi-energy pricing mechanism within PIES is proposed, effectively coordinating the interests of upper-level DSO and lower-level users, and safeguarding the common interests of DSOs and users.
• Utilizing the CET-GCT collaborative mechanism for simultaneous collaborative optimization of the source and load sides, results show that while optimizing the energy supply and demand balance, system carbon emissions are reduced and DSO revenue is increased.
• Based on the master-slave game within PIES and with users, a hybrid game model considering electricity-heat trading is constructed. Further consideration of carbon trading can further optimize system energy supply, increase system revenue, and promote the low-carbon economic operation of the system.

Generative AI statement

The author(s) declared that generative AI was not used in the creation of this manuscript.

Any alternative text (alt text) provided alongside figures in this article has been generated by Frontiers with the support of artificial intelligence and reasonable efforts have been made to ensure accuracy, including review by the authors wherever possible. If you identify any issues, please contact us.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Footnotes

Edited by: Hao Wang, Monash University, Australia

Reviewed by: Yongsheng Qiao, Taiyuan University, China

Boyang Qu, Liaoning Institute of Science and Technology, China