Coordinated reactive power planning in an industrial park integrated energy system using support vector regression

In this study, we develop and present a comprehensive multi-objective optimization model for reactive power planning in industrial park integrated energy systems (IPIES), addressing the challenges posed by extensive inductive loads. Our model prioritizes enhancing system resilience through three key objectives: minimizing total costs, improving static voltage stability margins, and enhancing voltage recovery capabilities. We employ support vector regression (SVR) to estimate voltage recovery metrics effectively, using reactive power capacity as input. The optimization process leverages the Non-dominated Sorting Genetic Algorithm-II (NSGA-II) for efficient solution finding. Validation on the IEEE-39 bus system demonstrates the model’s effectiveness, with SVR significantly reducing computational demands while maintaining satisfactory accuracy.


Introduction
Expanding upon the challenges and strategies within Industrial Park Integrated Energy Systems (IPIES), it's evident that managing reactive power is a nuanced task, particularly with the prevalence of inductive loads like induction motors, electric arc furnaces, and arc welding equipment.These loads are not merely passive elements; they actively shape the energy landscape of an industrial park by significantly increasing the demand for reactive power, leading to shortages that can disrupt operations [1].Disturbances exacerbate this situation, causing sudden surges in reactive power demand as induction motors restart, which can strain the system and hinder the smooth operation of critical machinery [2].This imbalance between reactive power supply and demand not only affects the operational efficiency but also poses economic risks due to potential downtime and reduced productivity.To counteract these challenges, a meticulous approach to reactive power management is essential, with a focus on strategically deploying resources to maintain balance and ensure continuous operation.Although solutions like parallel capacitor banks offer a cost-effective way to supply reactive power, their effectiveness is often limited by slow response times, which can be inadequate in scenarios requiring rapid voltage stabilization [3].This highlights the need for a comprehensive strategy that not only addresses the immediate reactive power needs but also considers the dynamic nature of industrial energy systems.
In the realm of dynamic reactive compensation systems, the distributed static synchronous compensator (DSTATCOM) is recognized for its swift responsiveness to system reactive demands, thereby enhancing voltage recovery capabilities [3].However, the effectiveness of deploying DSTATCOM to improve voltage recovery in a system is intricately linked to factors such as its installation location and capacity allocation [4] [5].Recent scholarly endeavors have delved extensively into the configuration aspects of DSTATCOM.In [6], Khan and colleagues applied an advanced algorithm based on bacterial foraging to determine the ideal size and location for DSTATCOM in the power distribution network, with goals of reducing losses and improving voltage stability.Further investigation in [7] examined how variations in wind energy and the power demand of loads affect the voltage levels in wind power systems, both with and without the integration of DSTATCOM.A unique strategy described in [8] employs a six-pulse controller in DSTATCOM for managing reactive power in wind farms, alongside strategic capacitor bank charging, to enhance voltage stability.The research methodology applied in [6][7] [8] involves time-domain simulations, which, while providing an accurate portrayal of voltage recovery ability.However, these methods are limited by their time-costed nature in a large system.
In addressing the prolonged computational demands inherent in time-domain simulations, machine learning emerges as a compelling alternative.Among the array of machine learning paradigms, Support Vector Machine (SVM), introduced by Vapnik and colleagues in the mid-1990s, has garnered significant attention for its promise and versatility [9].Grounded in statistical learning theory, SVM has achieved widespread acclaim for its outstanding performance in tasks involving classification and regression.In the context of power systems, SVM proves to be particularly pertinent.Notably, researchers in [10] harnessed Support Vector Regression (SVR), a variant of SVM, to forecast the static voltage stability of power systems.This involved utilizing training data derived from PV or QV curve data of the system.The trained SVR model adeptly evaluates the stability of hypothetical operating states within the system.Furthermore, An et al. [11] applied SVR for fault identification in power systems.Additionally, the work in [12] leveraged SVR to predict voltage deviations, demonstrating the versatility of SVM applications in enhancing the predictive capabilities and fault detection within power systems.
Thus, to address the above time-consuming problems, this paper introduces a method of coordinated reactive power planning in IPIES using SVR.The noteworthy contributions of this study are as follows: 1) A multi-objective optimization model for planning reactive power in an IPIES is proposed.Three indices are introduced, including the cost of reactive power source, static voltage margin, and voltage recovery ability.
2) The solution time can be significantly reduced for optimization by replacing it with SVR for traditional time-domain simulations.

