Dynamic Calculation Method for Zonal Carbon Emissions in Power Systems Based on the Theory of Production Simulation and Carbon Emission Flow Theory

Abstract

Power systems are the main source of carbon emissions. Currently, coordinated operation strategies of the source–grid–load–storage model considering carbon emissions is primarily expanded from the generation side. For practical power systems, where multiple types of generating units coexist at a single node, it is difficult to develop unit combination strategies that simultaneously consider carbon emission factors and power flow constraints. Therefore, a new power flow calculation method based on connectivity matrix theory was proposed, aiming to address the issues of existing approaches that are too coarse and unable to accurately represent the operating states of multiple units under each node. Furthermore, a new method for dynamic calculation of regional carbon emission based on connectivity matrix and carbon emissions flow was introduced to improve the accuracy of carbon emission measurements. Firstly, a simulation model for a coordinated optimization operation based on the minimum system cost for the source–grid–load was established and an optimal flow calculation method using a connectivity matrix was introduced. Second, a dynamic carbon emission calculation method, considering electricity sources, was developed by combining the results of the optimal power flow calculation with carbon emission flow theory. Finally, the effectiveness of the approach in this article was verified by the IEEE 14-bus system example and a provincial power grid, ensuring strict adherence to the conservation principle of carbon emissions between the supply and demand sides and satisfying power flow constraints.

1. Introduction

Most countries worldwide have reached a consensus on achieving ‘carbon peaking and carbon neutrality’ [1]. China has proposed accelerating the construction of a new type of power system to support the realization of the ‘carbon peaking and carbon neutrality’ goals. In the process of building the new power system, the focus is on increasing the proportion of clean energy and improving system efficiency to reduce carbon dioxide emissions from the source and process [2]. By the end of 2023, China’s cumulative installed capacity of grid-connected wind power and photovoltaic generation reached 1.05 billion kilowatts, accounting for 36% of the total installed capacity [3]. The scale of renewable energy generation has historically exceeded thermal power, with the annual increase exceeding half of the global total. However, with the promotion of electrification, it is estimated that by 2030, the proportion of load demand from end-use sectors such as industry, construction, and transportation through electrification will reach approximately 35% [4]. The power system load will continue to increase, and the task of achieving ‘carbon peak and carbon neutrality’ remains challenging [5].

The efficiency and accuracy of current carbon emission accounting techniques need improvement, and there is a lack of mature standards and unified measurement methods for carbon emission certification. Currently, carbon emission accounting methods mainly include macroscopic statistical methods [6], emission factor methods [7], mass balance methods [8], online monitoring methods [9], life cycle assessment methods [10], and carbon emission flow methods [11,12]. The carbon emission statistical accounting process highly relies on original production data and is currently limited to annual statistical calculations. The macroscopic statistical methods do not consider the specific technical processes of power system generation and require long-term fuel consumption statistics as data support, making it difficult to achieve detailed analysis of carbon emissions in power systems. The life cycle assessment methods analyze carbon emissions from a timeline perspective but lack consideration of the physical characteristics of power systems and cannot clarify the spatial–temporal transfer mechanisms of carbon emissions in the power system [13]. Emission factor methods, mass balance methods, and life cycle assessment methods are influenced by the statistical cycles of basic production data, such as fuel and raw materials, typically conducted on an annual basis. Online monitoring methods are currently the only means for continuous carbon emission monitoring, but due to economic and data reliability reasons, they have not been widely applied domestically. Carbon emission flow methods consider the transfer network of carbon emissions, bridging the gap between production and consumption, and can clarify carbon emission responsibilities [14,15]. However, they are still in the preliminary research stage and require further exploration. Additionally, the current carbon emission measurement methods cannot provide real-time guidance for the scheduling and operation of power plant units within the system.

To refine the strategy for a low-carbon target in the source–grid–load–storage system during its operational phase, the current simulation calculations that consider carbon emission flows [16] typically focus on several optimization objectives. These include minimizing system operation and maintenance costs [17,18], adhering to low carbon and system reliability constraints [19], reducing power losses [20], ensuring power balance [21], managing power flow [22], maintaining adequate system spinning reserves [23], optimizing power generation [24], and integrating energy storage systems [25]. However, these approaches often overlook the adequacy of system generation capacity. In addition, a coordinated low-carbon economic dispatch strategy is proposed in the article [26]. This strategy considers a dual-layer incentive-based carbon trading mechanism and provides a comprehensive energy optimization scheme for multi-agent systems with the park as the main agent. However, it does not guide the coordinated output of various types of power sources at individual nodes. Another paper [27] establishes a two-stage low-carbon dispatch model based on super-carbon demand response, which can effectively reduce system carbon emissions by shifting loads from high-carbon periods. However, it lacks optimization of the start-up combination of traditional units. Additionally, a power system carbon emission flow calculation method based on the power flow distribution matrix in the article [22] quantifies the carbon emissions sources and composition at each node. However, it lacks coordination optimization for the entire source–grid–load process, which prevents guidance on the coordinated output of various types of power sources. Paper [28] considers incentive loads, joint operation modes, and gas network characteristics to establish a low-carbon economic dispatch model for the electric–gas interconnected system, minimizing daily operating costs. While it provides the total output of multiple power sources at the system level. The feature of carbon emission calculation methods in other different literature is shown as

Table 1. Carbon emission calculation methods in different literature.

Literature	Simulation Object	Simulation Method	Merits	Deficiencies
[29]	Coupled transportation and power distribution networks	Multi-objective optimization	Achieve an optimal balance between economic efficiency and emission reduction	Without considering the carbon emission flow
[30]	Consumer-driven carbon intensities in power system	Carbon emission flow calculation method	Adaptive regression-based calculation framework	Unable to coordinate the output of various power sources
[31]	Carbon efficient power grid	The uplift payment scheme.	Enjoy the advantage of maximal flexibility	Lack of constraints on the actual output status of power generation units
[32]	Demand-side distribution networks	A spatial-temporal carbon response model based on geographically dispatchable loads	Emission target is embedded in the objective function via the exterior point method	Unable to coordinate the output of various power sources
[33]	Distribution networks considering the carbon emission allowance on the demand side	A multiagent-based bi-level operation framework	Transfer the responsibility for carbon emissions from generation to demand	Without considering the carbon emission flow
[34]	All datacenter clusters across Google’s fleet	A suite of analytical pipelines used to gather the next day’s carbon intensity forecasts	Minimizes electricity-based carbon footprint and by delaying temporally flexible workloads	Unable to coordinate the output of various power sources
[35]	A low-carbon optimal scheduling model with demand response-based carbon intensity control	The carbon emission flow model	The storage emission and relevant dynamic impacts on the system are fully considered.	Without considering the renewable energy
[36]	The direct current transmission system	Direct current dynamic optimal power flow with carbon emission trading	Presents a distributed alternating direction method of multipliers approach	Without considering the renewable energy
[37]	Combined renewable and coal power system	Multi-objective optimization framework	Consider the carbon capture systems	Without considering the carbon emission flow
[38]	Power system from source–grid–load–storage	Bi-level alternating optimal scheduling model	Integrate the carbon emission flow with the power flow	Unable to coordinate the output of various power sources
[39]	Power system and demand-side response	Two-stage low-carbon optimization scheduling model	Stimulate the change in electricity consumption behavior on the demand side to reduce the carbon emissions.	Unable to coordinate the output of various power sources
[40]	Power system and demand-side response	Load aggregators participate in the demand response market framework	Meet the target of minimum carbon emissions with the lowest cost of power grid operation and best sharpening effect	Without considering the carbon emission flow

