High-proportion wind power system transient voltage stability constraint optimal power flow dynamic dimension reduction method

Abstract

The invention belongs to the technical field of power systems of high-proportion wind power, and particularly relates to a high-proportion wind power system transient voltage stability constraint optimal power flow dynamic dimension reduction method, which comprises the following steps: S1, respectively carrying out collaborative dimension reduction on a bus level and a device level on the basis of an original high-order dynamic model to obtain a low-dimension equivalent model; s2, judging whether a preset transient voltage precision requirement is met or not; s3, equivalently decomposing the transient voltage stability constraint optimal power flow problem corresponding to the low-dimensional equivalent model meeting the precision requirement into two sub-models for coupling iteration solution; s4, performing coupling iteration solution on the two sub-models in the S3; and in the solving process, a self-adaptive correction mechanism is adopted for the transient voltage sensitivity matrix: according to the variable quantity of a system state variable between two continuous iterations, a tangent approximation mode or a tangent approximation mode is dynamically selected to update the transient voltage sensitivity matrix. According to the method, the reasonable voltage / reactive power level can be accurately and quickly evaluated.

Description

Technical Field

[0001] This invention belongs to the field of power system technology with high proportion of wind power, and particularly relates to a dynamic dimensionality reduction method for optimal power flow under transient voltage stability constraints in high proportion of wind power systems. Background Technology

[0002] With the large-scale integration of renewable energy sources such as wind power and power electronic interface equipment into modern power systems, the dynamic characteristics of the power grid have undergone profound changes. High wind power penetration significantly weakens the system's inertia support and voltage regulation capabilities, making transient voltage stability issues increasingly prominent. In the event of disturbances such as short circuits or line tripping, the local bus voltage may drop sharply within milliseconds. If reactive power support is insufficient or control strategies are lagging, this can easily trigger cascading trips or even large-scale power outages. Therefore, accurately simulating the transient process of the entire system during voltage / reactive power assessment is crucial.

[0003] However, power systems with a high proportion of wind power are typically characterized by their large scale, highly coupled differential-algebraic equations, dynamic interweaving across multiple time scales, and diverse fault scenarios. Their original dynamic models often contain thousands or even tens of thousands of state variables (such as the internal control loop of the wind turbine converter, the rotor motion of the synchronous machine, and the excitation system). When directly used to construct the Transient Voltage Stability-Constrained Optimal Power Flow (TVS-OPF) problem, they will face the "curse of dimensionality": on the one hand, the coupling between nonlinear differential equations and discrete control variables (such as capacitor switching and OLTC tap positions) makes the problem difficult to solve; on the other hand, the cost of detailed time-domain simulation calculations to meet engineering accuracy requirements is extremely high, making it difficult to meet the real-time requirements of online scheduling.

[0004] The existing technologies mainly suffer from the following bottlenecks: (1) The contradiction between model complexity and computational efficiency is prominent. Traditional methods either adopt simplified dynamic models (such as only retaining the dynamics of the synchronous machine) and ignore the fast response characteristics of key components such as wind power converters, resulting in overly optimistic safety assessment results; or directly use full-order models for joint simulation-optimization, which takes too long to calculate and cannot be put into practical use. (2) The handling of safety constraints is crude. Most studies simplify transient voltage constraints to static N-1 criteria or steady-state voltage limits, which fail to truly reflect the dynamic evolution of voltage after a fault, especially in weak power grids or high wind power penetration scenarios, which are prone to misjudgment. (3) Dimensionality reduction and optimization are disconnected. Although some dimension reduction methods can compress the model size, they do not consider the structural characteristics of subsequent optimization problems (such as mixed integer characteristics and sensitivity dependence), resulting in a "small" but "inaccurate" model after dimension reduction, or an inability to effectively embed it into the optimization framework.

[0005] Therefore, how to accurately and quickly assess reasonable voltage / reactive power levels has become an urgent problem to be solved. Summary of the Invention

[0006] To address the shortcomings of the existing technologies, this invention provides a dynamic dimensionality reduction method for optimal power flow under transient voltage stability constraints in high-proportion wind power systems, which can accurately and quickly assess reasonable voltage / reactive power levels.

[0007] To solve the above-mentioned technical problems, the present invention adopts the following technical solution:

[0008] A dynamic dimensionality reduction method for optimal power flow under transient voltage stability constraints in high-proportion wind power systems includes the following steps:

[0009] S1. For power systems with a high proportion of wind power, based on the original high-order dynamic model containing differential equations and network algebra equations of wind turbine generators and synchronous machines, synergistic dimensionality reduction is carried out at the bus level and device level respectively to obtain a low-dimensional equivalent model that retains the dynamic modes that play a dominant role in transient voltage stability.

[0010] At the bus level, voltage safety constraints are linearized based on trajectory sensitivity analysis, and redundant operating constraints are identified and eliminated based on the dynamic significance of bus voltage during transient processes, thereby achieving compression of the bus dimension.

[0011] At the device level, based on the empirical Gramian matrix decomposition method, the state and output data obtained by time-domain simulation are used to project the high-order differential equation model of wind turbine and synchronous machine to a low-dimensional subspace to generate an equivalent low-order dynamic model.

[0012] S2. Evaluate the dimensionality reduction error of the low-dimensional equivalent model obtained in S1 to determine whether it meets the preset transient voltage accuracy requirements. If it does not meet the requirements, adjust the dimensionality reduction parameters and return to S1 to re-execute the dimensionality reduction until a low-dimensional equivalent model that meets the accuracy requirements is obtained.

[0013] S3. The transient voltage stability constrained optimal power flow problem corresponding to the low-dimensional equivalent model that meets the accuracy requirements is equivalently decomposed into two coupled iteratively solved sub-models:

[0014] One sub-model is a mixed-integer linear programming model, which establishes a linearized dynamic relationship between state variables and continuous and discrete control variables through the transient voltage sensitivity matrix, so as to minimize system operating costs and optimize reactive power resource allocation;

[0015] Another sub-model is the combined correction model, which is based on transient simulation and steady-state power flow equations and is used to verify whether the control scheme output by the aforementioned mixed integer linear programming model meets the transient voltage safety constraints.

[0016] The two sub-models are coupled and solved through the iterative transfer of the sensitivity matrix and continuous and discrete control variables, forming a closed-loop optimization process;

[0017] S4. Solve the two sub-models of S3 in a coupled iterative manner. During the solution process, an adaptive correction mechanism is adopted for the transient voltage sensitivity matrix: based on the change of system state variables between two consecutive iterations, the transient voltage sensitivity matrix is ​​updated by dynamically selecting the secant approximation or tangent approximation method, thereby balancing computational efficiency and numerical accuracy under different operating conditions and achieving accelerated convergence.

[0018] Compared with the prior art, the present invention has the following beneficial effects:

[0019] 1. Significantly reduces model complexity while preserving key dynamic characteristics. By implementing collaborative dimensionality reduction at both the bus and device levels, redundant bus constraints with negligible impact on transient voltage stability are eliminated using trajectory sensitivity and convex polyhedron theory. Furthermore, the low-order dynamic modes dominating transient voltage responses in wind turbines and synchronous machines are accurately extracted based on the empirical Gramian method. Compared to existing technologies that only employ static simplification or full-order simulation, this method effectively ensures dynamic fidelity while significantly reducing model size.