Indices
This section introduces a variety of indices utilized to characterize power systems undergoing substantial integration of wind power.These indices encompass the investment cost index, static voltage stability margin index, and voltage recovery index.Within the scope of this discussion, these indices serve as pivotal metrics for evaluating and quantifying different aspects of power system performance in the context of increased wind power integration.The investment cost index addresses the economic considerations associated with the integration process, while the static voltage stability margin index gauges the resilience of the system under varying conditions.Additionally, the voltage recovery index provides insights into the system's ability to recover and stabilize voltage levels following perturbations.This comprehensive exploration of indices contributes to a nuanced understanding of the dynamic interactions within power systems augmented by significant wind power integration.generators within the system, thereby generating cost efficiencies.These reductions in generation expenses should be accounted for as a negative contribution to the reactive power source configuration cost within the overall investment cost.The methodology for calculating the total investment cost is delineated by Equation ( 1 where COST is the net expenditure, C inv is the overall expenditure, C save is the reduced operation cost by DSTATCOM, C ins is the fixed expenditure for installing a DSTATCOM, vector C cap s is the configuration cost for s th DSTATCOM, vector cap s is the configured reactive power capability of s th DSTATCOM, c unit is the unit reactive power cost for a DSTATCOM, C re-dis and C ori is the operation cost corresponding to system with DSTATCOMs, and operation cost corresponding to system without DSTATCOM, respectively, S is a set of DSTATCOMs.

Static voltage stability margin
The static voltage stability margin is defined as the distance between the current operating point of a system and its voltage stability boundary along a specific direction.This metric quantifies the system's resilience to external disturbances under a given operational state, providing a quantitative assessment of voltage stability.Specifically, it characterizes the relationship between the load at the system's steadystate operating point and its maximum load-carrying capacity, aiding in determining the maximum load the system can accommodate without compromising stability.For power systems, the static voltage stability margin holds significant implications.It furnishes crucial information for system operators and planners, facilitating an assessment of the voltage stability level under various operating conditions.Analyzing the static voltage stability margin of a power system allows the identification of the system's elasticity and robustness in the face of external disturbances.This, in turn, forms the foundation for devising effective power system planning and operational strategies.A comprehensive understanding and evaluation of this concept contribute to ensuring that power systems maintain reliable voltage levels across diverse operating conditions, thereby enhancing system robustness and sustainability.The index for the static voltage stability margin of the system is calculated as illustrated in Equation ( 2 where SVSM denotes the static voltage stability margin, λ denotes a load parameter, λ max denotes the maximum load power that the system can withstand, and λ op denotes the current load power that the system can handle.The λ max can be obtained via Equation ( where N is a set of buses, vector P b is a given increasing direction of active power, V is the bus voltage vector, P L and Q L are vectors for active and reactive power consumed by loads, P G and Q G are vectors for active and reactive power outputted by generators, Q S is the vector for output reactive power of DSTATCOMs, G ij is the conduction between bus i and bus j, B ij is the susceptance between bus i and bus j, θ ij is the angle equals to phase angle of bus i minus phase angle of bus j.It should be noted that Q Si ≡0 if there is no DSTATCOM for bus i.