. It can be seen that none of the mentioned studies offer guidance on the coordinated output combination for various types of units at individual nodes. Optimizing the coordinated output of diverse unit combinations directly impacts the size of carbon emission sources and is crucial for achieving safe and low-carbon operation in modern power systems. Therefore, analyzing the coordinated output combination of multiple types of units under each node is crucial.

In response to the above issue, this paper incorporates various types of power generation units under the nodes into actual operational simulations using a connectivity matrix. By combining optimal power flow and carbon emission flow theories, the paper provides recommendations for combining multiple types of units under a low-carbon objective. Firstly, an optimal power flow calculation method (OPFCM) based on the connectivity matrix is proposed, and a simulation model for coordinated optimization of source–grid–load–storage is established based on minimizing system costs. Secondly, a dynamic carbon emission calculation method is introduced based on carbon emission flow theory and OPFCM. Finally, an in-depth analysis is conducted on the accuracy of system flow, reasonableness of unit output, and carbon emission conservation in both the IEEE 30-bus system and a provincial-level power grid system, validating the effectiveness of the proposed methods.

2. Calculation Method of Power Flow and Carbon Emission Flow Based on Chain Matrix

The calculation method for carbon emission flow in the power system [16] provides a comprehensive approach to our study [16]. Through the operating simulation calculations, the output of various types of generating units at each moment can be obtained, resulting in the branch power flow distribution matrix $P_B=\left(P_{Bij}\right)N\times N$, the generator injection distribution matrix $P_G=\left(P_{Gkj}\right)K\times N$, the load distribution matrix $P_L=\left(P_{Lmj}\right)M\times N$, the demand response distribution matrix $P_O={\left(P_{Oij}\right)}_{M\times N}$, the load shedding distribution matrix $P_A={\left(P_{Aij}\right)}_{M\times N}$, and the renewable energy curtailment distribution matrix $P_R={\left(P_{Rlj}\right)}_{L\times N}$. To calculate the carbon flow and carbon potential introduced by each branch, the branch power flow distribution matrix undergoes element filtering, only retaining the power of the branches with power out flow. That is, when the power flow from node i to node j, P_Bij = p, P_Bji = 0.

The calculation of the node active flux matrix involves assessing the electric energy sources flowing into each node at various time intervals, along with their corresponding carbon emission factors. A node active power flux matrix is established to represent the sum of active power flowing into each node at each moment, which is an N-dimensional square matrix, denoted as $P_N=\left(P_{Nij}\right)N\times N$, specifically, it is a diagonal matrix expressed as $P_{Nij}=0,\left(i\ne j\right)$, whose diagonal elements calculation methods are as follows:

$$P_{Nii}={\displaystyle \sum _{n=1}^N{P}_{B_{n,i}}}+{\displaystyle \sum _{k=1}^K{P}_{G_{k,i}}}$$

(1)

The calculation of the generator carbon emission factor matrix: assuming that the carbon emission factor for the $k_{th}$ generator($k=1,2,\dots ,K$) is denoted as $e_{Gk}$, the generator carbon emission factor matrix can be represented as follows:

$$E_{\mathrm{G}}={\left[\begin{array}{cccc}{e}_{G1}& e_{G2}& \cdots & e_{Gk}\end{array}\right]}^{\mathrm{T}}$$

(2)

The calculation of the node carbon potential matrix: defining the carbon potential at each node as the average carbon emission obtained by distributing all carbon emissions flowing into that node over the total power flowing into the node. The node carbon potential matrix, denoted as $E_N={\left(e_{Ni}\right)}_{N\times 1}$, is an N × 1 matrix, and its elements are calculated as follows:

$$e_{Ni}={\displaystyle \frac{{\displaystyle \sum _{n=1}^N{p}_{B_{n,i}}{e}_{Nn}}+{\displaystyle \sum _{k=1}^K{p}_{G_{k,i}}{e}_{Gk,i}}}{P_{Nii}}}$$

(3)

where $p_{B_{n,i}}$ represents the flow from node n to node i, and $e_{Nn}$ denotes the average carbon potential at node n. Additionally, $p_{G_{k,i}}$ represents the $k_{th}$ generator connected to node i, and $e_{Gk,i}$ corresponds to the carbon emission factor associated with that generator.

The calculation of the node carbon emission matrix aims to characterize the carbon emissions corresponding to the electricity consumption of loads at each node. The node carbon emission matrix, denoted as $E_T={\left(e_{Ti}\right)}_{M\times 1}$, is an M × 1 matrix, expressed as follows:

$$E_T=P_L{E}_N$$

(4)

3. Simulation Models for New-Generation Power Systems Operation

Assuming a system with N nodes, including K power nodes (comprising generalized power sources such as coal-fired, gas-fired, wind, photovoltaic, and energy storage and so on), M load nodes, and L renewable energy nodes that experience power curtailment (primarily wind power and photovoltaic), the calculation flowchart for the carbon emission factors and carbon emissions in a multi-node system is illustrated as Figure 1:

3.1. Establishment of the Target Function

In the operational simulation system of the new-type power system, it is essential to consider the economy, low carbon, and reliability of the system’s operation. Therefore, considering the system’s maintenance costs, carbon emission costs, and the cost of power shortages, the objective function for the operational simulation model is established as follows [17]:

$$\min f=C^{Op}+C^{C{O}_2}+C^A$$

(5)

where $C^{Op},C^{C{O}_2},C^A$ donate to the system’s maintenance costs, system’s carbon emission costs, and the cost of power shortages.

3.1.1. System’s Maintenance Costs

The calculation methods of system’s maintenance costs are as follows [18]:

$$\begin{array}{cc}{C}^{Op}\hfill & =C^{\mathrm{OM}}+C^{\mathrm{F}\mathrm{U}\mathrm{E}\mathrm{L},G}+C^{\mathrm{S}\mathrm{U},G}+C^{OP,\mathrm{O}}+C^{OP,\mathrm{S}}\hfill \\ & ={\displaystyle \sum _{g\in \psi }{c}_g^{om}{C}_g}+{\displaystyle \sum _{o\in \Omega }{c}_o^{om}{C}_o}+\hfill \\ & {\displaystyle \sum _{s\in E}\left(c_s^{om,cap}{E}_s^{cap}+c_s^{\mathrm{o}\mathrm{m},\mathrm{e}\mathrm{n}}{E}_s^{en}\right)}+\hfill \\ & {\displaystyle \sum _{t\in T}{\omega }_t}\left[\begin{array}{l}{\displaystyle \sum _{g\in \psi }\left(a_g{P}_{g,t}^2+b_g{P}_{g,t}+c_g\right)U_{g,t}}+\\ \left(c_g^{on}{S}_{g,t}+c_g^{off}{D}_{g,t}\right)P_{g,t}\end{array}\right]+\hfill \\ & {\displaystyle \sum _{t\in T}{\omega }_t}{\displaystyle \sum _{s\in E}{c}_s^{\mathrm{p}\mathrm{r},\mathrm{o}}{p}_{o,t}}+{\displaystyle \sum _{t\in T}{\omega }_t}{\displaystyle \sum _{s\in E}{c}_s^{\mathrm{p}\mathrm{r},\mathrm{e}\mathrm{n}}{E}_{s,t}^{en}\left(1-{\eta }_s\right)}\hfill \end{array}$$

