A method for unit commitment based on cuppman partition considering regional inertia constraint
By identifying unit coherence relationships through Koopman decomposition and k-means clustering, regional inertia constraints are constructed, solving the problems of low accuracy and high computational cost in existing inertia constraint methods, and achieving a balance between frequency security coverage and economy.
Patent Information
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- SOUTHEAST UNIV
- Filing Date
- 2026-04-13
- Publication Date
- 2026-07-21
AI Technical Summary
Existing inertia-constrained unit combination methods cannot effectively reflect changes in the coherence relationship between units, resulting in uncovered frequency security risks or economic losses. Furthermore, they are computationally burdensome and difficult to efficiently couple with large-scale mixed-integer optimization.
By establishing a transient frequency response database, the angular velocity and phase angle deviation of the unit are extracted. Koopman decomposition and k-means clustering are used to identify the coherence relationship, construct regional inertia constraints, and embed them into the unit combination model to meet frequency safety requirements.
It improves frequency security coverage, reduces computational burden, enhances the engineering feasibility of unit combinations, and balances safety and economy.
Smart Images
Figure CN122437176A_ABST
Abstract
Description
Technical Field
[0001] This invention belongs to the technical field of power system operation optimization and frequency security constraint construction. Specifically, it relates to a unit combination method based on Koopman partitioning that takes into account regional inertia constraints. This method can be used to guarantee transient frequency security margin in the form of calculable algebraic constraints during the scheduling phase. Background Technology
[0002] Frequency security of power systems is one of the core issues for the safe and stable operation of the system. With the increasing penetration rate of power electronic power sources such as wind power and photovoltaics, the equivalent inertia of the system decreases and its spatial distribution is uneven. The traditional "system-level total inertia constraint" that relies on the synchronous machine aggregation assumption is difficult to characterize the real regional frequency dynamic differences: Under the same total inertia level, the geographical location and electrical coupling structure of different units may lead to more severe frequency change rate (RoCoF) and frequency minimum point (nadir) risks in local areas, which may trigger protection actions or cascading load shedding.
[0003] Existing inertia-constrained unit combination methods often derive lower bounds of inertia using system-level single-unit equivalent or network-wide aggregation models, and then use these as additional inequality constraints for unit combination. However, these methods typically implicitly assume "network-wide frequency consistency" or "uniform diffusion of disturbances," failing to reflect the changing characteristics of unit coherence relationships under different operating modes with variations in disturbances and operating conditions. This easily leads to two types of problems: first, overly conservative constraints result in economic losses; second, insufficient constraints leave local frequency risks uncovered. Furthermore, online safety verification based on detailed dynamic models often requires frequent calls to transient simulations or sensitivity calculations during iteration, resulting in heavy computational burdens and difficulty in efficiently coupling with large-scale mixed-integer optimization. Therefore, there is an urgent need for a method that can automatically extract "frequency coherence regions" from offline transient data and map regional frequency risks to regional inertia constraints that can be embedded in unit combination, thus balancing safety, computability, and economy. Summary of the Invention
[0004] This invention addresses the problems existing in the prior art by providing a unit combination method based on Koopman partitioning that considers regional inertia constraints. First, a transient frequency response database is established for different fault scenarios, and the angular velocity and phase angle deviation of each unit are extracted to form a frequency database. Then, Koopman spectrum analysis is performed on the system based on extended dynamic mode decomposition to obtain Koopman tuples containing Koopman eigenvalues, eigenfunctions, and Koopman modes, and the phase features of the dominant mode are further extracted. Based on the unit phase signature vectors under each fault scenario, k-means clustering is used to identify the coherence relationship between units, forming a candidate coherence partitioning library, and representative coherence region sets are obtained through coverage filtering. On this basis, offline fault samples are used to statistically analyze the power imbalance samples and their 95%... The perturbation magnitude is used to construct regional inertia constraints for the rate of frequency change and the minimum frequency point. Finally, these regional inertia constraints are embedded into a classic thermal power unit combination model, and together with operational constraints such as power balance, unit output limits, spinning reserve, ramp-up constraints, and minimum start-up / shutdown time, form a unit combination optimization model with regional inertia constraints. Solving this model yields the optimal unit combination scheme that meets regional frequency security requirements. This invention uses a data-driven Koopman operator method, which overcomes model inaccuracies while considering the nonlinear characteristics of the power system dynamic model, effectively solving the low accuracy problem of traditional linear methods and the excessive computational cost of nonlinear methods.
[0005] To achieve the above objectives, the technical solution adopted by the present invention is: a unit combination method based on Koopman partitioning that takes into account regional inertia constraints, comprising the following steps:
[0006] S1: Establish a transient frequency response database, and construct a frequency database based on the angular velocity and azimuth deviation of each unit under each fault condition;
[0007] S2: Establish an augmented Koupman operator, and extract spectral features of the unit's angular velocity and azimuth deviation under various fault conditions based on the transient frequency response database to obtain the Koupman tuples after decomposition of each variable; the Koupman tuples include Koupman eigenvalues, Koupman characteristic functions, and Koupman modes;
[0008] S3: Based on the Koopman tuples obtained in step S2, extract modal phase features, and form homology groups by k-means clustering based on the modal phase features. Further construct candidate homology partitions based on one or more of the homology groups. The homology group is a set of units with similar dynamic response characteristics under the same disturbance scenario or dominant mode. The candidate homology partition is a system region division scheme composed of one or more of the homology groups, used to characterize the regional homology structure that may occur in the power system under different fault scenarios.
[0009] S4: For any candidate homology partition, calculate its coverage and set a coverage threshold; determine the candidate homology partitions with coverage greater than or equal to the coverage threshold as the representative homology partition set;
[0010] S5: Based on partitioning, construct constraints, extract the sample set from the transient frequency response database and the representative homology region obtained in step S4, calculate the 95th quantile after taking the absolute value, obtain the regional disturbance magnitude, and construct constraints based on partitioning, including the frequency change rate constraint based on the regional disturbance magnitude and the maximum frequency deviation constraint based on the regional disturbance magnitude.
[0011] S6: Embed the frequency change rate constraint and maximum frequency deviation constraint constructed in step S5 into the classic thermal power unit combination model, and combine the minimization objective function and operation constraints of the unit combination to construct a unit combination optimization model with regional inertia constraints, and solve to obtain the optimal unit combination scheme that meets the regional inertia safety requirements.
[0012] As an improvement of the present invention, the transient frequency response database in step S1 includes at least generator active power, generator reactive power, load active power, load reactive power, voltage amplitude, and phase angle data, and extracts the angular velocity of each unit under each fault condition p. and angular deviation Build a frequency database :
[0013] .
[0014] As an improvement to the present invention, the method for obtaining the Koopman tuple in step S2 is specifically as follows:
[0015] From frequency database In, for any fault scenario Extracting data from each unit at adjacent sampling times and angular velocity at the bottom and phase angle deviation Construct state observation vectors respectively
[0016]
[0017] and
[0018]
[0019] in, This indicates the number of units; further, multiple current-moment state observation vectors are arranged column-wise to form a data matrix.
[0020]
[0021] The corresponding next-time state observation vectors are arranged into a data matrix by column.
[0022]
[0023] in, Indicates the number of sampling points;
[0024] Define vector-valued observation subfunction This is used to map the state observation vector to the lifting observation space; where
[0025]
[0026] Furthermore, for the data matrix respectively and Apply the vector-valued observation sub-function to each column of the state observation vector to construct the observable matrix.
[0027]
[0028] and
[0029]
[0030] Constructing the finite-dimensional approximation matrix of the Koopman operator
[0031] Its expression is
[0032]
[0033] in It is a Moore-Penrose pseudoreverse;
[0034] For the obtained finite-dimensional approximation matrix Perform eigenvalue decomposition, calculate the left eigenvectors, and store each left eigenvector in a matrix column-wise.
[0035]
[0036] in, Representation matrix Corresponding to the The left eigenvector of eigenvalues, This indicates the number of extracted feature pairs; further, based on the left feature vector matrix and the vector-valued observation sub-function, the Koopman feature function of the system is estimated.
[0037]
[0038] Based on Koupman characteristic function The evolution of the system state observations can be expressed in a finite-dimensional Koupman expansion:
[0039]
[0040] in, This represents the i-th Koopman eigenvalue. Let i represent the i-th Koopman characteristic function. This represents the i-th Koopman pattern;
[0041] Furthermore, based on the obtained Koopman eigenvalues, Koopman characteristic functions, and Koopman modes, the Koopman tuples of the system are constructed. The subsequent main use will be of the Koopman pattern in the Koopman tuple. Dominant mode screening is performed, and phase feature extraction is performed in conjunction with eigenvalue information.
[0042] As another improvement of the present invention, in step S3...
[0043] For any fault scenario Given the number of partitions ,by As the input for clustering, the K-means method is used to solve the following objective function to obtain the cluster centers. and the corresponding homology region
[0044]
[0045] Each failure scenario The clustering results below are recorded as a partition result:
[0046]
[0047] in, Indicates the fault scenario The next Each region is a homogeneous region. Furthermore, the single-partition results corresponding to all failure scenarios are aggregated to form a partition library:
[0048] .
[0049] As another improvement of the present invention, in step S4, for any candidate homology partition, its coverage rate is calculated and a coverage rate threshold is set; the candidate homology regions with coverage rates greater than or equal to the coverage rate threshold are determined as a representative homology region set.
[0050] For partitioned libraries The different cohomological regions appearing in the data are deduplicated and summarized to form a candidate cohomological region set. ;
[0051] For partitioned libraries For any candidate homology region appearing in the data, calculate its coverage:
[0052]
[0053] in, For partitioned libraries, For any candidate homology region, Indicates partitioned library The total number of candidate homology partitions, Indicates partitioned library The total number of results in the middle partition. This is an indicator function: it takes a value of 1 when the event is true, and a value of 0 otherwise;
[0054] Set coverage threshold , retain satisfaction The homology regions form a representative set of homology regions. The Used for the construction and embedding of subsequent regional inertia constraints.
[0055] Compared with existing technologies, the technical advantages and effects of this invention are as follows: This invention discloses a unit combination method based on Koopman partitioning that takes into account regional inertia constraints. It extracts the phase characteristics of the dominant dynamic modes of the power system in Koopman standard coordinates using an offline transient data-driven approach, thereby identifying the dynamic coherence relationships between units and constructing a representative set of coherent regions. Furthermore, it combines the power imbalance statistics and frequency safety indicators of each region under fault samples to form a lower bound constraint of regional inertia that can be directly embedded into the unit combination optimization model. This invention's scheme, through transient simulation, can reflect the spatial heterogeneity of inertia and the changes in coherence structure caused by operating modes, thereby improving frequency safety coverage. Simultaneously, the maximum frequency deviation constraint and the frequency change rate constraint constructed from the regional disturbance magnitude are expressed in a compact algebraic form, avoiding repeated transient simulations or complex sensitivity calculations during unit combination solution, significantly reducing the computational burden and improving engineering feasibility. This invention effectively solves the low accuracy problem of traditional linear methods and the excessive computational burden of nonlinear methods. Attached Figure Description
[0056] Figure 1 This is a flowchart of the steps of a unit combination method based on Koopman partitioning that takes into account regional inertia constraints according to the present invention. Detailed Implementation
[0057] The present invention will be further illustrated below with reference to the accompanying drawings and specific embodiments. It should be understood that the following specific embodiments are for illustrative purposes only and are not intended to limit the scope of the invention.
[0058] Example 1
[0059] A unit combination method based on Koopman partitioning that takes into account regional inertia constraints is proposed. The power system in this embodiment is a New England 39 bus system using a classic generator model, wherein a set of PMUs are placed at the generator terminals to provide... , , and The measured value. For example... Figure 1 As shown, this method specifically includes the following steps:
[0060] Step S1: Establish a transient frequency response database and fit the current operating scenario to obtain the frequency database.
[0061] An initial database containing a large amount of historical simulation results is established, including data such as generator active power, generator reactive power, load active power, load reactive power, voltage amplitude, and phase angle. The angular velocity of each unit under each fault condition p is extracted. and angular deviation Construct a frequency database, denoted as:
[0062]
[0063] In this embodiment, short-circuit faults and tangential disturbances are applied to 10 relative positions within each transmission line, and the position parameters are taken as follows: The simulation duration was 20 seconds with intervals of 0.01 seconds; the angular velocity of each unit under each fault condition p was extracted. and angular deviation Construct a frequency database, denoted as:
[0064] .
[0065] Step S2: Establish the augmented Koopman operator, and extract spectral features of the unit's angular velocity and azimuth deviation under various fault conditions based on the frequency database to obtain the Koopman tuples after decomposition of each variable; the Koopman tuples include Koopman eigenvalues, Koopman characteristic functions, and Koopman modes.
[0066] S21: For frequency databases For any fault scenario p, extract the angular velocity of each unit at discrete time t. Phase angle deviation Construct state observation vector And construct a set of observable functions based on the state observation vector. ;
[0067] S22: Based on the set of observable functions A Koopman operator is established to characterize the dynamic mapping relationship between the unit's angular velocity and phase angle deviation over time under various fault scenarios in the frequency database;
[0068] S23: Construct a mapping function based on the selected angular velocity and phase angle observation set:
[0069]
[0070] in, This represents the state observation vector at time t. This represents the mapping function from the state observation vector to the observation space. Indicates the dimension of the observation space. This represents the dimension of the state observation vector; it is used for subsequent estimation of the Koopman eigenvalues and Koopman characteristic functions of the system.
[0071] S24: Based on the state observation vector and mapping function constructed in steps S21-S23, the extended dynamic mode decomposition method is used to approximate the Koopman operator in finite dimensions, and the Koopman eigenvalues, Koopman characteristic functions and Koopman modes of the system are obtained, thereby obtaining the Koopman tuples.
[0072] The extended dynamic mode decomposition method in step S24 specifically includes the following steps:
[0073] S241, from frequency database In, for any fault scenario Extracting data from each unit at adjacent sampling times and angular velocity at the bottom and phase angle deviation Construct state observation vectors respectively
[0074]
[0075] and
[0076]
[0077] in, Indicates the number of generating units;
[0078] Furthermore, multiple current state observation vectors are arranged into a data matrix column by column.
[0079]
[0080] The corresponding next-time state observation vectors are arranged into a data matrix by column.
[0081]
[0082] in, Indicates the number of sampling points;
[0083] S242. Based on the state observation vector constructed in step S241, define the vector value observation sub-function. This is used to map the state observation vector to the lifting observation space; where
[0084]
[0085] Furthermore, for the data matrix respectively and Apply the vector-valued observation sub-function to each column of the state observation vector to construct the observable matrix.
[0086]
[0087] and
[0088]
[0089] S243, The observable matrix constructed based on step S242 and Establish the finite-dimensional approximation matrix of the Koopman operator.
[0090] Its expression is
[0091]
[0092] in It is a Moore-Penrose pseudoreverse;
[0093] S244. The finite-dimensional approximation matrix obtained in step S243 Perform eigenvalue decomposition, calculate the left eigenvectors, and store each left eigenvector in a matrix column-wise.
[0094]
[0095] in, Representation matrix Corresponding to the The left eigenvector of eigenvalues, This indicates the number of extracted feature pairs; further, based on the left feature vector matrix and the vector-valued observation sub-function, the Koopman feature function of the system is estimated.
[0096]
[0097] S245. Based on the Koopman characteristic function obtained in step S244 The evolution of the system state observations is expressed as a finite-dimensional Koupman expansion.
[0098] in, This represents the i-th Koopman eigenvalue. Let i represent the i-th Koopman characteristic function. This represents the i-th Koopman pattern;
[0099] S246. Based on the obtained Koopman eigenvalues, Koopman characteristic functions, and Koopman patterns, construct the Koopman tuples of the system. The subsequent main use is of the Koopman pattern in the Koopman tuple. Dominant mode screening is performed, and phase feature extraction is performed in conjunction with eigenvalue information.
[0100] Step S3: Extract modal phase features based on the obtained Koopman tuples, and form homology groups based on the modal phase features through k-means clustering. Further construct candidate homology partitions based on one or more of the homology groups. The homology group is a set of units with similar dynamic response characteristics under the same disturbance scenario or dominant mode. The candidate homology partition is a system region division scheme composed of one or more of the homology groups, used to characterize the regional homology structure that may occur in the power system under different fault scenarios.
[0101] S31, Based on the obtained The Koopman tuple will be the first... Each Koopman pattern is represented as a vector:
[0102]
[0103] in To measure the number of units, Indicates the first Measurement locations in the mode Modal components below;
[0104] S32. Calculate the Euclidean norm of each modal vector as a contribution index.
[0105]
[0106] in Used for quantizing modes The overall contribution to the measured dynamic data;
[0107] S33, sort all candidate modes by Sort the modes from largest to smallest, and select the mode with the largest norm to form the target mode set, denoted as:
[0108]
[0109] S34. Perform amplitude-phase decomposition on the selected dominant mode: for any and any The modal components are written in amplitude-phase form:
[0110]
[0111] in For amplitude, This is the initial phase, used to characterize the phase consistency of different units in this mode;
[0112] S35. The above-mentioned Koopman mode extraction, dominant mode selection, and amplitude-phase decomposition processes are all performed separately for each fault scenario. Therefore, when constructing the unit phase characteristics, a fault scenario index is introduced. To characterize the phase response features of the unit under different disturbance scenarios;
[0113] Phase characteristics of each unit based on the dominant mode Construct the first Phase characteristic vector of each unit:
[0114]
[0115] in, Indicating a fault scenario Lower unit In the selected dominant mode set Phase signature vector on; Indicating a fault scenario Lower unit In the Initial phase on each dominant mode;
[0116] S36. For any fault scenario Given the number of partitions ,by As the input for clustering, the K-means method is used to solve the following objective function to obtain the cluster centers. and the corresponding homology region
[0117]
[0118] in, Indicating a fault scenario Next The set of samples corresponding to each cluster Indicating a fault scenario Lower unit The phase eigenvector, Indicates the fault scenario Next The center vector of each cluster; Represents the L2 norm, Indicates sample and its cluster center The square of the Euclidean distance between them;
[0119] S37. Each fault scenario... The clustering results below are recorded as a partition result:
[0120]
[0121] in, Indicates the fault scenario The next One homology region;
[0122] S38. Summarize the single-partition results corresponding to all failure scenarios to form a partition library:
[0123] .
[0124] Step S4: For any candidate homology partition, calculate its coverage rate and set a coverage rate threshold; determine the candidate homology regions with coverage rates greater than or equal to the coverage rate threshold as the representative homology region set.
[0125] S41, Partition database The different cohomological regions appearing in the data are deduplicated and summarized to form a candidate cohomological region set. ;
[0126] S42, Partitioning the database For any candidate homology region appearing in the data, calculate its coverage:
[0127]
[0128] in, For partitioned libraries, For any candidate homology region, Indicates partitioned library The total number of candidate homology partitions, Indicates partitioned library The total number of results in the middle partition. This is an indicator function: it takes the value 1 when the event is true, and 0 otherwise.
[0129] S43. Set coverage threshold , retain satisfaction The homology regions form a representative set of homology regions. The Used for the construction and embedding of subsequent regional inertia constraints.
[0130] Step S5: Extract the regional power imbalance value and construct constraints based on the partitioning situation:
[0131] S51, Based on representative homology region sets Obtain the magnitude of regional disturbance Specifically, it includes:
[0132] S511, For any representative homology region Construct regional power imbalance samples based on fault samples:
[0133]
[0134] in Indicates the first In a fault scenario, the synchronization region Regional power loss measurement; Indicates the first Unit under various failure scenarios Mechanical input power; Indicates the first Unit under various failure scenarios Electromagnetic output power.
[0135] S512. Calculate the 95th quantile after taking the absolute value of the sample set to obtain the magnitude of the regional disturbance:
[0136]
[0137] in This represents the 95th percentile of the sample set;
[0138] S52. Constructing a rate of change of frequency (RoCoF) constraint based on the magnitude of regional perturbations:
[0139]
[0140] in Indicates the unit During the period The power-on status, The system's rated frequency, This is the upper limit of the allowable rate of change of frequency. For the region The power base value, Obtained from offline perturbation statistics Percentage power imbalance value; Indicates the unit During the period The startup status variable, when the unit During the period The value is 1 when the device is powered on, otherwise it is 0. Indicates the unit inertia constant; A set of representative homology regions;
[0141] S53. Constructing the maximum frequency deviation (Nadir) constraint based on the magnitude of regional perturbation:
[0142]
[0143] in This is the load damping coefficient. Obtained from offline perturbation statistics quantile power imbalance value This is the equivalent droop parameter for primary frequency modulation. The equivalent coefficient parameters for a primary frequency modulation channel are: The equivalent time constant parameter for primary frequency modulation; and These represent the damping ratio and natural angular frequency of the equivalent second-order frequency response model, respectively. The time when the frequency is at its lowest point. It is the exponential factor of the second-order response decay term. Indicates the unit During the period The startup status variable, when the unit During the period The value is 1 when the device is powered on, otherwise it is 0. Indicates the unit The inertia constant; It is a set of representative homology regions.
[0144] Step S6: Embed the regional inertia constraint constructed in Step S5 into the classic thermal power unit combination model, and combine it with the minimization objective function and operating constraints of the unit combination to construct a unit combination optimization model with regional inertia constraints, and solve to obtain the optimal unit combination scheme that meets the regional inertia safety requirements.
[0145] S61. Construct the minimization objective function for the thermal power unit combination problem:
[0146]
[0147] constant The coefficients of the cost function for the k-th unit are represented. This indicates the on / off state of the k-th unit during time period t. The value is 1 if the unit is online, and 0 otherwise. This represents the power generation of the k-th generating unit during time period t;
[0148] S62. The constraints for the thermal power unit combination problem are mainly: power generation constraints, power balance constraints, spinning reserve constraints, minimum start-up / shutdown time constraints, start-up cost constraints, ramp rate constraints, and the previously constructed minimum inertia constraints.
[0149] The power generation constraint formula is:
[0150]
[0151] in Indicates the unit Minimum effort, Indicates the unit Maximum output, Indicates the unit exist The output status at any given moment;
[0152] The power balance constraint formula is:
[0153]
[0154] in Indicates time period The total system load requirement;
[0155] The formula for the rotational reserve constraint is:
[0156]
[0157] This represents the rotating reserve requirement for time period t, where Indicates time period The total system load requirement;
[0158] The minimum power-on / power-off time constraint formula is:
[0159]
[0160] Indicates the unit Minimum continuous power-on duration; Indicates the unit During the period The startup status variable, if the unit is in the time period If enabled, set the value to 1; otherwise, set the value to 0. Indicates the unit The minimum continuous downtime; Indicates the unit During the period The shutdown status variable, if the unit is in the time period If the machine is stopped, use 1; otherwise, use 0.
[0161] The start-up cost constraint formula is:
[0162]
[0163] Indicates the unit Start-up and shutdown costs during time period t; Indicates the unit Startup costs; Indicates the unit Downtime costs; Indicates the unit The start indicator variable for time period t is set to 1 when start is initiated, and 0 otherwise. Indicates the unit The shutdown indicator variable in time period t is set to 1 when shutdown occurs, and 0 otherwise.
[0164] The formula for the climbing rate constraint is:
[0165]
[0166] Indicates the unit The upper limit of the upward climbing rate, that is, the maximum output that can be increased per unit time; Indicates the unit The upper limit of the downward ramp rate, that is, the maximum output that can be reduced per unit time;
[0167] S63. Add the regional inertia constraints constructed in step S5 to obtain the unit combination optimization model containing partitioned inertia constraints:
[0168]
[0169]
[0170] in Indicates the unit During the period The power-on status, The system's rated frequency, This is the upper limit of the allowable rate of change of frequency. For the region The power base value, Obtained from offline perturbation statistics Percentage power imbalance value; This is the load damping coefficient. This is the equivalent droop parameter for primary frequency modulation. The equivalent coefficient parameters for a primary frequency modulation channel are: The equivalent time constant parameter for primary frequency modulation; and These represent the damping ratio and natural angular frequency of the equivalent second-order frequency response model, respectively. The time when the frequency is at its lowest point. It is the exponential factor of the second-order response decay term. Indicates the unit During the period The startup status variable, when the unit During the period The value is 1 when the device is powered on, otherwise it is 0. Indicates the unit The inertia constant; It is a set of representative homology regions.
[0171] In the application of this invention, the Gurobi solver was called on the MATLAB 2022b platform to optimize and solve the classical thermal power unit combination model and the unit combination optimization model with partitioned inertia constraints, respectively. In the simulation, taking the IEEE-39 Node system as an example, the 24-hour output and start-up / shutdown results of each generator unit in the unit combination optimization model with partitioned inertia constraints were obtained. Then, the proposed unit combination optimization model with partitioned inertia constraints was compared with the classical thermal power unit combination model. Taking the IEEE-39 Node system as an example, the day-ahead scheduling simulation was performed and compared using the Gurobi solver on the MATLAB 2022b platform. The data is shown in Table 1 below:
[0172] Table 1. Statistical Comparison of 24-Hour Scheduling Economy and Commitment for Different Inertia Constraint Strategies
[0173]
[0174] As shown in Table 1, compared with the classic thermal power unit combination model, the total cost of the unit combination optimization model with zonal inertia constraints increased from $48,849.2 to $50,370.5. Simultaneously, the average number of online units increased from 7.17 to 9.33, indicating that the system improves its overall safety factor by increasing online synchronous support resources to meet regional inertia safety requirements. The results in Table 1 demonstrate that the zonal inertia constraints proposed in this invention can effectively influence unit combination decisions, ensuring that scheduling results no longer solely aim at minimizing operating costs, but also consider regional frequency safety requirements while maintaining economic efficiency. Although the system operating cost increases, the increased number of online units provides the system with more sufficient inertia support and frequency regulation capabilities, thereby enhancing the system's frequency safety margin in response to local disturbances. This indicates that the method of this invention can effectively transmit the regional inertia requirements identified offline to the scheduling layer in the form of optimization constraints, demonstrating good safety improvement effects and engineering application value.
[0175] Example 2
[0176] The power system in this embodiment is an IEEE 118 bus system using a classic generator model. A set of PMUs is placed at the generator terminal to provide power. , , and Based on the measured values and simulation parameters of Example 1, day-ahead scheduling simulations were performed on both the classic thermal power unit combination model and the unit combination optimization model with partition inertia constraints. The results are compared in Table 2.
[0177] Table 2. Statistical Comparison of 24-Hour Scheduling Economy and Commitment for Different Inertia Constraint Strategies
[0178]
[0179] As shown in Table 2, compared with the classic thermal power unit combination model, the total cost of the unit combination optimization model with zonal inertia constraints increased from $266,187 to $275,596. Meanwhile, the average number of online units increased from 18.08 to 39.79, indicating that the system meets regional inertia safety requirements by increasing online synchronous support resources, thus improving the overall safety factor. This demonstrates that in larger-scale systems, the method of this invention can also effectively influence the unit combination decision-making results, enabling the scheduling scheme to further meet regional frequency safety requirements while considering operational economy. Although the system operating cost has increased, keeping more synchronous units online provides the system with more sufficient inertia support and frequency regulation capabilities, thereby enhancing the system's frequency safety margin in response to local disturbances. This shows that the method of this invention still has good effectiveness and engineering applicability in large-scale power systems.
[0180] It should be noted that the above content merely illustrates the technical concept of the present invention and should not be construed as limiting the scope of protection of the present invention. For those skilled in the art, various improvements and modifications can be made without departing from the principle of the present invention, and all such improvements and modifications fall within the scope of protection of the claims of the present invention.
Claims
1. A unit combination method based on Koopman partitioning that takes into account regional inertia constraints, characterized in that, The method includes the following steps: It is characterized by including the following steps: S1: Establish a transient frequency response database, and construct a frequency database based on the angular velocity and azimuth deviation of each unit under each fault condition; S2: Establish an augmented Koopman operator. Based on the transient frequency response database established in step S1, extract spectral features from the angular velocity and azimuth deviation of the unit under various fault conditions to obtain the Koopman tuples after decomposition of each variable. The Koopman tuples include Koopman eigenvalues, Koopman characteristic functions, and Koopman modes. S3: Extract modal phase features based on the Koopman tuples obtained in step S2, and form homology groups based on the modal phase features through k-means clustering. Construct candidate homology partitions based on one or more homology groups. The homology group is a set of units with similar dynamic response characteristics under the same disturbance scenario or dominant mode. The candidate homology partition is a system region division scheme composed of one or more of the homology groups. S4: For any candidate cohomology region obtained in step S3, calculate its coverage rate and set a coverage rate threshold; determine the candidate cohomology regions with coverage rates greater than or equal to the coverage rate threshold as the representative cohomology region set; S5: Based on partitioning, construct constraints, extract the sample set from the transient frequency response database and the representative homology region obtained in step S4, calculate the 95th quantile after taking the absolute value, obtain the regional disturbance magnitude, and construct constraints based on partitioning, including the frequency change rate constraint based on the regional disturbance magnitude and the maximum frequency deviation constraint based on the regional disturbance magnitude. S6: Embed the frequency change rate constraint and maximum frequency deviation constraint constructed in step S5 into the classic thermal power unit combination model. Based on the minimization objective function and operating constraints of the unit combination, construct a unit combination optimization model with partition inertia constraints, and solve it to obtain the optimal unit combination scheme.
2. The unit combination method based on Koopman partitioning considering regional inertia constraints as described in claim 1, characterized in that: The transient frequency response database in step S1 includes at least generator active power, generator reactive power, load active power, load reactive power, voltage amplitude, and phase angle data. The angular velocity of each unit under each fault condition p is extracted from this database. and angular deviation Build a frequency database : 。 3. The unit combination method based on Koopman partitioning considering regional inertia constraints as described in claim 2, characterized in that: Step S2 specifically includes the following steps: S21: For frequency databases Any of the following fault scenarios Extract the angular velocity of each unit at discrete time t. Phase angle deviation Construct state observation vector And construct a set of observable functions based on the state observation vector. ; S22: Based on the aforementioned set of observable functions A Koopman operator is established to characterize the dynamic mapping relationship between the unit's angular velocity and phase angle deviation over time under various fault scenarios in the frequency database; S23: Constructing a mapping function based on angular velocity and phase angle observation sets: ； in, This represents the state observation vector at time t. This represents the mapping function from the state observation vector to the observation space. Indicates the dimension of the observation space. This represents the dimension of the state observation vector; S24: Based on the state observation vector and mapping function constructed in steps S21-S23, the extended dynamic mode decomposition method is used to perform a finite-dimensional approximation of the Koopman operator, and the Koopman eigenvalues, Koopman characteristic functions and Koopman modes of the system are obtained to obtain the Koopman tuples.
4. The unit combination method based on Koopman partitioning considering regional inertia constraints as described in claim 3, characterized in that: The extended dynamic mode decomposition method in step S24 specifically includes the following steps: S241, from frequency database In, for any fault scenario Extracting data from each unit at adjacent sampling times and angular velocity at the bottom and phase angle deviation Construct state observation vectors respectively: ； and ； in, Indicates the number of generating units; Arrange multiple current state observation vectors into a data matrix by column: ； The corresponding next-time state observation vectors are arranged into a data matrix by column. ； in, Indicates the number of sampling points; S242. Based on the state observation vector constructed in step S241, define the vector value observation sub-function. : ； For the data matrix respectively and Apply the vector-valued observation sub-function to each column of the state observation vector to construct the observable matrix: ； and ； S243, The observable matrix constructed based on step S242 and Establish the finite-dimensional approximation matrix of the Koopman operator: ； Its expression is: ； in It is a Moore-Penrose pseudoreverse; S244. The finite-dimensional approximation matrix obtained in step S243 Perform eigenvalue decomposition, calculate the left eigenvectors, and store each left eigenvector in a matrix column-wise. middle: ； in, Representation matrix Corresponding to the The left eigenvector of eigenvalues, This indicates the number of feature pairs extracted; Based on the left eigenvector matrix and the vector-valued observation subfunction, estimate the Koopman eigenfunction of the system: ； S245. Based on the Koopman characteristic function obtained in step S244 The evolution of the system state observations can be expressed in a finite-dimensional Koupman expansion: ； in, This represents the i-th Koopman eigenvalue. Let i represent the i-th Koupman characteristic function. Let i represent the i-th Koopman pattern; the Koopman tuple of the system is .
5. The unit combination method based on Koopman partitioning considering regional inertia constraints as described in claim 1, characterized in that: In step S3, the method for extracting modal phase features based on Koopman tuples specifically involves: extracting the Koopman modes from the Koopman tuples. Representing modes as vectors, calculate the Euclidean norm of each mode vector, sort all candidate modes in descending order of their Euclidean norms, and select the mode with the largest norm to form the target mode set. For the target mode set Dominant mode in Amplitude-phase decomposition was performed, and the phase characteristics of each unit in the dominant mode were considered. Construct the first Phase characteristic vector of each unit: ； in, Indicating a fault scenario Lower unit In the selected dominant mode set Phase signature vector on; Indicating a fault scenario Lower unit In the The initial phase on each dominant mode.
6. The unit combination method based on Koopman partitioning considering regional inertia constraints as described in claim 5, characterized in that: The method for forming homogeneous groups and constructing candidate homogeneous partitions through k-means clustering in step S3 is as follows: For any fault scenario Given the number of partitions ,by As the input for clustering, the K-means method is used to solve the following objective function to obtain the cluster centers. and the corresponding homology region: ； in, Indicating a fault scenario Next The set of samples corresponding to each cluster Indicating a fault scenario Lower unit The phase eigenvector, Indicates the fault scenario Next The center vector of each cluster; Represents the L2 norm, Indicates sample and its cluster center The square of the Euclidean distance between them; Each failure scenario The clustering results below are recorded as a partition result: ； in, Indicates the fault scenario The next One homogeneous region; The single-partition results corresponding to all failure scenarios are summarized to form a partition library: 。 7. The unit combination method based on Koopman partitioning considering regional inertia constraints as described in claim 1, characterized in that: The specific method for calculating the coverage rate in step S4 is as follows: ； in, Indicates candidate homology region coverage, For partitioned libraries, For any candidate homology region, Indicates partitioned library The total number of results in the middle partition. This is an indicator function: it takes the value 1 when the event is true, and 0 otherwise.
8. The unit combination method based on Koopman partitioning considering regional inertia constraints according to claim 1, characterized in that: In step S5, the frequency change rate constraint is constructed based on the magnitude of the regional disturbance as follows: ； in Indicates the unit During the period The power-on status, The system's rated frequency, This is the upper limit of the allowable rate of change of frequency. For the region The power base value, Obtained from offline perturbation statistics Percentage power imbalance value; Indicates the unit During the period The startup status variable, when the unit During the period The value is 1 when the device is powered on, otherwise it is 0. Indicates the unit inertia constant; A set of representative homology regions; The specific method for constructing the maximum frequency deviation constraint based on the magnitude of regional disturbances is as follows: ； in This is the load damping coefficient. Obtained from offline perturbation statistics quantile power imbalance value This is the equivalent droop parameter for primary frequency modulation. The equivalent coefficient parameters for a primary frequency modulation channel are: The equivalent time constant parameter for primary frequency modulation; and These represent the damping ratio and natural angular frequency of the equivalent second-order frequency response model, respectively. The time when the frequency is at its lowest point. It is the exponential factor of the second-order response decay term. Indicates the unit During the period The startup status variable, when the unit During the period The value is 1 when the device is powered on, otherwise it is 0. Indicates the unit The inertia constant.
9. The unit combination method based on Koopman partitioning considering regional inertia constraints according to claim 1, characterized in that: The objective function to be minimized in step S6 is specifically: ； in, The coefficients of the cost function for the k-th unit are represented. This indicates the on / off state of the k-th unit during time period t. This represents the power generation of the k-th generating unit during time period t; Classic thermal power unit operating constraints include power balance constraints, spinning reserve constraints, minimum start-up / shutdown time constraints, start-up cost constraints, and ramp rate constraints. The power balance constraint is specifically as follows: ； in Indicates time period The total system load requirement; The rotational spare constraint is specifically as follows: ； in, This represents the rotating reserve requirement during time period t. This represents the maximum output of unit k; The minimum power-on / power-off time constraint is specifically as follows: ； in, Indicates the unit Minimum continuous power-on duration; Indicates the unit During the period The startup status variable, if the unit is in the time period If enabled, set the value to 1; otherwise, set the value to 0. Indicates the unit The minimum continuous downtime; Indicates the unit During the period The shutdown status variable, if the unit is in the time period If the machine is stopped, use 1; otherwise, use 0. The specific startup cost constraint is as follows: ； in, Indicates the unit Start-up and shutdown costs during time period t; Indicates the unit Startup costs; Indicates the unit Downtime costs; Indicates the unit The start indicator variable for time period t is set to 1 when start is initiated, and 0 otherwise. Indicates the unit The shutdown indicator variable in time period t is set to 1 when shutdown occurs, and 0 otherwise. The climbing rate constraint is specifically as follows: ； in, Indicates the unit The upper limit of the upward climbing rate; Indicates the unit The upper limit of the downward climbing rate.