Voltage recovery ability
The index for voltage recovery ability comprises two components: 1) The first part characterizes the rate of voltage amplitude recovery after fault clearance, defined as VRSI, depicted by the blue-shaded region in Figure 1.In cases of substantial disturbances such as short circuits, the voltage may sharply drop below the allowable minimum amplitude.By calculating the area S 1 between the voltage amplitude curve and the set voltage lower limit, the extent of voltage declines and the rate of voltage recovery after fault clearance can be effectively assessed.Furthermore, considering the presence of overshoot in the reactive power control system due to inherent dynamics, the area S 2 corresponding to the voltage overshoot beyond the upper limit must also be accounted for.
2) The second part is concerned with characterizing the time it takes for the system voltage to reach its final steady-state amplitude, defined as VRTI, represented by the red dashed line in Figure 1.Within this paper's context, the voltage load is considered to have reached a steady state when its amplitude falls within a ±2% range of the final steadystate value.The VRI, which is the index for describing voltage recovery represented by both VRSI and VRTI in a weighted manner, is formulated as shown in Equation ( 4), (1 ) where α denotes a weight parameter.Figure 1 illustrates voltage dynamics, with V indicating post-event voltage levels and V 0 as the preevent benchmark.The term t trip marks the fault clearance instant, leading to a recovery phase where t 1 is the point at which voltage climbs back to 90% of V0.Subsequently, t 2 and t 3 are moments when voltage equals to 110% of V 0 , during its ascent and descent, respectively.V s symbolizes the eventual voltage stabilization value, with t s pinpointing the time when voltage steadies within a narrow band of 0.98V s to 1.02V s .

Theory of SVR
The realm of Support Vector Regression (SVR) is intricately divided into two fundamental divisions, namely linear regression and nonlinear regression.A more detailed exploration reveals that linear regression within the SVR framework focuses on modeling and analyzing relationships between variables that display a linear pattern.On the other hand, nonlinear regression involves the investigation of connections between variables characterized by a more intricate, non-linear structure.Essentially, SVR encapsulates the dual characteristics of simplicity found in linear regression and the intricacies inherent in nonlinear regression.This amalgamation renders SVR a comprehensive and versatile framework for predictive modeling, well-suited for addressing a myriad of scenarios and complexities.

Linear regression
If the function that needs to be regressed is a linear function, it can be represented using the general equation, as shown in Equation ( 5), T ( ) where ω T and b are regression parameters that are to be determined.For a given set of k samples (x 1 , y 1 ), …, (x m, y m ), …, (x k , y k ), the loss function L ε (x m, y m ) is defined as shown in Equation ( 6). Figure 2 illustrates the conceptual diagram of the loss function.When the error between the actual value and the regression value is less than the given ε, it is considered accurate at that point.If the error is greater than the ε, it is included in the error function.In order to enhance the generalization ability of the regression function, relaxation factors ξ and ξ* can be introduced.In this case, the regression function transforms into a constrained optimization problem, as shown in Equation ( 7), ( ) ( ) , 0 where the purpose of the first term is to enhance the generalization ability of the regression function, while the purpose of the second term is to reduce regression error.C is the penalty factor.The larger the regression error is, the larger the penalty factor will be.
Based on Equation (7), the Lagrange equation is constructed as shown in Equation (8).
where, α m , α * m , λ m , and λ * m are all Lagrange multipliers.When the partial derivatives with respect to the parameters ω T , b, ξ, and ξ* are all equal to zero, Equation (8) achieves its minimum value, at which point the relationships depicted in Equations.( 9) to (12) exist.