[0020] 2. Achieve controllable closed-loop feedback for dimensionality reduction accuracy. A dimensionality reduction evaluation mechanism based on output response error is introduced, and the dimensionality reduction parameters are automatically adjusted and remodeled when accuracy requirements are not met, forming a closed-loop process of "dimensionality reduction—evaluation—correction". Compared to traditional one-time dimensionality reduction methods, this mechanism avoids safety misjudgments caused by oversimplification, ensuring that the low-dimensional model always meets the engineering accuracy requirements for transient voltage safety analysis.

[0021] 3. Effective integration of mixed-integer optimization and transient safety verification. The original nonlinear, non-convex transient voltage stability constrained optimal power flow problem is decomposed into an efficient solvable mixed-integer linear programming (MILP) sub-model and a high-fidelity combined correction sub-model, which are iteratively coupled through control variables and sensitivity matrices. This architecture overcomes the difficulty of existing methods in simultaneously handling discrete control actions and continuous transient simulations, balancing the solvability of optimization with the accuracy of safety verification.

[0022] 4. Improve iterative convergence efficiency and numerical robustness. An adaptive sensitivity correction mechanism is introduced into the coupled solution process, dynamically switching between secant and tangent approximation update strategies based on changes in state variables. Compared to traditional methods that use a single fixed approximation method, this mechanism accelerates convergence during large step phases and ensures accuracy during fine-tuning phases, significantly enhancing the algorithm's adaptability and stability under different operating conditions.

[0023] In summary, this method can accurately and quickly assess reasonable voltage / reactive power levels.

[0024] Preferably, in S1, at the bus level, the process of linearizing the voltage safety constraint based on trajectory sensitivity analysis includes:

[0025] Calculate the sensitivity of the bus voltage to the reactive power injection from the generator at time t under fault scenario s.

[0026] in, The change in bus voltage caused by the change in reactive power output of the generator in the i-th interval; Let be the change in reactive power output of the generator in the i-th interval. This represents the lower limit or reference starting value of the generator reactive power output in the i-th interval. This is the reactive power setting value for the (i+1)th interval;

[0027] This sensitivity is used to transform the nonlinear voltage safety constraint into a linear form:

[0028]

[0029] In the formula, Let be the predicted voltage value of the (i+1)th interval at time t in fault scenario s. This represents the current voltage value of the i-th interval within the same scenario; V represents the actual reactive power output of the generator in the i-th iteration; max and V min These represent the upper and lower limits of the voltage value, respectively.

[0030] This setup allows for two key advantages: 1) Achieving a dynamic linear expression of transient voltage constraints. Traditional methods often treat voltage safety constraints as static or simplify them to steady-state limits, failing to reflect the true evolution of voltage over time after a fault. This approach transforms the originally highly nonlinear time-varying voltage constraints into a form that can be handled within a linear programming framework by calculating the voltage's sensitivity to reactive power changes at each moment and in each iteration. Compared to existing static or coarse approximation methods, this technique significantly improves the modeling accuracy and dynamic response matching capability of voltage safety constraints.

[0031] 2. Effectively identify and eliminate redundant operational constraints to reduce problem dimensionality. In transient analyses across multiple scenarios and time steps, many buses may experience small voltage fluctuations under most operating conditions, and their constraints have a limited impact on the overall system safety. This solution utilizes sensitivity indicators to assess the "dynamic significance" of each bus voltage, and then combines this with convex polyhedron theory to determine which constraints can be safely ignored. This mechanism avoids full simulation and constraint modeling for all buses, thereby significantly reducing the number of optimization variables and computational complexity.

[0032] Preferably, in S1, the process of identifying and eliminating redundant operating constraints based on the dynamic saliency of the bus voltage during the transient process, combined with convex polyhedral theory, includes:

[0033] Based on the convex polyhedron theory, a feasible region for the control input is constructed. By analyzing the dynamic significance of the bus voltage during the transient process, redundant operational constraints unrelated to transient voltage stability are identified and eliminated. The dynamic significance is determined by the voltage change at adjacent sampling times, and its constraint condition is as follows:

[0034]

[0035] In the formula, Let be the voltage response value of the original high-order dynamic model of bus b at time t under fault scenario s. To predict the voltage value for the same fault scenario s at the same time t using the currently obtained low-dimensional equivalent model; ε is the preset voltage fluctuation allowable error; t cl t is the fault clearing time. f N is the simulation termination time; s This represents the number of fault scenarios.

[0036] If a busbar b satisfies the above inequality in all scenarios and transient periods, its voltage dynamics are deemed negligible, and the corresponding safety operation constraints are eliminated.

[0037] This setup allows for accurate identification of non-dominant dynamic buses, simplifying constraints. In high-proportion wind power systems, many buses experience minimal voltage fluctuations during transient processes, and their safety constraints contribute little to overall system stability. Traditional methods often include all buses in the constraint set, leading to an inflated problem scale. This approach, by quantitatively analyzing the dynamic changes in bus voltage under different fault scenarios, retains only buses with significant dynamic behavior for modeling, effectively reducing the number of redundant constraints and improving the model's simplicity and tractability.

[0038] Preferably, in S1, generating an equivalent low-order dynamic model at the device level based on the empirical Gramian matrix decomposition method includes:

[0039] Trajectory data of the system state variable x(t) and output variable y(t) were collected through time-domain simulation under various input excitations, and empirical Gramian matrix P was constructed using these data. e ;

[0040] Then, the original high-order model is projected onto a low-dimensional subspace to obtain the dimensionality-reduced system matrix:

[0041]

[0042] In the formula, P is the projection matrix; T is the coordinate transformation matrix; f() and h() are the right-hand side terms and output equations of the original system in the original high-order dynamic model, respectively; These are the output variables after dimensionality reduction; n represents the state variables after dimensionality reduction. x denoted by , where is the order of the low-dimensional equivalent model; u(t) is the external controllable input vector applied to the system at time t.

[0043] This setup preserves the dominant dynamic modes and improves the physical fidelity of the low-dimensional model. The empirical Gramian method automatically identifies the dynamic modes (such as fast transient responses and oscillating modes) that dominate the system's behavior by analyzing the system's energy response to input excitation. Compared to traditional simplified models (such as ignoring the internal control loop of the wind turbine), this method effectively preserves the key dynamic processes that affect transient voltage stability, ensuring that the model after dimensionality reduction still accurately reflects the dynamic characteristics of the original system.

[0044] Preferably, in S2, the dimensionality reduction error evaluation includes:

[0045] Construct the output response error function between the original high-order dynamic model and the low-dimensional equivalent model:

[0046]

[0047] in,

[0048] In the formula, y s (t) represents the output response of the original high-order dynamic model under fault s. Output response for the low-dimensional equivalent model; W y The weighted matrix for the output variables; N s t represents the number of fault scenarios. f n is the simulation termination time; x N represents the order of the low-dimensional equivalent model. x u is the order of the original high-order dynamic model; s , Let s be the control input vector applied to the system under fault scenario s (this input acts on both the original high-order dynamic model and the low-dimensional equivalent model); y These are the upper and lower limits of the output, respectively; x These are the upper and lower limits of the state, respectively; the objective is to minimize the integral of the error.