(6)

where C^OM represents the fixed operation and maintenance costs required by the equipment loss and normal operation of the generating units, demand-side response, and energy storage; $C^{\mathrm{F}\mathrm{U}\mathrm{E}\mathrm{L},G}$ denotes the generation cost of the units; $C^{\mathrm{S}\mathrm{U},G}$ indicates the start-up and shut-down costs of the units; $C^{OP,\mathrm{O}}$ is the operational cost for demand-side response; and $C^{OP,\mathrm{S}}$ signifies the operational cost of the energy storage. $c_g^{om}$, $c_o^{om}$, $c_s^{om,cap}$, and $c_s^{\mathrm{o}\mathrm{m},\mathrm{e}\mathrm{n}}$, respectively, represent the operation and maintenance costs corresponding to the unit capacity of conventional units, the unit capacity of demand-side response, the unit power, and the unit capacity for energy storage. ${\omega }_t$ denotes the operating duration of the corresponding entity. $a_g$, $b_g$, and $c_g$ are the generation cost coefficients for conventional units. $P_{g,t}$ is the power generated by conventional units at time t, $U_{g,t}$ represents the start-up and shut-down state of the units, where $U_{g,t}$ = 1 indicates that the unit is operational at time t, and $U_{g,t}$ = 0 indicates that the unit is shut down at time t. $c_g^{on}$ and $c_g^{off}$, respectively, represent the start-up and shut-down costs per unit capacity for conventional units. $S_{g,t}$ and $D_{g,t}$, respectively, denote the start-up and shut-down actions of the units at time t, where $S_{g,t}$ = 1 indicates that the unit starts up at time t. $D_{g,t}$ = 1 indicates that the unit shuts down at time t, and $S_{g,t}$ = 0 or $D_{g,t}$ = 0 indicates that there is no change in the state of the unit at time t. $c_s^{\mathrm{p}\mathrm{r},\mathrm{o}}$ represents the unit power operation and maintenance cost for demand-side response at time t, which is mainly reflected in the unit power opportunity cost loss due to the inability to participate in other auxiliary service markets after participating in the electricity quantity market balance. $p_{o,t}$ is the adjustment power of demand-side response at time t, and $c_s^{\mathrm{p}\mathrm{r},\mathrm{e}\mathrm{n}}$ denotes the unit power cost of energy storage at time t, primarily reflecting the unit power cost caused by charge–discharge efficiency loss during energy storage discharge, which can be characterized by the average electricity price during energy storage charging. $E_{s,t}^{en}$ is the discharge power of energy storage at time t, and ${\eta }_s$ represents the charge–discharge efficiency of energy storage, which can generally be taken 86% as a certain value [18].

3.1.2. Annual Carbon Emission Cost of System

The annual carbon emission cost of the system is mainly related to the power generation of conventional units such as coal power and gas power in the system. The calculation method for the annual carbon emission cost of the system is as follows [19]:

$$C^{C{O}_2}=c^{C{O}_2}{\displaystyle \sum _{t\in T}{\omega }_{g,t}}{\displaystyle \sum _{g\in {\psi }_T}{e}_g{P}_{g,t}}$$

(7)

where $c^{C{O}_2}$ represents the penalty cost corresponding to the unit carbon emission, which can be characterized by the carbon trading price. For reference, the closing price of the national carbon market on 20 June 2022, was CNY 59.25 per ton CO₂. ${\omega }_{g,t}$ denotes the operating duration of unit g, and $e_g$ represents the carbon emission factor of unit g. The typical power generation carbon emission factors for various types of units are shown in

Table 2. The typical carbon emission factors for various types of power generation units.

Types of Power Generation	(tCO2/MWh)
Ultra-Supercritical 1000 MW Coal-Fired Power Unit	0.7938
Ultra-Supercritical 600 MW Coal-Fired Power Unit	0.8067
Supercritical 600 MW Coal-Fired Power Unit	0.8385
Subcritical 300 MW Coal-Fired Power Unit	0.8875
Subcritical 300 MW Circulating Fluidized Bed Coal-Fired Power Unit	0.9171
Ultra-High Voltage and Below Units	0.9363
Gas Turbine Unit	0.3789
Biogas Power Generation	0.23
Hydro Power Generating Units	0.0021
Nuclear Power Units	0.0065

[41,42].

3.1.3. The Cost of Power Shortage Loss

The cost of power shortage loss is mainly caused by the system’s inability to guarantee the basic electricity demand of users due to emergencies. The calculation method is as follows [20]:

$$C^A=c^A{\displaystyle \sum _{i\in N}{\displaystyle \sum _{t\in T}{\displaystyle \sum _{m\in M}{\omega }_{A_{\mathrm{m},i,t}}{P}_{A_{\mathrm{m},i,t}}}}}$$

(8)

where $c^A$ represents the penalty cost per unit of unserved energy. It can be initially characterized by the degree of electricity GDP, which is obtained by combining the annual GDP of the region and the total electricity consumption of the whole society. In this context, it is taken as 16.41 CNY/kWh. ${\omega }_{A_{\mathrm{m},i,t}}$ represents the duration of power outage for the load m connected through node i. $P_{A_{\mathrm{m},i,t}}$ represents the power outage of the load m connected through node i at time t.

3.2. Constraint Conditions

3.2.1. Power Balance Constraint

Considering that all nodes need to satisfy the power balance constraint at every moment, the power balance constraint of each node can be described as follows: at any moment, the power flowing into node i + the power generation of the generator connected to node i + the load response power of node i + the power shortage of the load at node i = the load power connected to node i + the power flowing out from node i + the discarded power of the renewable energy connected to node i. The mathematical expression is as follows [21]:

$${\displaystyle \sum _{n=1}^N{P}_{B_{n,i}}}+{\displaystyle \sum _{k=1}^K{\eta }_{G_{k,i}}{P}_{G_{k,i}}}+{\displaystyle \sum _{m=1}^M{P}_{O_{m,i}}}+{\displaystyle \sum _{m=1}^M{P}_{A_{m,i}}}={\displaystyle \sum _{m=1}^M{P}_{L_{m,i}}}+{\displaystyle \sum _{l=1}^L{P}_{R_{l,i}}}$$

(9)

where $P_{B_{n,i}}$ represents the power flowing into node i from node n. $P_{G_{k,i}}$ donates the power generated by generator k connected at node i. ${\eta }_{G_{k,i}}$ refers to the adjustable depth of generator k connected at node i. According to the “Notice on Carrying out the Upgrading of Coal-fired Power Units Nationwide” issued by the National Development and Reform Commission and the National Energy Administration in 2021, the adjustable depth of coal-fired power units is taken as 0~35%; $P_{O_{m,i}}$ represents the response power of demand response m connected at node i, $P_{A_{m,i}}$ donates to the power shortage of load m connected at node i, $P_{L_{m,i}}$ refers to the power consumption of load m connected at node i, $P_{B_{i,n}}$ represents the power flowing out of node i to node n, and $P_{R_{l,i}}$ donates to the wasted power of renewable energy l connected at node i. The above system balance matrix can be further expanded and applied to calculate the power balance of each node in a multi-node system. For a multi-node system, a branch power flow distribution matrix, a unit injection distribution matrix, a load distribution matrix, a demand-side response distribution matrix, a load power shortage distribution matrix, and a renewable energy waste distribution matrix can be established, as shown below [15]:

(1)

Branch Power Flow Distribution Matrix

First, establish the branch power flow distribution matrix, which is an N-dimensional square matrix denoted as $P_B={\left(P_{Bij}\right)}_{N\times N}$; since the power flow between non-adjacent nodes is 0, the optimization result cannot consider the connection relationship between nodes, and there may be power flow between non-adjacent nodes. Therefore, a connectivity matrix needs to be defined and multiplied with the pre-optimized branch power flow distribution matrix to eliminate the power flow between non-adjacent nodes, and only the power flow between adjacent nodes is retained for the model optimization solution. The branch power flow connectivity matrix is also an N-order square matrix, denoted as $A_B={\left(A_{Bij}\right)}_{N\times N}$, if there is a branch connected between node i and node j (i, j = 1, 2, …, N); hence, $A_{Bij}$ = 1, $A_{Bji}$ = −1, and the branch power flow distribution matrix is equal to the product of the branch current connectivity matrix multiplied by the branch current pre-distribution matrix. The mathematical expression is as follows [16]:

$$P_B=A_B\odot P_{preB}$$

(10)

where $P_{preB}$ is the N × N branch power flow pre-distribution matrix obtained by system optimization, and $P_B$ is the branch power flow distribution matrix. If there is a branch connected between node i and node j (i, j = 1, 2, …, N), and the positive active current p flows from node i to node j through this branch, then $P_{Bij}$ = p, $P_{Bji}$ = −p; if the active current p flowing through this branch is a reverse current, then $P_{Bij}$ = −p, $P_{Bji}$ = p; in other cases, $P_{Bij}$ = $P_{Bji}$ = 0, specifically, for all diagonal elements, $P_{Bii}$ = 0(i, j = 1, 2, …, N).

(2)

Unit Injection Distribution Matrix

The unit injection distribution matrix is a K × N order matrix, denoted by $P_G={\left(P_{Gkj}\right)}_{K\times N}$. If the k_th (k = 1, 2, …, K) power generating unit is connected to node j, and the active current injected from the k_th node containing a generator to node j is p, then $P_{Gkj}$ = p; otherwise, $P_{Gkj}$ = 0. Since the unit injection pre-distribution matrix obtained by optimization may yield output values at nodes where generator unit k is not connected, it is necessary to correct the unit injection pre-distribution matrix with the unit connectivity matrix. The unit connectivity matrix is denoted by $A_G={\left(A_{Gkj}\right)}_{K\times N}$, and if generator unit k is connected to node j (j = 1, 2, …, N), then $A_{Gkj}$ = 1. The relationship between the unit injection distribution matrix and the connectivity matrix is as follows:

$$P_G=A_G\odot P_{preG}$$

(11)

Since one unit can only inject into one node, there is the following constraint for any row of the unit injection distribution matrix: only one element in $P_{Gij}$(j = 1, 2, …, N) is non-zero; that is, the element corresponding to the column of the access point of unit i is non-zero, and the elements of the other N-1 columns are zero. The mathematical expression is as follows:

$${\displaystyle \sum _{j=1}^N{A}_{Gkj}}=1,\forall k$$

(12)

where $A_{Gkj}$ represents the element in the k_th row and j_th column of the unit connectivity matrix; if generator unit k is connected to node j (j = 1, 2, …, N), then $A_{Gkj}$ = 1; otherwise, it is 0.

(3)

Load Distribution Matrix

The load distribution matrix is an M × N order matrix denoted by $P_L={\left(P_{Lmj}\right)}_{M\times N}$ If node j is the m_th (m = 1, 2, …, M) node with load, and the active load is p, then $P_{Lmj}$ = p; otherwise, $P_{Lmj}$ = 0. Since the load data are an input boundary condition, there is no need to constrain them with a connectivity matrix. As one load can only be injected into one node, there is the following constraint for any row in the load distribution matrix, which means that only one element in $P_{Lij}$(j = 1, 2, …, N) is non-zero, the element at the access point of load i is non-zero, and the other N-1 elements are all zero.

(4)

Demand-side Response Distribution Matrix

The demand-side response distribution matrix is an M × N order matrix denoted by $P_O={\left(P_{Oij}\right)}_{M\times N}$. If node j is the m_th (m = 1, 2, …, M) node with demand-side response, and the demand-side response adjustment power is p, then $P_{Amj}=p$; otherwise, it is $P_{Amj}=0$. Since the demand-side response pre-distribution matrix obtained through optimization may result in output values in nodes where demand-side response m is not connected, it is necessary to correct the demand-side response pre-distribution matrix with the load connectivity matrix. The load connectivity matrix is denoted by $A_O={\left(A_{Omj}\right)}_{M\times N}$, and if demand-side response m is connected to node j (j = 1, 2, …, N), then $A_{Omj}$ = 1. The relationship between the demand-side response distribution matrix and the connectivity matrix is as follows:

$$P_O=A_O\odot P_{preO}$$

(13)

Since one demand-side response can only be connected to one node, there is the following constraint for any row of the demand-side response distribution matrix: only one element in $P_{Omj}(j=1,2,\cdots N)$ is non-zero; that is, the element corresponding to the column of the access point of demand-side response m is non-zero, and the elements of the other N-1 columns are all zero. The mathematical expression is:

$${\displaystyle \sum _{j=1}^N{A}_{Omj}}=1,\forall m$$

(14)

where $A_{Omj}$ represents the connection relationship between demand-side response m and node j (j = 1, 2, …, N); if demand-side response m is connected to node j (j = 1, 2, …, N), then $A_{Omj}$ = 1; otherwise, it is 0.

(5)

Load Power Shortage Distribution Matrix

The load power shortage distribution matrix is an M × N order matrix denoted by $P_A={\left(P_{Aij}\right)}_{M\times N}$. If node j is the m_th (m = 1, 2, …, M) node experiencing a power shortage, and the power shortage is p, then $P_{Amj}=p$; otherwise, it is $P_{Amj}=0$. Due to the possibility of obtaining power outage values in nodes that are not connected to load m in the pre-distribution matrix of load power outage obtained through optimization, it is necessary to correct the load power shortage pre-distribution matrix using the load connectivity matrix. The relationship between the load power shortage distribution matrix and the connectivity matrix is as follows:

$$P_A=A_O\odot P_{preA}$$

(15)

(6)

Renewable energy curtailment distribution matrix

The renewable energy curtailment distribution matrix is an L × N order matrix denoted by $P_R={\left(P_{Rlj}\right)}_{L\times N}$. If node j is the l_th (l = 1, 2, …, L) node with renewable energy curtailment, and the curtailment power is p, then $P_{Rlj}=p$; otherwise, it is $P_{Rlj}=0$. Since the renewable energy curtailment pre-distribution matrix obtained by optimization may yield curtailment values at nodes where renewable energy l is not connected, it is necessary to correct the renewable energy curtailment pre-distribution matrix with the renewable energy connectivity matrix. The renewable energy connectivity matrix is denoted by $A_R={\left(A_{Rlj}\right)}_{L\times N}$, and if renewable energy l is connected to node j (j = 1, 2, …, N), then $A_{Rlj}$ = 1. The relationship between the demand-side response distribution matrix and the connectivity matrix is as follows:

$$P_R=A_R\odot P_{preR}$$

(16)