(
)( The Lagrange multipliers α m and α * m can be obtained by solving Equation.(13), and in turn, we can obtain the parameter ω T as depicted in Equation.( 14).
( ) Finally, the variable b can be obtained as Equation.( 16) according to the KKT condition.

Non-linear regression
In the domain of non-linear regression problems, unraveling the intricacies of intricate relationships demands a paradigm shift.This entails transcending the constraints of the original vector by orchestrating a non-linear transformation, thereby navigating the data landscape in a higher-dimensional space.Through this transformative leap, the inherent complexities and non-linear patterns latent within the data emerge in a more analytically accessible form.Within this expanded higher-dimensional realm, the strategic application of linear regression techniques proves to be an effective means of modeling and comprehending the underlying relationships.The resulting outcome, denoted as z, encapsulates the transformed manifestation of vector x, symbolizing the enriched understanding gained through this augmented mapping process.The vector z satisfies the relationship depicted in Equation.( 17), ( , ) where the K(x m , x n ) is referred to as the kernel function, a non-negative definite symmetric function that satisfies the Mercer condition.Equation.(18) can be obtained by substituting z into Equation.(8).
After the nonlinear transformation, the quadratic programming problem is transformed into an optimization problem, as depicted in Equation.( 19).
The final regression function is shown in Equation.(20).

Training procedure
In this paper, the elements of vector x represent the capacity of DSTATCOM, and the predictive target is the index VRI.The process of training the SVR regression model using samples is illustrated in the flowchart in Figure 3.After obtaining the well-trained SVR model, the DSTATCOM capacity vectors are simply input into the model to calculate the VRI index, as shown in Equation.(21). ) ) , ( 0 , ) ( , ( , where F represent the comprehensive objective function, x is a set of system parameters, y is a set of state variables, z is a set of control variables, respectively.Additionally, E and D are constraints in terms of equality and inequality, respectively.

Constraints
The planning model incorporates a comprehensive set of constraints, encompassing both steady-state and dynamic aspects.Steady-state constraints encompass power balance constraints and system security constraints.These constraints are formally expressed as Equations.( 23)-( ] , , where P Gi and Q Gi are the active and reactive power generation in the IPIES, respectively, the generators include thermal power plants and combined heat and power plants, P fedi and Q fedi denote the power feed by grid, G and B are sets of generators and branches, L is a vector for the branch current.
In terms of dynamic constraints, the system encompasses a spectrum of considerations, including power angle difference constraints for synchronous generators, constraints for maximum recovery time, and constraints for load shedding.These constraints are mathematically formulated and can be expressed generically as presented in Equations.( 27 where the symbol vector  represents the power angle of synchronous generators.The parameter  max signifies the permissible maximum disparity in power angles.T r i denotes the moment when the voltage amplitude returns to the acceptable range, while T C represents the time point when the fault is rectified.T max corresponds to the maximum time threshold.Vectors P LS and P LS,max stand for the active load shedding magnitude and the allowable maximum active load shedding, respectively.The NSGA-II (Non-dominated Sorting Genetic Algorithm II) is employed as the optimization algorithm to address the planning model.The associated procedural diagram for solving the planning model is visually represented in Figure 4.In this iterative process, NSGA-II systematically explores and refines potential solutions, leveraging genetic principles and non-dominated sorting techniques to efficiently navigate the solution space.Expanding upon this, the algorithmic steps involve the initialization of a population of potential solutions, the evaluation of their fitness based on predefined objectives, and the application of genetic operators such as crossover and mutation to generate new candidate solutions.The non-dominated sorting process then categorizes solutions into Pareto fronts, emphasizing trade-offs between conflicting objectives.Through successive generations, NSGA-II converges towards a diverse set of well-balanced solutions, offering decision-makers a spectrum of trade-off options for the planning model.Figure 4 serves as a visual guide to elucidate the intricate steps and interactions within the NSGA-II-based optimization process.

Case study
The IEEE 39-bus test system, also known as the 10-machine New England System, is a widely recognized benchmark for power system dynamic stability and power flow studies.It simulates a smallscale, yet complex, electrical grid comprising 39 buses, 10 generators, and 19 loads, along with a variety of transformers and reactive power compensation devices.Designed to mimic a portion of the power network in the Northeastern U.S., the IEEE 39-bus system serves as a practical platform for analyzing power system behaviors, stability studies, fault simulations, and testing various control strategies under different operational conditions.This system enables researchers and engineers to assess grid performance in terms of voltage stability, frequency response, and dynamic behavior following system disturbances.The IEEE 39-bus test system is adopted for the case analysis, and the topology is illustrated in Figure 5.In Figure 5, the generators at bus 30, bus 32, and bus 34 are combined heat and power plants, and others are thermal plants.The expected contingencies encompass three-phase short-circuit faults occurring between line segments 2-3, 6-7, and 16-24, each with a duration of 0.1 seconds.The capacity of DSTATCOM is within a range of 17 to 100 Mvar.The fixed expenditure for installing a DSTATCOM is set at 1.5 million dollars per DSTATCOM, and the unit expenditure for configured reactive power capability is set at 50 thousand dollars per Mvar.
To enhance the simulation of the voltage amplitude curve and accurately compute the VRI, the DSTATCOM's configuration capacity is randomly produced within a predetermined range, ensuring a broad spectrum of scenarios.This generated capacity information is then integrated into sophisticated time-domain simulation tools to model the system's behavior accurately.The simulations yield VRI data, which is instrumental in training the SVR model, refining its predictive accuracy.The model's effectiveness is further validated through a meticulous comparison between its forecasts and actual system data, with discrepancies illustrated in a detailed error analysis graphically represented in Figure 6, showcasing the model's precision and areas for enhancement.where result time and result SVR are the values obtained via time domain simulation and SVR, respectively.Based on Figure 6, it's observed that a significant portion, precisely 584 data sets, exhibited errors beyond the 7% threshold, constituting around 6% of the total.A closer examination of this subset revealed that these data points generally represent larger values when contrasted with those exhibiting lesser errors.In the context of employing the NSGA-II algorithm for problem-solving, it becomes apparent that data points marred by larger errors are systematically excluded due to their elevated regression figures.This selective exclusion mechanism inherently diminishes the influence of regression inaccuracies within the SVR model, thereby enhancing the precision of the derived solutions.This refined approach improves the effectiveness of leveraging advanced algorithms to mitigate error impacts, ensuring more reliable outcomes in computational models.
The trained SVR model is brought into the solving process, which is shown in Figure 4.After solving it, a Pareto front is obtained, as shown in Figure 7.In Figure 7, two methods for calculating the VRI are applied.A comparison of the Pareto front obtained through NSGA-II optimization is illustrated.As seen in Figure 7, the Pareto front obtained using time-domain simulation and SVR model largely overlap.Decision-makers can choose the DSTATCOM planning solution from the Pareto front based on their specific priorities.The results in Figure 7 validate the feasibility of the method proposed in this paper.The comparison of solution time for using time domain simulation and SVR is listed in Table 1.  1, it can be observed that the time required for time-domain simulations is 26.35 hours, whereas the time spent using the SVR surrogate model is only one-tenth of that (0.27 hours).Therefore, employing the SVR surrogate model for calculating this multi-objective optimization model can result in significant time savings.

Conclusion
This paper introduces an innovative approach to coordinated reactive power planning for IIPIES, leveraging SVR as a key tool in the process.In this context, the DSTATCOM has been chosen as the primary system element for the power supply.The results obtained through extensive numerical analysis are particularly promising, showcasing a significant reduction in solution time.This enhancement not only streamlines the process of coordinated reactive power planning for IPIES but also presents a remarkable step forward in its efficiency.One noteworthy outcome of our study is the emergence of a diverse Pareto Front, encompassing a wide spectrum of potential solutions.This expanded range of choices empowers planners with greater flexibility and adaptability, enabling them to tailor their

Figure 4 .
Figure 4. Solution flowchart of the model.

Figure 6 .
Figure 6.Error analysis for SVR.In Figure6, a total of 10, 000 sets of comparative data are shown, and the error function is defined as shown in Equation.(30).

Figure 7 .
Figure 7. Pareto front for using time domain simulation and SVR.

Table 1 .
Comparison for solution time.