[0049] And based on this error function, it is determined whether the low-dimensional equivalent model meets the preset transient voltage accuracy requirements. If L(n) x If the value is greater than a preset threshold, then the order n of the low-dimensional equivalent model will be reduced. xIncrease by 1 and return to S1 to reduce the dimension again.

[0050] This setup allows for: 1) quantitative evaluation and closed-loop control of dimensionality reduction accuracy. Traditional dimensionality reduction methods often rely on empirically set orders or subjective judgments, making it difficult to guarantee the model's accuracy during critical transient processes. This scheme constructs an integral error function L(n... x This method unifies and quantifies the output differences between the low-dimensional model and the original model under multi-fault scenarios, thereby objectively evaluating the quality of dimensionality reduction. Compared with qualitative analysis, this method has stronger scientific rigor and reproducibility.

[0051] 2. Ensure the low-dimensional model meets the engineering accuracy requirements for transient voltage safety analysis. The error function focuses on key electrical quantities such as voltage and covers multiple fault scenarios and the entire simulation period, comprehensively capturing the deviations of the low-dimensional model in dynamic response. When L(n x When the threshold is exceeded, the system automatically increases the order and remodels, avoiding the risk of safety misjudgment caused by oversimplification and ensuring the reliability of constraints in subsequent optimal power flow problems.

[0052] Preferably, in S2, the implicit trapezoidal integral method is used to perform numerical simulations on the original high-order dynamic model and the low-dimensional equivalent model, and its discretization form is as follows:

[0053] [x(t+1)-x(t)] / h=[f(x(t+1),u(t+1))+f(x(t),u(t))] / 2;

[0054] 0=[y(t+1)-h(x(t+1))+y(t)-h(x(t))] / 2,t=1,2,...,t f ;

[0055] Where h is the integration step size; x(t) is the state variable vector of the power system at time t; u(t) is the control input vector of the system at time t; f(x(t),u(t)) is the right-hand side term of the system differential equation, representing the time rate of change of the state variables; y(t) is the output variable vector of the original high-order model at time t; h(x(t)) is the output equation, which maps the state variable x(t) to the output y(t).

[0056] This setup, by uniformly employing the same numerical integration method, allows for a fair comparison of the output responses of the original high-order model and the low-dimensional equivalent model. This avoids misjudgments caused by differences in integration algorithms, ensuring that the subsequently constructed error function L(n) is accurate and reliable. x It is more objective and comparable, enhancing the scientific basis for dimensionality reduction evaluation.

[0057] Preferably, in S3, the mixed-integer linear programming model aims to minimize the system's reactive power scheduling cost:

[0058]

[0059] In the formula, u is a continuous control variable, u k+1 Let r(u) be the value of this variable in the (k+1)th iteration. k+1 () is the running cost function associated with the control input of this iteration;

[0060] The constraints of the mixed-integer linear programming model include:

[0061]

[0062] In the formula, the subscript k / k+1 indicates the k / k+1th iteration; x k c is a state variable. k For discrete control variables; and These are the sensitivity matrices for the state with respect to continuous and discrete control variables, respectively, used to linearize the dynamic relationship.

[0063] This setup offers several advantages: 1. It enables linear modeling of transient dynamic processes, improving the solvability of the optimization problem. The original transient processes of power systems are described by complex nonlinear differential equations, making them difficult to directly embed into the optimization framework. This scheme introduces a transient voltage sensitivity matrix... and By representing the changes in state variables as a linear combination of changes in continuous and discrete control variables, this method successfully transforms high-dimensional nonlinear dynamic relationships into low-order linear expressions. Compared to traditional static or simplified models, this method significantly improves the clarity of the mathematical structure and the feasibility of solving the optimization problem while preserving its dynamic characteristics.

[0064] 2. In actual power grids, reactive power regulation resources include both continuously adjustable equipment (such as generators and STATCOMs) and discrete-action equipment (such as capacitor bank switching and transformer tap adjustment). This model considers both types of variables simultaneously and uses sensitivity matrices to characterize their impact on the system state, achieving unified modeling and joint optimization of multiple types of regulation methods, thus enhancing the flexibility and practicality of scheduling strategies.

[0065] Preferably, in S3, the mathematical expression of the combined correction model is:

[0066]

[0067] In the formula, F() is the implicit nonlinear mapping function of the combined correction model; and Let u be the initial state variable and algebraic variable for the (k+1)th iteration, respectively; k+1 c is the continuous control variable vector for the (k+1)th iteration;k+1 h is the vector of discrete control variables in the (k+1)th iteration. eq () = 0 is the steady-state equilibrium equation, used to determine the initial point before the fault; the superscript s indicates the fault s; f s () represents the right-hand side term of the transient differential equation; g s ()=0 is an algebraic constraint equation describing the relationship between the algebraic variables and state variables of the system; h s ()≤0 is a transient safety constraint; t∈(t cl ,t f [N] represents the time interval from when the fault is cleared to when the simulation ends; s This represents the number of fault scenarios.

[0068] Compared to traditional static or simplified models that rely solely on steady-state voltage limits to determine safety, this approach solves a complete system of nonlinear differential-algebraic equations (DAEs) to accurately simulate the entire system response from fault occurrence to clearance under multi-fault scenarios. This allows the calibration model to realistically reflect key dynamic characteristics such as voltage dip depth and recovery speed, significantly improving the accuracy and reliability of safety assessments.

[0069] Preferably, in S4, the adaptive correction mechanism used for the transient voltage sensitivity matrix is ​​as follows:

[0070]

[0071] In the formula, and The transient voltage sensitivity matrices for the k-th and (k+1)-th iterations, respectively; x k x k+1 x k-1 These are the state variables for the k-th, k+1-th, and k-1-th iterations, respectively; u k and u k+1 , respectively, are the control variable vectors for the k-th and k+1-th iterations; F() is the implicit nonlinear mapping function of the combined correction model, representing the response from control input to state output; Δu is the disturbance of the control variable; ε is a preset small positive threshold used to distinguish between small step size and large step size regions; This represents the Euclidean norm; the superscript i is the index. This represents the i-th control variable in the k-th iteration; the above segmentation conditions are based on the state change. With control change The relative size of the values ​​enables automatic switching between the first row secant approximation and the second row tangent approximation.

[0072] This setup achieves two key advantages: 1. Intelligent switching of sensitivity update methods, improving algorithm robustness. During optimization iterations, the system may undergo large step size adjustments (e.g., in the initial stage) and small step size fine-tuning (e.g., in the convergence stage). Traditional methods consistently use a single approximation method (e.g., always using the secant method), which is prone to numerical noise or divergence risks in the small step size region. This scheme automatically determines whether the system is in a "large step size" or "small step size" region by comparing the changes in the state and control variables, and dynamically switches between secant and tangent approximations, significantly enhancing the algorithm's adaptability and stability in different operating stages.

[0073] 2. Balancing computational efficiency and numerical accuracy. The secant approximation relies solely on existing data, requiring no additional simulations, resulting in low computational costs and suitability for rapid iteration; while the tangent approximation, although requiring the introduction of perturbations for simulation, offers higher accuracy and is suitable for highly sensitive regions. This mechanism achieves an optimal balance between efficiency and accuracy by ensuring accuracy in critical regions while avoiding high-cost global computations, outperforming fixed-update strategies in single-mode approaches. Attached Figure Description

[0074] To make the objectives, technical solutions, and advantages of the invention clearer, the invention will now be described in further detail with reference to the accompanying drawings, wherein:

[0075] Figure 1 This is a flowchart of the method;

[0076] Figure 2 This is a schematic diagram showing the voltage and reactive power of synchronous generator 1 in the IEEE-9 bus system of Example 2.

[0077] Figure 3 This is a schematic diagram comparing the order reduction error of the IEEE-9 bus system in Example 2;

[0078] Figure 4 This is a schematic diagram of the voltage of wind turbine generator 1 in the IEEE-118 bus system under different control modes in Example 2;

[0079] Figure 5 This is a schematic diagram comparing the order reduction errors of the four test systems in Example 2;

[0080] Figure 6 This is a schematic diagram illustrating the convergence of the two methods in Example 2 under different test systems;

[0081] Figure 7 This is a schematic diagram of the dynamic reactive power reserve margin assessment results of the IEEE-39 bus system in Example 2.

[0082] Figure 8 This is a schematic diagram showing the capacitor bank optimization results of two methods under the IEEE-39 bus system in Example 2;

[0083] Figure 9 This is a simulation verification diagram of different methods under different conditions in Example 2;

[0084] Figure 10 This is a schematic diagram of the dynamic reactive power reserve assessment results of different systems under different wind power penetration rates in Example 2;

[0085] Figure 11 This is a schematic diagram of the simulation verification of the IEEE-118 bus system with a wind power penetration rate of 40% in Example 2. Detailed Implementation

[0086] The following detailed explanation illustrates the specific implementation methods:

[0087] Example 1

[0088] like Figure 1 As shown in the figure, this embodiment discloses a dynamic dimensionality reduction method for optimal power flow under transient voltage stability constraints in high-proportion wind power systems, including the following steps:

[0089] S1. For power systems with a high proportion of wind power, based on the original high-order dynamic model containing differential equations and network algebra equations of wind turbine generators and synchronous machines, synergistic dimensionality reduction is carried out at the bus level and device level respectively to obtain a low-dimensional equivalent model that retains the dynamic modes that play a dominant role in transient voltage stability.

[0090] At the bus level, voltage safety constraints are linearized based on trajectory sensitivity analysis, and redundant operating constraints are identified and eliminated based on the dynamic significance of bus voltage during transient processes, thereby achieving compression of the bus dimension.

[0091] At the device level, based on the empirical Gramian matrix decomposition method, the high-order differential equation model of the wind turbine and the synchronous machine is projected to a low-dimensional subspace using the state and output data obtained from time-domain simulation, generating an equivalent low-order dynamic model.

[0092] It should be noted that the nonlinear constraints can be simplified to a general linear programming form, the mathematical expression of which is:

[0093]

[0094] Ω O Ω represents the original set of nonlinear constraints. L S is the set of constraints after linearization. i Let u be the sensitivity matrix. 0 and u 1 These are the initial and updated control input vectors, respectively; C1, C2, and C3 are the steady-state operation constraints, network topology and power flow balance constraints, and transient safety constraints, respectively; x0 x 1 The initial control inputs u are respectively 0 The corresponding system state vector, and the system state vector corresponding to the updated control input u1; u i u i+1 These are the control input vectors used for optimization in the i-th and (i+1)-th iterations, respectively;

[0095] d1, d2, d3, and d4 are auxiliary decision variables used to distinguish different constraint types and processing methods. d1 = d2 = d3 = 0 indicates that the first three types of constraints (such as steady state, power flow, and transient) are not introduced as slack variables in the current linearization stage and are kept strictly satisfied. d4 ≥ 0 are optional slack variables used to handle soft constraints that are difficult to fully satisfy (such as some non-critical voltage fluctuations), allowing a certain deviation to improve the feasibility of the solution.

[0096] In practice, at the bus level, the process of linearizing voltage safety constraints based on trajectory sensitivity analysis includes:

[0097] Calculate the sensitivity of the bus voltage to the reactive power injection from the generator at time t under fault scenario s.

[0098] in, The change in bus voltage caused by the change in reactive power output of the generator in the i-th interval; Let be the change in reactive power output of the generator in the i-th interval. This represents the lower limit or reference starting value of the reactive power output of the i-th generator. This is the reactive power setting value for the (i+1)th interval;

[0099] This sensitivity is used to transform the nonlinear voltage safety constraint into a linear form:

[0100]

[0101] In the formula, Let be the predicted voltage value of the (i+1)th interval at time t in fault scenario s. This represents the current voltage value of the i-th interval within the same scenario; V represents the actual reactive power output of the generator in the i-th iteration; max and V min These represent the upper and lower limits of the voltage value, respectively.

[0102] Traditional methods often treat voltage safety constraints as static or simplify them to steady-state limits, failing to reflect the true evolution of voltage over time after a fault. This scheme transforms the originally highly nonlinear time-varying voltage constraints into a form that can be handled within a linear programming framework by calculating the voltage's sensitivity to reactive power changes at each moment and in each iteration. Compared to existing static or coarse approximation methods, this technique significantly improves the modeling accuracy and dynamic response matching capability of voltage safety constraints. Furthermore, in transient analyses across multiple scenarios and time steps, many buses may experience small voltage fluctuations under most operating conditions, limiting the impact of their constraints on the overall system safety. This scheme utilizes sensitivity indices to assess the "dynamic significance" of each bus voltage, and then combines this with convex polyhedral theory to determine which constraints can be safely ignored. This mechanism avoids full simulation and constraint modeling of all buses, thereby significantly reducing the number of optimization variables and computational complexity.

[0103] In specific implementation, the process of identifying and eliminating redundant operating constraints based on the dynamic saliency of the bus voltage during the transient process, combined with convex polyhedral theory, includes:

[0104] Based on the convex polyhedron theory, a feasible region for the control input is constructed. By analyzing the dynamic significance of the bus voltage during the transient process, redundant operational constraints unrelated to transient voltage stability are identified and eliminated. The dynamic significance is determined by the voltage change at adjacent sampling times, and its constraint condition is as follows:

[0105]

[0106] In the formula, Let be the voltage response value of the original high-order dynamic model of bus b at time t under fault scenario s. To predict the voltage value for the same fault scenario s at the same time t using the currently obtained low-dimensional equivalent model; ε is the preset voltage fluctuation allowable error; t cl t is the fault clearing time. f N is the simulation termination time; s This represents the number of fault scenarios.

[0107] If a busbar b satisfies the above inequality in all scenarios and transient periods, its voltage dynamics are deemed negligible, and the corresponding safety operation constraints are eliminated.

[0108] In high-proportion wind power systems, a large number of buses exhibit minimal voltage fluctuations during transient processes, and their safety constraints contribute little to the overall system stability. Traditional methods often include all buses in the constraint set, leading to an inflated problem size. This approach quantifies the dynamic changes in bus voltage under different fault scenarios, retaining only buses with significant dynamic behavior for modeling. This effectively reduces the number of redundant constraints and improves the model's simplicity and tractability.

[0109] In practical implementation, generating an equivalent low-order dynamic model at the device level based on the empirical Gramian matrix decomposition method includes:

[0110] Trajectory data of the system state variable x(t) and output variable y(t) were collected through time-domain simulation under various input excitations, and empirical Gramian matrix P was constructed using these data. e ;

[0111] Then, the original high-order model is projected onto a low-dimensional subspace to obtain the dimensionality-reduced system matrix:

[0112]

[0113] In the formula, P is the projection matrix; T is the coordinate transformation matrix; f() and h() are the right-hand side terms and output equations of the original system in the original high-order dynamic model, respectively; These are the output variables after dimensionality reduction; n represents the state variables after dimensionality reduction. x denoted by , where is the order of the low-dimensional equivalent model; u(t) is the external controllable input vector applied to the system at time t.

[0114] In this way, the empirical Gramian method automatically identifies the dynamic modes (such as fast transient response, oscillation mode, etc.) that dominate the system behavior by analyzing the system's energy response under input excitation. Compared with traditional simplified models (such as ignoring the internal control loop of the wind turbine), this method can effectively retain the key dynamic processes that affect the stability of transient voltage, ensuring that the model after dimensionality reduction can still accurately reflect the dynamic characteristics of the original system.

[0115] In practice, the dimensionality-reduced state variables can be divided into primary modes. and residual modes N x Let be the dimension of the original high-order dynamic model, and its dynamic equation is:

[0116]

[0117]

[0118] Where P and Q are the projection matrices of the principal mode and the residual mode, respectively.

[0119] S2. Evaluate the dimensionality reduction error of the low-dimensional equivalent model obtained in S1 to determine whether it meets the preset transient voltage accuracy requirements. If it does not meet the requirements, adjust the dimensionality reduction parameters and return to S1 to re-execute the dimensionality reduction until a low-dimensional equivalent model that meets the accuracy requirements is obtained.

[0120] In specific implementation, the dimensionality reduction error assessment includes:

[0121] Construct the output response error function between the original high-order dynamic model and the low-dimensional equivalent model:

[0122]

[0123] in,

[0124] In the formula, y s (t) represents the output response of the original high-order dynamic model under fault s. Output response for the low-dimensional equivalent model; W y The weighted matrix for the output variables; N s t represents the number of fault scenarios. f n is the simulation termination time; x N represents the order of the low-dimensional equivalent model. x u is the order of the original high-order dynamic model; s , Let s be the control input vector applied to the system under fault scenario s (this input acts on both the original high-order dynamic model and the low-dimensional equivalent model); y These are the upper and lower limits of the output, respectively; x These are the upper and lower limits of the state, respectively; the objective is to minimize the integral of the error.

[0125] And based on this error function, it is determined whether the low-dimensional equivalent model meets the preset transient voltage accuracy requirements. If L(n) x If the value is greater than a preset threshold, then the order n of the low-dimensional equivalent model will be reduced. x Increase by 1 and return to S1 to reduce the dimension again.

[0126] Traditional order reduction methods often rely on empirically set orders or subjective judgments, making it difficult to guarantee the accuracy of the model during critical transient processes. This scheme addresses this by constructing an integral form error function L(n) x This method unifies and quantifies the output differences between the low-dimensional model and the original model under multiple fault scenarios, thereby objectively evaluating the quality of dimensionality reduction. Compared to qualitative analysis, this method has stronger scientific rigor and repeatability. Furthermore, the error function focuses on key electrical quantities such as voltage and covers multiple fault scenarios and the entire simulation period, comprehensively capturing the deviations of the low-dimensional model in dynamic response. When L(n x When the threshold is exceeded, the system automatically increases the order and remodels, avoiding the risk of safety misjudgment caused by oversimplification and ensuring the reliability of constraints in subsequent optimal power flow problems.

[0127] In practice, the implicit trapezoidal integral method is used to perform numerical simulations on the original high-order dynamic model and the low-dimensional equivalent model. Its discretization form is as follows:

[0128] [x(t+1)-x(t)] / h=[f(x(t+1),u(t+1))+f(x(t),u(t))] / 2;

[0129] 0=[y(t+1)-h(x(t+1))+y(t)-h(x(t))] / 2,t=1,2,...,t f ;

[0130] Where h is the integration step size; x(t) is the state variable vector of the power system at time t; u(t) is the control input vector of the system at time t; f(x(t),u(t)) is the right-hand side term of the system differential equation, representing the time rate of change of the state variables; y(t) is the output variable vector of the original high-order model at time t; h(x(t)) is the output equation, which maps the state variable x(t) to the output y(t).

[0131] By uniformly employing the same numerical integration method, the output responses of the original high-order model and the low-dimensional equivalent model can be compared fairly. This avoids misjudgment due to differences in integration algorithms, ensuring that the subsequently constructed error function L(n) is accurate. x It is more objective and comparable, enhancing the scientific basis for dimensionality reduction evaluation.

[0132] S3. The transient voltage stability constrained optimal power flow problem corresponding to the low-dimensional equivalent model that meets the accuracy requirements is equivalently decomposed into two coupled iteratively solved sub-models:

[0133] One sub-model is a mixed-integer linear programming model, which establishes a linearized dynamic relationship between state variables and continuous and discrete control variables through the transient voltage sensitivity matrix, so as to minimize system operating costs and optimize reactive power resource allocation;

[0134] Another sub-model is the combined correction model, which is based on transient simulation and steady-state power flow equations and is used to verify whether the control scheme output by the aforementioned mixed integer linear programming model meets the transient voltage safety constraints.

[0135] The two sub-models are coupled and solved through the sensitivity matrix and the iterative transfer of continuous and discrete control variables, forming a closed-loop optimization process.

[0136] The transient voltage sensitivity matrix refers to the partial derivative matrix of system state variables (such as bus voltage) with respect to control variables (such as generator reactive power output and capacitor connection) during the critical period after a fault. It is used to characterize the linear influence of control actions on the transient voltage response. This matrix is ​​updated in each iteration by the simulation results of the combined correction model, forming a closed loop of "simulation → correction → optimization".

[0137] In practical implementation, the mixed-integer linear programming model aims to minimize the system's reactive power scheduling cost.

[0138]

[0139] In the formula, u is a continuous control variable, u k+1 Let r(u) be the value of this variable in the (k+1)th iteration. k+1 () is the running cost function associated with the control input of this iteration;

[0140] The constraints of the mixed-integer linear programming model include:

[0141]

[0142] In the formula, the subscript k / k+1 indicates the k / k+1th iteration; x k c is a state variable. k For discrete control variables; and These are the sensitivity matrices for the state with respect to continuous and discrete control variables, respectively, used to linearize the dynamic relationship.

[0143] The transient processes of the original power system are described by complex nonlinear differential equations, making them difficult to directly embed into an optimization framework. This scheme introduces a transient voltage sensitivity matrix. and By representing the changes in state variables as a linear combination of changes in continuous and discrete control variables, this method successfully transforms high-dimensional nonlinear dynamic relationships into low-order linear expressions. Compared to traditional static or simplified models, this method significantly improves the clarity of the mathematical structure and the feasibility of solving the optimization problem while preserving dynamic characteristics. Furthermore, in actual power grids, reactive power regulation resources include both continuously adjustable equipment (such as generators and STATCOMs) and discrete-action equipment (such as capacitor bank switching and transformer tap adjustment). This model considers both types of variables simultaneously and characterizes their impact on the system state through sensitivity matrices, achieving unified modeling and joint optimization of multiple types of control measures, thus enhancing the flexibility and practicality of scheduling strategies.

[0144] In practical implementation, the mathematical expression of the combined correction model is:

[0145]

[0146] In the formula, F() is the implicit nonlinear mapping function of the combined correction model; and Let u be the initial state variable and algebraic variable for the (k+1)th iteration, respectively; k+1 c is the continuous control variable vector for the (k+1)th iteration; k+1h is the vector of discrete control variables in the (k+1)th iteration. eq () = 0 is the steady-state equilibrium equation, used to determine the initial point before the fault; the superscript s indicates the fault s; f s () represents the right-hand side term of the transient differential equation; g s ()=0 is an algebraic constraint equation describing the relationship between the algebraic variables and state variables of the system; h s ()≤0 is a transient safety constraint; t∈(t cl ,t f [N] represents the time interval from when the fault is cleared to when the simulation ends; s This represents the number of fault scenarios.

[0147] Compared to traditional static or simplified models that judge safety solely based on steady-state voltage limits, this approach solves a complete system of nonlinear differential-algebraic equations (DAEs) to accurately simulate the entire system response from fault occurrence to clearance under multi-fault scenarios. This enables the calibration model to realistically reflect key dynamic characteristics such as voltage dip depth and recovery speed, significantly improving the accuracy and reliability of safety assessments.

[0148] S4. Solve the two sub-models of S3 in a coupled iterative manner. During the solution process, an adaptive correction mechanism is adopted for the transient voltage sensitivity matrix: based on the change of system state variables between two consecutive iterations, the transient voltage sensitivity matrix is ​​updated by dynamically selecting the secant approximation or tangent approximation method, thereby balancing computational efficiency and numerical accuracy under different operating conditions and achieving accelerated convergence.

[0149] In specific implementation, the adaptive correction mechanism used for the transient voltage sensitivity matrix is ​​as follows:

[0150]

[0151] In the formula, and The transient voltage sensitivity matrices for the k-th and (k+1)-th iterations, respectively; x k x k+1 x k-1 These are the state variables for the k-th, k+1-th, and k-1-th iterations, respectively; u k and u k+1 , respectively, are the control variable vectors for the k-th and k+1-th iterations; F() is the implicit nonlinear mapping function of the combined correction model, representing the response from control input to state output; Δu is the disturbance of the control variable; ε is a preset small positive threshold used to distinguish between small step size and large step size regions; This represents the Euclidean norm; the superscript i is the index. This represents the i-th control variable in the k-th iteration; the above segmentation conditions are based on the state change. With control change The relative size of the values ​​enables automatic switching between the first row secant approximation and the second row tangent approximation.

[0152] During the optimization iteration process, the system may undergo large-step adjustments (such as in the initial stage) and small-step fine-tuning (such as in the convergence stage). Traditional methods use a single approximation method (such as always using the secant method), which is prone to numerical noise or divergence risks in the small-step region. This scheme automatically determines whether the current state is in a "large-step" or "small-step" region by comparing the changes in the state and control variables, and dynamically switches between the secant and tangent approximations, significantly enhancing the algorithm's adaptability and stability in different operating stages. Furthermore, the secant approximation only relies on existing data, requires no additional simulation, has low computational cost, and is suitable for rapid iteration; while the tangent approximation, although requiring the introduction of perturbations for simulation, has higher accuracy and is suitable for highly sensitive regions. This mechanism avoids high-cost global computation while ensuring accuracy in critical regions, achieving an optimal balance between efficiency and accuracy, which is superior to the fixed update strategy in a single mode.

[0153] By implementing collaborative dimensionality reduction at both the bus level and the device level, this method utilizes trajectory sensitivity and convex polyhedron theory to eliminate redundant bus constraints that have a negligible impact on transient voltage stability. Furthermore, it accurately extracts the low-order dynamic modes that dominate the transient voltage response in wind turbines and synchronous machines based on the empirical Gramian method. Compared to existing technologies that only employ static simplification or full-order simulation, this method significantly reduces the model size while effectively ensuring dynamic fidelity. In addition, a dimensionality reduction evaluation mechanism based on output response error is introduced, automatically adjusting the dimensionality reduction parameters and remodeling when accuracy requirements are not met, forming a closed-loop process of "dimensionality reduction—evaluation—correction." Compared to traditional one-time dimensionality reduction methods, this mechanism avoids safety misjudgments caused by oversimplification, ensuring that the low-dimensional model always meets the engineering accuracy requirements for transient voltage safety analysis.

[0154] In addition, the original nonlinear, nonconvex transient voltage stability constrained optimal power flow problem is decomposed into an efficient solvable mixed-integer linear programming (MILP) sub-model and a high-fidelity combined correction sub-model, which are iteratively coupled through control variables and sensitivity matrices. This architecture overcomes the difficulty of existing methods in simultaneously handling discrete control actions and continuous transient simulations, balancing the solvability of optimization with the accuracy of safety verification. Furthermore, an adaptive sensitivity correction mechanism is introduced during the coupled solution process, dynamically switching between secant and tangent approximation update strategies based on changes in state variables. Compared to traditional methods that use a fixed single approximation method, this mechanism accelerates convergence in the large step phase and ensures accuracy in the fine-tuning phase, significantly enhancing the algorithm's adaptability and stability under different operating conditions.

[0155] This method can accurately and quickly assess reasonable voltage / reactive power levels.

[0156] Example 2

[0157] To facilitate a better understanding of this method by those skilled in the art, the following explanation is provided.

[0158] The computational efficiency and accuracy of this method were tested in four different scale systems (i.e., IEEE-9, -39, -118 and -300 buses).

[0159] Wind power penetration level is the proportion of total wind power generation to the total power generation in the area where the system is located, initially set at 40% for the four test systems. Next, taking the IEEE-9 bus system as an example, the device dimensionality reduction process will be analyzed in detail.

[0160] Figure 2 The figure shows the voltage and reactive power of synchronous generator 1 in the IEEE-9 bus system at different orders. The fault is located at bus 6, and the original model is order 15. The results show that the 10th-order reduced model can effectively represent the original model, especially for the terminal voltage of synchronous generator 1, with an error of 0.16%. However, the 7th-order reduced model introduces larger errors, with errors of 1.39% for voltage and 7.79% for reactive power. The simulation of other synchronous generator state variables and output variables is similar to the figure above.

[0161] To demonstrate the effectiveness of the order reduction error assessment, Figure 3 The order reduction errors of the IEEE-9-bus system at different orders are presented. It can be concluded that the introduced order reduction error gradually decreases as the order of the reduced model increases. However, in some cases, the error of a higher-order reduced model may be greater than that of a lower-order model. For example, the voltage error of the 8th-order synchronous generator 1 is 2.37%, higher than the 1.51% voltage error of the 7th-order synchronous generator 1. Therefore, selecting an appropriate order reduction model is crucial for calculating errors and subsequent optimization solutions. Here, the error threshold is set to 0.05, and using a 10th-order model to represent the original model satisfies the accuracy requirements well.

[0162] Similarly, the reduced-order model of the wind turbine was validated in the IEEE-118 bus system. The original IEEE-118 bus system model had an order of 270; reduced-order models of orders 175, 165, and 160 were selected to capture the voltage of the original wind turbine model in the system, such as... Figure 4 As shown, with the increase of the value of k, the constructed continuous function gets closer and closer to the original control switching function. Therefore, the value of k in the continuous function is set to 40. The 160th and 165th order models are used to characterize the wind turbine model under constant power control, with errors of 4.05% and 3.40%, respectively. This shows that under constant voltage control and constant power control modes, the 160th and 175th order models can effectively capture the voltage of wind turbine 1.

[0163] Figure 5 The errors corresponding to the optimal reduced-order models for four test systems, obtained based on the error evaluation model, are presented. It can be concluded that the order of the reduced-order models for the four test systems is approximately 2 / 3 of the original model's order, and the reduction errors can all be controlled within the error threshold of 0.05. This ensures that the optimal low-order model, with the dominant module as the primary element, is retained to replace the original model, fully guaranteeing the accuracy of the reduction.

[0164] To verify the impact of device dimensionality reduction on the sensitivity of secant / tangent correction, the convergence of the secant approximation method and the matrix adaptive correction dynamic dimensionality reduction method of the original decomposition model were tested using four test systems of different scales.

[0165] Figure 6 The iteration counts for the two methods under different test systems are presented. In all four test systems, the matrix adaptive correction dynamic dimensionality reduction method requires significantly fewer iterations than the original method. The convergence speed of the matrix adaptive correction dynamic dimensionality reduction algorithm is significantly faster than the original algorithm. Particularly in the large IEEE-300 system, the matrix adaptive correction dynamic dimensionality reduction method requires only 14 iterations, significantly lower than the 25 iterations required by the original model, thus demonstrating the necessity of dimensionality reduction.

[0166] Figure 7 The results show the dynamic reactive power reserve evaluation of the IEEE-39 bus system before and after dimensionality reduction. The red line represents the rated maximum reactive power of the synchronous generator, and the blue line represents the maximum allowable reactive power of the synchronous generator while ensuring system transient voltage safety. The area enclosed by these two lines represents the dynamic reactive power reserve boundary. Furthermore, the maximum reactive power output of the generator before and after dimensionality reduction was calculated, with the calculation error controlled within 3%.

[0167] To verify the effectiveness of the matrix adaptive correction dynamic dimensionality reduction method in computing test systems with capacitor banks, Figure 8 This paper presents a capacitor switching scenario for the IEEE-39 bus system over a 24-hour period. In this scenario, capacitor banks C1 and C2 are installed on buses 10 and 15, each consisting of three 20Mvar capacitors. The capacitor switching behavior determined by the matrix adaptive correction dynamic dimensionality reduction method is consistent with the results of the original model, confirming the method's effectiveness in solving problems involving integer variables.

[0168] Figure 9 Simulation verification of the three methods under different test cases is presented. Since the matrix adaptive correction dynamic dimensionality reduction method does not simplify transient simulations, its time-domain voltage simulation results meet the conditions for safe voltage operation. However, compared to the other two methods, the transient voltage exceeds the safety limit. For example, for... Figure 9As shown in method B (b), the voltage of bus 6 falls below the lower limit at approximately 0.4 s. For Figure 9 As shown in (d) of Method A, the voltage on bus 52 exceeds the upper limit at approximately 0.3 s. This demonstrates that only the matrix adaptive correction dynamic dimensionality reduction algorithm can guarantee the transient voltage safety of the system.

[0169] To further verify the computational practicality of the matrix adaptive correction dynamic dimensionality reduction method under different wind power penetration rates. Figure 10 The results of total dynamic reactive power reserve assessments for different systems under varying wind power penetration rates are presented. It can be concluded that as the wind power penetration rate increases from 20% to 80%, the required dynamic reactive power reserve capacity of the system also increases. Specifically, in the IEEE-39 bus system, the penetration rate ranges from 27.78% to 41.58% to 43.71%. Due to the relatively weak grid structure and high load operation level of the IEEE-39 system, when the wind power grid-connected capacity is large, a significant amount of dynamic reactive power reserve capacity needs to be reserved from synchronous generators to ensure the safety of system transient voltage.

[0170] Furthermore, the effectiveness of the dynamic reactive power reserve was verified using the IEEE-118 system as an example. Wind turbines 1, 2, 3, and 4 of the wind farm are connected to buses 4, 12, 25, and 31 of the system, respectively. The system's dynamic reactive power reserve was calculated using a matrix adaptive correction dynamic dimensionality reduction method, assuming a wind power penetration rate of 40% in the system. Figure 11 The green bar chart is shown. The yellow bar represents the maximum reactive power output allowed by the synchronous generator during steady-state operation. Scenario A (dynamic reactive power reserve capacity of 3052.90 Mvar) and Scenario B (dynamic reactive power reserve capacity of 4052.90 Mvar) were selected to verify the transient voltage safety of the wind turbines. Scenario A: Insufficient dynamic reactive power reserve caused voltage anomalies in wind turbines 1, 2, and 3 for 5.5s to 6s. Scenario B shows that the system has sufficient dynamic reactive power reserve, and the transient voltages of all wind turbines are within the safe range under this scenario, further verifying the effectiveness of the matrix adaptive correction dynamic dimensionality reduction algorithm.

[0171] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention and not to limit the technical solutions. Those skilled in the art should understand that any modifications or equivalent substitutions to the technical solutions of the present invention without departing from the spirit and scope of the present invention should be covered within the scope of the claims of the present invention.

Claims

1. A dynamic dimension reduction method for optimal power flow with transient voltage stability constraints for high penetration wind power systems, characterized in that, The method comprises the following steps: S1, for a power system containing a high proportion of wind power, based on the original high-order dynamic model containing the differential equations of wind turbine generators and synchronous machines and network algebraic equations, dimensionality reduction is performed at the bus level and the device level to obtain a low-dimensional equivalent model that retains the dynamic mode that plays a leading role in transient voltage stability; wherein, at the bus level, the voltage safety constraint is linearized based on trajectory sensitivity analysis, and combined with the convex polyhedron theory, the dynamic significance of the bus voltage in the transient process is identified and redundant operating constraints are removed to compress the dimension of the bus; at the device level, based on the empirical Gramian matrix decomposition method, the high-order differential equation model of the wind turbine generator and the synchronous machine is projected into a low-dimensional subspace using the state and output data obtained by time domain simulation to generate an equivalent low-order dynamic model; S2, the low-dimensional equivalent model obtained in S1 is subjected to dimensionality reduction error evaluation to determine whether it meets the preset transient voltage accuracy requirement; if not, the dimensionality reduction parameters are adjusted and the dimensionality reduction is re-executed in S1 until a low-dimensional equivalent model that meets the accuracy requirement is obtained; S3, the transient voltage stability constraint optimal power flow problem corresponding to the low-dimensional equivalent model that meets the accuracy requirement is equivalent to two sub-models that are coupled and iteratively solved: one sub-model is a mixed integer linear programming model that establishes a linearized dynamic relationship between state variables and continuous control variables and discrete control variables through a transient voltage sensitivity matrix to minimize system operation cost and optimize reactive power resource allocation; the other sub-model is a combined correction model based on transient simulation and steady-state power flow equations for checking whether the control scheme output by the aforementioned mixed integer linear programming model meets the transient voltage safety constraint; the two sub-models are coupled and iteratively solved through iterative transmission of the sensitivity matrix and continuous control variables and discrete control variables to form a closed-loop optimization process; S4, the two sub-models in S3 are coupled and iteratively solved; in the solving process, an adaptive correction mechanism is used for the transient voltage sensitivity matrix: according to the variation of the system state variables between two consecutive iterations, the transient voltage sensitivity matrix is updated in a secant approximation or tangent approximation mode to balance the calculation efficiency and numerical accuracy under different operating conditions, thereby realizing convergence acceleration.

2. The high penetration wind power system transient voltage stability constrained optimal power flow dynamic dimension reduction method of claim 1, wherein: In S1, at the bus level, the process of linearizing the voltage safety constraint based on trajectory sensitivity analysis includes: Calculating the sensitivity of the bus voltage to the generator reactive power injection at time t under fault scenario s wherein, is the bus voltage variation caused by the reactive power output variation of the generator in the ith interval; is the variation of the reactive power output of the generator in the ith interval, is the lower limit value or the reference starting value of the reactive power output of the generator in the ith interval, is the reactive power setting value of the ith+1 interval; using the sensitivity to convert the nonlinear voltage safety constraint into a linear form: wherein, is the predicted voltage value for the (i+1)th interval at time t for the fault scenario s, is the current voltage value for the ith interval under the same scenario; is the actual generator reactive power output for the ith iteration; V max and V min are the upper and lower limits of the voltage value, respectively.

3. The high penetration wind power system transient voltage stability constrained optimal power flow dynamic dimension reduction method of claim 2, wherein: In S1, the process of identifying and removing redundant operating constraints based on the dynamic significance of the bus voltage in the transient process in combination with the convex polyhedron theory includes: based on the convex polyhedron theory, the control input feasible region is constructed, the dynamic significance of the bus voltage in the transient process is analyzed, and redundant operating constraints irrelevant to transient voltage stability are identified and removed; the dynamic significance is determined by the voltage variation at adjacent sampling times, and the constraint condition is: In the formula, is the voltage response value of the original high-order dynamic model of the bus b at the time t under the fault scenario s, is the voltage prediction value of the same fault scenario s and the same time t using the currently obtained low-dimensional equivalent model; ε is a preset voltage fluctuation allowable error; t cl is the fault clearing time, t f is the simulation termination time; N s is the number of fault scenarios; if a bus b satisfies the above inequality in all scenarios and transient periods, it is determined that the voltage dynamic is negligible, and the corresponding safety operating constraint is removed.

4. The high penetration wind power system transient voltage stability constrained optimal power flow dynamic dimension reduction method of claim 3, wherein: In S1, the equivalent low-order dynamic model is generated based on an empirical Gramian matrix decomposition method at the device level, including: Through time-domain simulation, the trajectory data of system state variable x(t) and output variable y(t) are collected under various input excitations, and the empirical Gramian matrix P is constructed using the data e ; The original high-order model is projected to a low-dimensional subspace to obtain a reduced system matrix: where P is the projection matrix; T is the coordinate transformation matrix; f() and h() are the right end item and output equation of the original high-order dynamic model, respectively; is the output variable after dimension reduction; is the state variable after dimension reduction; n x is the order of the low-dimensional equivalent model; u(t) is the external controllable input vector applied to the system at time t.

5. The method of claim 4, wherein the high wind power penetration system transient voltage stability constrained optimal power flow dynamic dimension reduction method is characterized by: In S2, the dimension reduction error evaluation includes: An output response error function is constructed between the original high-order dynamic model and the low-dimensional equivalent model: wherein n x ∈ [1, N x ], where y s (t) is the original high-order dynamic model output response under fault s, is the low-dimensional equivalent model output response; W y is the weighting matrix of output variables; N s is the number of fault scenarios; t f is the simulation termination time; n x is the order of the low-dimensional equivalent model; N x is the order of the original high-order dynamic model; u s , is the control input vector applied to the system under fault scenario s; y are the upper and lower limit values of the output, respectively; x are the upper and lower limit values of the state, respectively; the objective is to minimize the error integral; And based on the error function, it is judged whether the low-dimensional equivalent model meets the preset transient voltage accuracy requirement. If L(n x ) is greater than a preset threshold, the order n x of the low-dimensional equivalent model is increased by 1, and S1 is returned to reduce the dimension again.

6. The high penetration wind power system transient voltage stability constrained optimal power flow dynamic dimension reduction method of claim 5, wherein: In S2, the implicit trapezoidal integration method is used for numerical simulation of the original high-order dynamic model and the low-dimensional equivalent model, and the discrete form is: [x(t+1)-x(t)] / h=[f(x(t+1),u(t+1))+f(x(t),u(t))] / 2; 0 = [y(t+1) - h(x(t+1)) + y(t) - h(x(t))] / 2, t = 1, 2,..., t f ; Where h is the integration step size; x(t) is the state variable vector of the power system at time t; u(t) is the control input vector of the system at time t; f(x(t), u(t)) is the right-hand side of the system differential equation, representing the time rate of change of the state variable; y(t) is the output variable vector of the original high-order model at time t; h(x(t)) is the output equation that maps the state variable x(t) to the output y(t).

7. The high penetration wind power system transient voltage stability constrained optimal power flow dynamic dimension reduction method of claim 6, wherein: In S3, the mixed integer linear programming model aims to minimize the reactive power dispatching cost: where u is the continuous control variable, u k+1 is the value of the variable in the k+1 iteration, r(u k+1 ) is the operating cost function associated with the control input of this iteration. The constraints of the mixed integer linear programming model include: where the subscript k / k+1 denotes the k / k+1 iteration; x k is the state variable; c k is the discrete control variable; and are the sensitivity matrices of the state with respect to the continuous and discrete control variables, respectively, used to linearize the dynamic relationship. 8.The high penetration wind power system transient voltage stability constrained optimal power flow dynamic dimension reduction method of claim 7, wherein: In S3, the mathematical expression of the combined correction model is: where F() is the implicit nonlinear mapping function of the combined correction model; and are the initial state and algebraic variables of the (k+1)th iteration, respectively; u k+1 is the continuous control variable vector of the (k+1)th iteration; c k+1 is the discrete control variable vector of the (k+1)th iteration; h eq () = 0 is the steady-state equilibrium equation for determining the initial point before the fault; the superscript s denotes the fault s; f s () is the transient differential equation right-hand side term; g s () = 0 is the algebraic constraint equation describing the relationship between the system algebraic variables and the state variables; h s () ≤ 0 is the transient security constraint; t e (t cl , t f ] is the time interval from the fault clearing to the simulation termination; N s is the number of fault scenarios.

9. The high penetration wind power system transient voltage stability constrained optimal power flow dynamic dimension reduction method of claim 8, wherein: In S4, the adaptive correction mechanism adopted for the transient voltage sensitivity matrix is: wherein, and are the transient voltage sensitivity matrices of the kth and k+1th iterations, respectively; x k , x k+1 , x k-1 are the state variables of the kth, k+1th and k-1th iterations, respectively; u k and u k+1 are the control variable vectors of the kth and k+1th iterations, respectively; F() is an implicit nonlinear mapping function of the combined correction model, representing the response from control input to state output; Δu is the perturbation of the control variable; ε is a preset small positive threshold value for distinguishing small step and large step regions; denotes the Euclidean norm; the superscript i is an index, denotes the ith control variable of the kth iteration; the above piecewise condition realizes automatic switching of the first row secant approximation and the second row tangent approximation based on the relative size of the state change and the control change .

Images