Since a renewable energy plant can only curtail at the node it is connected to, there is the following constraint for any row of the renewable energy curtailment distribution matrix: only one element in $P_{Rlj}$(j = 1, 2, …, N) is non-zero, meaning the element corresponding to the column of the access point of renewable energy l is non-zero, and the elements of the other N − 1 columns are zero. The mathematical expression is as follows:

$${\displaystyle \sum _{j=1}^N{A}_{Rlj}}=1,\forall l$$

(17)

where $A_{Rlj}$ represents the connection relationship between the renewable energy source l and node j (j = 1, 2, …, N); if the renewable energy source l with curtailment is connected to node j (j = 1, 2, …, N), then $A_{Rlj}$ = 1; otherwise, it is 0.

3.2.2. System Power Flow Constraints

To further introduce the direct current flow to describe the flow constraints between nodes [22], the system node matrix ${\theta }_L={\left({\theta }_{Li}\right)}_{N\times 1}$ and the system impedance matrix $B_L={\left(B_{Lij}\right)}_{N\times N}$ are established. ${\theta }_{Li}$ represents the phase angle deviation of node i relative to the balance node, and $B_{Lij}$ characterizes the reciprocal of the admittance between node i and j. Then, the system flow constraints are represented as follows:

$$P_B=B_L\times {\theta }_L$$

(18)

3.2.3. System Spinning Reserve Constraints

Considering the uncertainty of wind power and photovoltaic output, the system needs to consider satisfying the load demand with non-renewable energy sources during unit dispatch optimization. Therefore, the startup of non-renewable energy sources in the system need to meet the following constraint [23]:

$${\displaystyle \sum _{k=1}^K{P}_{G_{k,i}}}\ge \left(1+{\eta }_L\right)\left({\displaystyle \sum _{m=1}^M{P}_{L_{m,i}}}+{\displaystyle \sum _{l=1}^L{P}_{R_{l,i}}}\right)$$

(19)

where $P_{G_{k,i}}$ represents the maximum power generation capacity of generator unit k connected at node i, ${\eta }_L$ refers to the system load reserve rate, which is generally set at 2% to 3%. $P_{L_{m,i}}$ donates the power consumption of load m connected at node i, and $P_{R_{l,i}}$ represents the curtailed power of renewable energy source l connected at node i.

3.2.4. Constraints on Power Generation Technology

The constraints on power generation technology mainly include the output constraints of renewable energy units, the technical constraints of the conventional units’ output, and so on. For renewable energy units (such as wind power and photovoltaic), their operational constraints primarily depend on the availability of resources, and their maximum dispatchable capacity is represented as the product of the maximum capacity factor and the installed capacity. The above constraints are expressed as follows [24]:

$${\gamma }_g^{\min }{U}_{g,t}{C}_g\le P_{g,t}\le {\gamma }_g^{\max }{U}_{g,t}{C}_g,\forall g\in {\Psi }_T,t$$

(20)

$$\left\{\begin{array}{l}{P}_{g,t}-P_{g,t-1}\le {\rho }_g^{RU}{C}_g{U}_{g,t}+\\ \max ({\gamma }_g^{\min }{C}_g,{\rho }_g^{RU}{C}_g)S_{g,t},\forall g\in {\Psi }_T,t\\ P_{g,t-1}-P_{g,t}\le {\rho }_g^{RD}{C}_g{U}_{g,t}+\\ \max ({\gamma }_g^{\max }{C}_g,{\rho }_g^{RU}{C}_g){\mathrm{D}}_{g,t},\forall g\in {\Psi }_T,t\end{array}\right.$$

(21)

$$U_{g,t}-U_{g,t-1}=S_{g,t}-{\mathrm{D}}_{g,t},\forall g\in {\Psi }_T,t$$

(22)

$$U_{g,t},S_{g,t},{\mathrm{D}}_{g,t}\in \left\{0,1\right\},\forall g\in {\Psi }_T,t$$

(23)

$$U_{g,t}\ge {\displaystyle \sum _{\tau =t-m_g^{on}}^t{S}_{g,\tau }},\forall g\in {\Psi }_T,t$$

(24)

$$1-U_{g,t}\ge {\displaystyle \sum _{\tau =t-m_g^{off}}^t{\mathrm{D}}_{g,\tau }},\forall g\in {\Psi }_T,t$$

(25)

$$0\le P_{g,t}\le {\alpha }_{g,t}{C}_g,\forall g\in {\Psi }_R,t$$

(26)

where ${\Psi }_T$ represents the set of conventional units, including coal-fired and gas-fired units, and ${\Psi }_R$ refers to the set of renewable energy units, including wind power and photovoltaic generation units; ${\gamma }_g^{\min }$ and ${\gamma }_g^{\max }$, respectively, donate the minimum and maximum technical output levels of the units; ${\rho }_g^{RU}$ and ${\rho }_g^{RD}$, respectively, represent the maximum upward and downward ramping capabilities of the units, excluding the minimum output level required when the unit transitions from shutdown state to startup state; $m_g^{on}$ and $m_g^{off}$, respectively, refer to the minimum continuous operating and shutdown time of the units; ${\alpha }_{g,t}$ is the predicted energy factor of the renewable energy units. Constraint (20) represents the upper and lower limits of the output of conventional units; constraint (21) represents the ramping constraints of the units, with the downward ramping constraint being active when the unit transitions from a startup state to a shutdown state, and vice versa for the upward ramping constraint; constraints (22) and (23) indicate that the operating state and state transition flag of the units are both 0–1 variables, and restrict the units to have only one operating state at any given time; constraint (24) limits the minimum continuous startup time of the units; constraint (25) limits the minimum continuous shutdown time of the units; and constraint (26) limits the maximum and minimum output of renewable energy units.

3.2.5. Constraints of an Energy Storage System

The charging and discharging process of the energy storage system (ESS) must consider physical constraints such as energy storage charging and discharging efficiency, power constraints, and rated capacity constraints [25].

$$S_{s,t+1}=S_{s,t}+\Delta t{P}_{s,t}^c{\eta }_s^c-{\displaystyle \frac{\Delta t{P}_{s,t}^d}{{\eta }_s^d}}$$

(27)

$$0\le P_{s,t}^c,P_{s,t}^d\le E_s^{cap}$$

(28)

$$E_{s,\min }^{en}\le S_{s,t}\le E_{s,\max }^{en}$$

(29)

$$S_{s,1}=S_{s,\mathrm{T}}$$

(30)

where $S_{s,t}$ represents the energy state of the energy storage; $\Delta t$ refers to the time interval. Constraint (27) represents the energy evolution process of the ESS. Constraint (28) represents the power constraints for charging and discharging of the ESS. Constraint (29) represents the upper and lower capacity constraints of the ESS. Constraint (30) ensures that the initial and final energy states of the ESS are equal over the optimization period [1,T].

3.2.6. System Generation Adequacy Constraint

The adequacy of a power system refers to its ability to meet user demand in terms of power generation, transmission, and supply after considering factors such as the outage or underperformance of some equipment, which characterizes the steady-state performance of the system. In the planning of traditional power sources, the planning target is usually set at a certain percentage higher than the annual peak load [43]. However, due to the volatility and randomness of renewable energy output, it does not possess the strong controllability of traditional thermal power units, and its generation adequacy differs significantly from conventional thermal power units. Traditional power system planning often neglects or simplifies the capacity value of renewable energy, leading to underinvestment or overinvestment in scenarios with a high proportion of renewable energy. Therefore, scholars have proposed the concept of capacity credit [44] (expressed as a percentage), which represents the proportion of wind and photovoltaic units that can be considered equivalent to conventional units under the premise of equal reliability, enabling the comparison of uncontrollable, volatile, and random generation technologies with conventional units on the same level. The system generation adequacy constraint is expressed as follows [45]:

$${\displaystyle \sum _{g\in \Psi }{c}_g^{ad,G}}{C}_g+{\displaystyle \sum _{s\in \mathrm{E}}{c}_s^{ad,S}}{E}_s^{cap}\ge (1+{\alpha }_R)\max \left({\displaystyle \sum _{z\in \mathrm{Z}}{D}_{z,t}}\right)$$

(31)

In the formula, $c_g^{ad,G}$ is the capacity credit coefficient of generators; $c_s^{ad,S}$ is the capacity credit coefficient of energy storage, where the credible capacity of energy storage is related to the rated discharge duration. Energy storage can be used to replace traditional peaking power units, and its capacity credit is related to the maximum discharge duration of the energy storage. According to years of historical operational data from the California electricity market in the United States, the capacity credit of energy storage batteries with a discharge duration of 4 h is approximately about 75% [46]. ${\alpha }_R$ is the system generation adequacy coefficient. Considering the possible faults or maintenance of the generator units in the system, the actual planned generation capacity should exceed a certain proportion of the system’s peak load. $D_{z,t}$ represents the power demand of important load z at time t, and Z represents the set of important loads.

3.2.7. System Reliability Constraint

The system reliability constraint is described through the load capacity shortage rate, ensuring that in emergency scenarios such as isolated operation of local microgrids, emergency power support can be provided to critical loads through internal distributed energy storage [19]. The constraint is expressed as follows:

$${\displaystyle \frac{E_{s,t}^{cap}{\eta }_s}{P_{L_{imp,i}}}}\ge {\beta }_m$$

(32)

where $E_{s,t}^{cap}$ represents the remaining energy of storage s at time t; ${\eta }_s$ denotes the charging and discharging efficiency of the storage, which is generally around 86%; $P_{L_{imp,i}}$ indicates the important load amount at node I; ${\beta }_m$ is the power supply guarantee coefficient for important loads, meaning the proportion of the important load that the remaining energy of the storage at the node should cover at any time, which is generally around 5%.

4. Instance Application Results

4.1. IEEE 14-Bus System

4.1.1. Verification of Optimal Power Flow Accuracy Based on Node Connectivity Matrix

To verify the accurate calculation of system power flow and carbon potential using the algorithm proposed, a comparative analysis was conducted with the IEEE 14-bus system and the literature [16], and present the calculation deviation percentage of key parameters in flowing table. The steady-state active power flow of the IEEE 14-bus system is shown in Appendix A, and the related system parameters are shown in

Table 3. Power flow comparison of the IEEE 14-bus system.

Start Node i	Termination Node j	Branch Real Power Flow/MW
		Results	References	Deviation/%
1	2	77.949	77.949	0.00
2	3	32.982	32.980	0.01
2	4	35.852	35.848	0.01
3	4	−1.217	−1.219	0.23
1	5	42.051	42.051	0.00
2	5	27.414	27.409	0.02
4	5	−36.919	−36.919	0.00
5	6	24.946	24.941	0.02
4	7	12.422	12.420	0.02
7	8	−20.000	−20.000	0.00
4	9	11.332	11.330	0.02
7	9	32.422	32.419	0.01
9	10	5.074	5.078	0.07
6	11	7.426	7.422	0.05
10	11	−3.926	−3.930	0.09
6	12	7.710	7.710	0.00
6	13	17.610	17.610	0.00
12	13	1.610	1.610	0.01
9	14	9.180	9.18	0.00
13	14	5.720	5.720	0.00

.

Based on the above analysis, it can be seen that the error of the active power flow of each branch obtained through the proposed optimal power flow calculation method based on the node connectivity matrix does not exceed 0.23%. This indicates that the proposed method can accurately calculate the system branch power flow, and the results are consistent with the AC power flow calculation results from Matpower, verifying the accuracy of the proposed algorithm.

4.1.2. Verification of Carbon Emission Accuracy Based on Node Linkage Matrix

By combining the active power flux of each node obtained from the optimal power flow calculation with the theoretical carbon emission flow, the carbon potential of each node can be obtained, as shown in

Table 4. Nodal active power flux and nodal carbon intensity of IEEE 14-bus system.

Node	Active Power Flux			Nodal Carbon Intensity
	Results	References	Deviation/%	Results	References	Deviation/%
1	120.000	120.000	0	875.000	875.000	0.00
2	117.949	117.949	0	756.305	756.305	0.00
3	94.200	94.200	0	275.053	275.053	0.00
4	72.771	72.771	0	792.758	792.758	0.00
5	69.465	69.465	0	828.157	828.157	0.00
6	43.946	43.946	0	694.926	694.926	0.00
7	32.422	32.422	0	303.731	303.700	0.01
8	20.000	20.000	0	0	0	0.00
9	43.754	43.754	0	430.386	430.343	0.01
10	9.000	9.000	0	545.794	545.903	0.02
11	7.426	7.426	0	694.926	694.926	0.00
12	7.710	7.710	0	694.926	694.926	0.00
13	19.220	19.220	0	694.926	694.926	0.00
14	14.900	14.900	0	531.938	531.938	0.00

.

Through the proposed carbon potential calculation based on the linkage matrix, the active power flux of each node is consistent with the case study, with a carbon potential error of no more than 0.02%. This not only accurately reflects the active power flowing into each node of the system but also accurately reflects the carbon emission factors corresponding to the electricity consumption of each node, which is conducive to reasonably allocating the carbon emissions from the power supply side to each node.

4.1.3. Verification of Carbon Emission Conservation Based on Node Linkage Matrix

Further, by combining the carbon potential of each node with the power generation and load electricity consumption, the comparison of the load carbon flow rate of each node can be obtained, as shown in

Table 5. Load carbon emission rate of IEEE 14-bus system.

Node	Load Carbon Emission Rate/tCO2/h
	Results	References	Deviation/%
2	16.412	16.410	0.01
3	25.910	25.910	0.00
4	37.894	37.890	0.01
5	6.294	6.290	0.06
6	7.783	7.782	0.04
9	12.696	12.700	−0.03
10	4.912	4.910	0.04
11	2.432	2.430	0.09
12	4.239	4.240	−0.02
13	9.381	9.379	0.02
14	7.926	7.930	−0.05
Total	135.879	135.871	0.01

.

Based on the above analysis, it is known that the total carbon emissions from the load side obtained through the proposed dynamic calculation method for regional carbon emissions based on the linkage matrix are consistent with the total carbon emissions from the source side, indicating the proposed method can effectively measure carbon emissions from the load side. Furthermore, as can be seen from the table above, the traditional method of calculating carbon emissions from the source side based on carbon emission factors cannot account for the carbon emissions corresponding to different electricity consumption. For example, the electricity consumption of Node 3 is significantly higher than other nodes, but the carbon emissions calculated from the source side cannot reflect the actual carbon emissions caused by the real electricity consumption at this node. Therefore, it is necessary to carry out regional carbon emission calculations based on the connectivity matrix to calculate and verify the actual carbon emissions corresponding to the load electricity consumption of each node.

4.2. Actual Provincial System

To verify the rationality of the coordinated output of various types of generating units under each node under the premise of satisfying the safe and stable operation of the power system, a 21-node system composed of various cities in a province was selected as a case study. Its topological structure and power source composition are shown in Figure 2. The bold horizontal lines in the figure represent the busbars of each node, and the numbers next to the busbars represent the node numbers. The numbers in the colored circles correspond to the nodes where each unit is located. The output of wind power and photovoltaic at each node, and the system parameters are shown in Appendix A.

4.2.1. Analysis of Generating Unit Output Combination Characteristics

Through the new power system operation simulation based on the linkage matrix, the output of various types of power generating units in the system can be obtained as shown in Figure 3 and

Table 6. Carbon potential, carbon emissions, and energy composition at each node.

System Boundary	Node	1	2	3	4	5	6	7	8	9	10	11	12	13	14	15	16	17	18	19	20	21	Total
	Total Load/10,000 MWh	3.24	2.42	2.79	1.66	3.77	11.81	5.13	17.43	4.62	2.94	7.70	2.36	2.38	3.57	4.53	2.33	1.38	12.05	4.64	3.56	16.08	116.38
	Photovoltaic Output/10,000 MWh	0.00	0.00	0.00	0.01	0.00	0.00	0.31	0.00	0.09	0.10	0.00	0.00	0.00	0.08	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.60
	Wind Power Output/10,000 MWh	1.13	0.23	1.72	0.06	0.09	0.00	0.16	0.00	0.67	0.42	0.16	0.07	0.04	0.46	0.23	0.14	0.71	0.00	0.00	0.23	0.00	6.51
	Imported Electricity/10,000 MWh	0.00	2.08	0.00	0.00	1.65	0.58	3.58	1.39	0.81	0.00	5.55	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	15.63
	Renewable energy Output Ratio	34.81%	9.46%	61.75%	3.92%	2.29%	0.01%	9.16%	0.00%	16.44%	17.88%	2.11%	3.18%	1.52%	15.17%	4.97%	5.99%	51.40%	0.00%	0.00%	6.47%	0.00%	6.10%
Excluding Renewable energy Sources	Coal Power Output/10,000 MWh	5.56	3.52	3.47	2.06	0.80	3.60	4.82	4.77	0.00	4.37	2.18	2.28	3.41	2.28	4.34	2.37	3.19	3.55	0.79	1.91	3.61	62.87
	Gas Power Output/10,000 MWh	0.00	0.17	0.00	0.00	0.75	1.21	0.63	2.39	0.42	0.00	1.29	0.00	0.00	0.20	0.00	0.13	0.00	3.31	2.15	1.46	2.88	17.01
	Hydropower Output/10,000 MWh	0.00	0.02	0.10	0.00	0.14	0.00	0.01	0.03	0.35	0.32	0.06	0.38	0.36	0.02	0.00	0.04	0.03	0.00	0.00	0.00	0.00	1.86
	Nuclear Power Output/10,000 MWh	0.00	0.00	15.64	0.00	0.00	0.00	8.40	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	14.69	38.73
	Waste-to-Energy Power Output/10,000 MWh	0.39	0.20	0.08	0.00	0.10	0.56	0.03	1.94	0.17	0.78	0.50	0.00	0.04	0.07	0.30	0.14	0.17	0.46	0.42	0.20	1.13	7.68
	Energy Storage Discharge/10,000 MWh	0.02	0.00	0.08	0.00	0.00	0.00	0.00	0.10	0.05	0.00	0.10	0.00	0.08	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.11	0.57
	Energy Storage Charge/10,000 MWh	0.01	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.01
	Wasted Electricity/10,000 WMh	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00
	Carbon Emissions/10,000 tons	3.41	1.80	0.85	1.47	1.95	4.51	1.67	8.45	1.72	2.90	3.04	1.72	2.29	3.89	4.17	2.56	1.72	7.86	2.29	2.80	3.45	64.52
	Carbon Emissions per Unit Load (t/MWh)	1.05	0.74	0.31	0.88	0.52	0.38	0.33	0.49	0.37	0.99	0.40	0.73	0.96	1.09	0.92	1.10	1.25	0.65	0.49	0.79	0.21	0.55
	Inter-node Incoming Electricity Volume (10,000 MWh)	−0.78	−2.29	−15.36	0.32	0.87	7.55	−10.44	9.31	2.96	−1.00	−0.76	0.50	−0.32	1.87	1.41	0.52	−0.89	7.12	2.31	1.15	−4.07	0.00
	Power Supply-Side Calculated Carbon Emissions/10,000 tons	4.94	3.49	3.08	1.83	1.23	3.74	5.04	5.34	0.28	3.88	3.23	2.02	3.03	2.10	3.85	2.15	2.83	4.40	1.52	2.25	4.29	64.52
Including Renewable energy Sources	Coal Power Output/10,000 MWh	5.59	3.51	3.17	2.05	0.81	3.63	4.86	4.84	0.00	4.42	2.19	2.28	3.42	2.22	4.12	2.31	3.16	3.55	0.81	1.92	3.56	62.42
	Gas Power Output/10,000 MWh	0.00	0.18	0.00	0.00	0.51	1.04	0.49	1.75	0.26	0.00	0.88	0.00	0.00	0.20	0.00	0.12	0.00	2.53	1.55	1.06	2.29	12.85
	Hydropower Output/10,000 MWh	0.00	0.02	0.10	0.00	0.14	0.00	0.01	0.03	0.35	0.32	0.06	0.38	0.36	0.02	0.00	0.04	0.03	0.00	0.00	0.00	0.00	1.86
	Nuclear Power Output/10,000 MWh	0.00	0.00	15.64	0.00	0.00	0.00	8.40	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	14.69	38.73
	Waste-to-Energy Power Output/10,000 MWh	0.17	0.10	0.03	0.00	0.05	0.31	0.01	1.13	0.08	0.45	0.28	0.00	0.02	0.03	0.17	0.07	0.09	0.25	0.23	0.11	0.66	4.24
	Energy Storage Discharge/10,000 MWh	0.03	0.00	0.20	0.00	0.00	0.01	0.00	0.30	0.16	0.01	0.22	0.00	0.11	0.01	0.01	0.00	0.01	0.00	0.00	0.00	0.17	1.24
	Energy Storage Charge/10,000 MWh	0.02	0.00	0.14	0.00	0.00	0.01	0.00	0.28	0.15	0.01	0.16	0.00	0.08	0.01	0.01	0.00	0.01	0.00	0.00	0.00	0.14	1.01
	Wasted Electricity/10,000 WMh	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.34	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.00	0.34
Comparison	Carbon Emissions/10,000 tons	2.96	1.75	0.80	1.41	1.91	4.29	1.63	8.45	1.40	2.79	3.20	1.87	2.44	3.49	4.07	2.49	1.35	7.92	2.22	2.70	3.40	62.55
	Carbon Emissions per Unit Load (t/MWh)	0.91	0.72	0.29	0.85	0.51	0.36	0.32	0.48	0.30	0.95	0.42	0.79	1.02	0.98	0.90	1.07	0.99	0.66	0.48	0.76	0.21	0.54
	Inter-node Incoming Electricity Volume (10,000 MWh)	−1.70	−2.41	−16.83	0.27	0.98	7.89	−10.82	10.57	2.44	−1.23	−0.40	0.76	−0.29	1.40	1.46	0.49	−1.50	7.84	2.88	1.28	−3.10	0.00
	Degree of Reduction in Carbon Emissions per Unit	13.27%	2.80%	5.67%	4.03%	2.13%	4.98%	2.49%	0.02%	18.68%	3.73%	−5.24%	−9.00%	−6.63%	10.27%	2.34%	2.89%	21.23%	−0.77%	2.97%	3.68%	1.33%	3.06%
	Power Supply Side Calculated Carbon Emissions/10,000 tons	4.96	3.49	2.81	1.82	1.15	3.70	5.02	5.16	0.21	3.92	3.08	2.03	3.03	2.05	3.66	2.10	2.80	4.11	1.30	2.10	4.03	62.55
	Power Supply Side Degree of Reduction in Carbon Emissions per Unit	−0.53%	0.17%	8.67%	0.70%	6.21%	1.18%	0.48%	3.36%	22.62%	−1.25%	4.45%	−0.13%	−0.25%	2.36%	5.00%	2.54%	0.95%	6.64%	14.20%	6.57%	6.18%	3.06%

.

Through the system operation simulation, it is known that the system load is mainly met by coal power, nuclear power, imported electricity, and renewable energy. During the day from 9:00 to 12:00, due to the significant increase in imported electricity from outside the region, the energy storage begins to charge; to ensure that the amount of electricity at the starting and ending times is equal, the energy storage discharges from 22:00 to 24:00 to meet the load electricity demand. The operating status of other types of power sources also meets the corresponding constraints and is consistent with the actual operating status. To verify the rationality of the established generating unit output combination, the output status of various types of units is represented graphically, as shown in Figure 4.

In Figure 4, the horizontal axis represents the numbers of different types of generating units, and the vertical axis indicates the operating status of the units at each moment (one point every 15 min). The gradual change in colors in the figure represents the output status of the units, with dark blue indicating a shutdown state, and yellow indicating a state of full output. As can be seen from Figure 4a, by analyzing the output status of each unit at each moment, it is known that coal-fired power units have low operating costs, but due to their high carbon emission factors, there are fewer units operating at full capacity for extended periods. The reason why most units operate at a lower status for most of the time is due to line transmission constraints, thus requiring some units to maintain a lower output status. Further attention to the shutdown and continuous operation time of the units reveals that they all meet the start–stop constraints and ramping constraints of the units, with no frequent switching between startup and shutdown states.

From Figure 4b, although the gas-fired power units have lower carbon emission factors, their operating costs are higher than those of coal-fired units. Therefore, most units operate at low load for extended periods, and there is no frequent start–stop phenomenon, which is consistent with actual operating experience. This verifies the rationality of the new power system operation strategy based on the linkage matrix proposed in this paper, which coordinates the operation of source–grid–load–storage.

4.2.2. Impact of Renewable Energy on Regional Carbon Emissions

To analyze the impact of different renewable energy penetration rates on the calculation of node carbon potential and carbon emissions, this paper designed a scenario without renewable energy as a control group, referred to as Case 1, and a scenario with renewable energy as Case 2. Through simulation analysis, the carbon potential and carbon emissions of each node in the system before and after the addition of renewable energy are shown in Figure 5.

From the carbon potential and emissions of each node, it is evident that the system’s carbon emissions trend is consistent with the load electricity consumption. With the introduction of renewable energy, the carbon potential and emissions of each node have decreased. It is particularly noteworthy that the carbon emissions of nodes 6, 8, 18, and 21 at each moment are significantly higher than those of other nodes. Combined with Figure A4d, it is known that the load levels of nodes 6, 8, 18, and 21 are significantly higher than those of other nodes, resulting in significantly higher carbon emissions attributed to the load side. To verify the correlation between the carbon potential of each node and the composition of the power generation sources of the node, a detailed analysis of the output of various types of power sources under each node at 14:00 is conducted, as shown in Figure 6.

Figure 7 shows that nodes 1 and 17 have significantly reduced carbon potential after the introduction of renewable energy due to their high proportion of renewable energy output. However, since the load of these two nodes is relatively low, the reduction in carbon emissions is not significant. In addition, nodes 3, 7, and 21 are connected to nuclear power units. During production operation, the output proportion of nuclear power is high, and the carbon emission factor of nuclear power is much lower than that of coal-fired and gas-fired power. Therefore, the carbon emissions of these three nodes are always lower than those of other nodes with equivalent loads. Further sorting out the comparison of the hourly carbon emission curves of each node under the two scenarios is shown in Figure 7.

Comparing the carbon emissions of each node at each moment, it can be seen that node 3, with a renewable energy proportion of 61.75%, has significantly lower carbon emissions than other nodes. Node 7 sees a significant increase in photovoltaic output in the morning, resulting in a clear downward trend in its carbon emissions. A detailed breakdown of the carbon emissions and energy composition of each node is further shown in Table 6.

Through the comparative analysis of carbon emissions in the provincial system under scenarios with and without renewable energy, it is observed that the output of gas-fired power significantly decreases under the scenario with renewable energy. However, since the proportion of gas-fired power output is small and the carbon emission factor of gas-fired power units is 0.3789, the degree of reduction in carbon emissions per unit is only 3.06% when the proportion of renewable energy is 6.10%.

4.2.3. Analysis of the Rationality of Carbon Emissions Calculated on the Load Side

To analyze the correlation between the proportion of renewable energy output at each node and the degree of reduction in carbon potential, the relationship between the proportion of renewable energy output at each node and the degree of reduction in carbon emissions is organized, as shown in Figure 8.

From Figure 8, it can be seen that the degree of carbon emissions reduction in various nodes is generally positively correlated with the proportion of renewable energy. However, due to the existence of branch links between cities, carbon emissions are influenced by adjacent cities.

Taking Node 3 as an example, although it has a high proportion of renewable energy, the degree of carbon emission reduction calculated from the source side is relatively low. This is because some of the renewable energy is provided to the adjacent Node 4, which significantly reduces the degree of carbon emission reduction at Node 4.

A further comparison of the carbon emissions of system nodes calculated from the source side and the load side is shown in Figure 9.

Figure 9 shows that the traditional method of calculating carbon emissions from the power supply side is related to the power source structure, while the calculation of carbon emissions from the load side is related to the node carbon potential during the load electricity consumption phase, which can more accurately reflect the actual carbon emissions of the region.

5. Conclusions

This paper establishes a coordinated optimization operation simulation model for source–grid–load–storage systems with the objective of minimizing system costs. It introduces a power flow calculation method and a dynamic zonal carbon emission calculation method, both based on a linkage matrix and the theory of carbon emission flows. The power flow calculation method replaces traditional nonlinear constraints with linear ones, thereby avoiding the computational burden of invoking Matpower and not only reducing the complexity of model solving but also ensuring the accuracy of model solutions. The dynamic zonal carbon emission calculation method addresses the shortcomings of existing methods, which are too coarse-grained and unable to accurately determine the operational states of multiple units at a node. This method enhances the precision of carbon emission metering and clarifies the improvement effect of renewable energy output on the carbon emissions of adjacent regions. The effectiveness and practicality of the proposed methods are validated through case studies based on the IEEE 14-bus system and actual operational data from a provincial power grid. The results of this paper not only provide theoretical support for the coordinated optimization of power system operation and carbon emission management but also offer significant technical references and practical guidance for the sustainable development of future power systems. Although the methods presented in this paper demonstrate good feasibility and guiding value in practical applications, further validation and optimization are needed under broader system scales and more complex operational conditions. Future research could consider incorporating more practical factors, such as market mechanisms and policy elements, into the model to enhance its adaptability and predictive accuracy.