Journal of Modern Power Systems and Clean Energy

ISSN 2196-5625 CN 32-1884/TK

网刊加载中。。。

使用Chrome浏览器效果最佳,继续浏览,你可能不会看到最佳的展示效果,

确定继续浏览么?

复制成功,请在其他浏览器进行阅读

Large-signal Analysis and Controller Synthesis of Droop-based DC Power System with Saturation Constraints  PDF

  • Jinghan Zhao
  • Keting Wan
  • Yongpan Chen
  • Miao Yu (Member, IEEE)
College of Electrical Engineering, Zhejiang University, Hangzhou 310027, China

Updated:2025-05-21

DOI:10.35833/MPCE.2024.000164

  • Full Text
  • Figs & Tabs
  • References
  • Authors
  • About
CITE
OUTLINE

Abstract

In DC power systems dominated by power electronic devices, constant power loads (CPLs) and saturation components significantly impact large-signal stability. During the large-signal stability analysis process, the presence of multiple state variables and high-order system poses substantial challenges. To address this, considering the complete control dynamics, this paper proposes an equivalent single-machine (ESM) model of the droop-based DC power systems to reduce the complexity of the large-signal analysis. Building on the proposed ESM model, considering the dynamics of CPL and saturation constraints, a region of attraction (ROA) estimation algorithm based on sum of squares (SOS) programming is proposed, which significantly reduces the conservativeness compared with other existing methods. Furthermore, a control parameter optimization algorithm based on SOS programming is proposed with the aim of expanding the ROA. Furthermgre, with the aim of expanding the ROA, controller sythesis is conducted with proposed control parameter optimization algorithm based on SOS programming. Ultimately, simulation experiments validate the accuracy of the proposed ESM model and the proposed ROA estimation algorithm, as well as the effectiveness of the control parameter optimization algorithm.

I. Introduction

iL,i Inductance current of the ith converter
is Input current of ESM system
io Output current of ESM system
io,i Output current of the ith converter
is,i Input current of the ith converter
L Equivalent filter inductance
Li Filter inductance of the ith converter
n Number of converters
P Constant power load (CPL)
R Resistive load
vbus Bus voltage of load side
vo Output voltage of ESM system
vo,i Output voltage of the ith converter
vs Bus voltage of source side
vs,i Input voltage of the ith converter
B. Control Parameters
α Intermediate optimization parameter
d Duty cycle of ESM system
di Duty cycle of the ith converter
iLref Current reference of ESM system
iL,iref Current reference for the ith converter
k1-k6 Coefficients of monomials of optimization result
Kp, Ki Equivalent proportional and integral coefficients of ESM system
Kp,i, Ki,i Equivalent proportional and integral coefficients of the ith converter
Kcp Proportional coefficient of current loop control for ESM system
Kcp,i Proportional coefficient of current loop control for the ith converter
Kvi Integral coefficient of voltage loop control for ESM system
Kvi,i Integral coefficient of voltage loop control for the ith converter
Kvp Proportional coefficient of voltage loop control for ESM system
Kvp,i Proportional coefficient of voltage loop control for the ith converter
pi Current sharing coefficient of the ith converter
r Droop coefficient of ESM system
ri Droop coefficient of the ith converter
ulow, uup Lower bound and upper bound of saturation part
voref Voltage reference of ESM system
vo,iref Voltage reference of the ith converter
vref Rated output voltage reference

THE DC power system attracts widespread attention benefiting from its friendly integration of renewable energy sources (RESs) [

1]. In addition, it also possesses advantages such as high efficiency, high reliability, and simple control [2]. In engineering practice, due to the limitations imposed by the power rating of a single converter and economic considerations, the DC power system is mostly implemented by employing multiple power electronic converters in parallel through droop control [3]. This renders the DC power system essentially a power electronics-dominated system, exhibiting multiple state variables, high system order, strong nonlinear, low inertia, and weak damping characteristics. As a result, the system faces the risk of transient instability, and large-signal stability issues become prominent.

Firstly, because of the large number of state variables and high-order dynamics, the full-order model of droop-based DC power system is difficult to directly investigate, necessitating the application of order reduction techniques. In the droop-based DC power system, most studies only focus on the droop control loop, and represent the voltage droop control loop as a line impedance. Then, the whole system is further simplified as an equivalent single-machine (ESM) RLC circuit system [

4]. This simplification process neglects the voltage and current control loops, which is proven to be of significant importance during large disturbance transients [5]. Alternatively, from a frequency-domain perspective, we can identify the dominant dynamics components of the original system dynamics while neglecting other dynamics, thereby obtaining a corresponding reduced-order model [6]. However, the selection of dominant dynamics lacks specific guidance, and the neglect of other dynamics also introduces a certain degree of conservatism.

Furthermore, the current large-signal stability analysis in DC power systems predominantly relies on nonlinear system theories such as Lyapunov’s direct method, Takagi-Sugano (T-S) fuzzy method, and mixed potential theory (MPT) method, and tackles with nonlinear elements such as constant power loads (CPLs), saturation, and bilinear terms. However, when applying the Lyapunov’s direct method, the comprehensive guidelines for constructing Lyapunov functions are lacked. In traditional power system studies, it is common to directly select the energy function of the system as the Lyapunov function [

7], which overlooks the influence of controller and saturation on stability and is challenging to adapt to power electronics dominated power systems. Besides, inappropriately chosen Lyapunov functions can yield excessively conservative analysis results.

Based on fuzzy rules, the nonlinear state-space matrix of complex DC power systems can be approximated by a weighted sum of a finite number of linearized steady-state subsystems, namely T-S fuzzy modeling. This approximation facilitates the application of the Lyapunov’s indirect method for analysis and the construction of linear matrix inequality (LMI) problems to solve for a quadratic Lyapunov function [

8]. It is able to handle the system that includes saturation, rational, and bilinear components. Nevertheless, T-S fuzzy method and LMI method tend to produce conservative outcomes from both modeling and stability analysis steps. For high-order systems accompanied with intricate nonlinearities, the exponential growth of fuzzy rules poses challenges in T-S fuzzy modeling, problem solving, and stability analysis.

MPT method constructs a scalar function resembling an energy function based on the circuit structure. When applied to DC power systems, the MPT method requires the equivalence of power electronic devices to controlled sources [

9]. The key challenge lies in preserving control strategies and system nonlinear characteristics to the greatest extent possible during the equivalence process, meanwhile without violating the completeness of topology [10]. Obviously, it fails to account for the complete control loop and cannot handle the saturation parts well.

Table I shows the feature comparison of region of attraction (ROA) estimation methods.

TABLE I  Feature Comparison of ROA Estimation Methods
MethodApplicable systemCan nonlinearity be handled?Can control loop be fully considered?Preprocessing costComputational burdenConservatism
T-S fuzzy [8] AC and DC power systems Fractional, saturated, trigonometric, polynomial High Medium High
MPT [9], [10] Circuit system with complete topology - × Medium Low High
Lyapunov’s direct [7] Traditional AC power system - × Medium Low High
Monte Carlo All - Low Very high Very low
This paper AC and DC power systems Fractional, saturated, trigonometric, polynomial Low High Low

The methods mentioned above for constructing Lyapunov functions can only be applied to the systems with specific structures. The nonlinearities they can handle are limited. The applications of these methods do require a certain level of researcher knowledgement, resulting in a certain degree of preprocessing cost. The obtained Lyapunov functions are fixed and cannot fully incorporate the system control elements into the construction process. Apart from these methods, although the Monte Carlo method is accurate, it imposes a significant computational burden. Without relying on the researcher knowledgement, directly searching for the candidate of the Lyapunov function that satisfies the non-negativity conditions in the Lyapunov stability theorem is an NP-hard problem [

11]. Notice that sum of squares (SOS) polynomial is semi-positive definite everywhere. Inspired by this property, [12] replaces the non-negativity conditions in the Lyapunov stability theorem with SOS polynomials, which makes the solving process computationally tractable through SOS programming. Consequently, this approach is applied to AC power systems for large-signal stability analysis such as ROA estimation [13] and low-voltage ride-through analysis [14]. In the case of DC power systems, the ROA of an islanded DC microgrid is estimated based on SOS polynomial [15]. The method presented in [15] can be applied to common power systems. It can handle various types of nonlinearities to achieve the results with lower conservativeness at an acceptable computational burden. However, the neglect of inner control components and the approximation of non-polynomial parts in these previous studies may introduce conservatism into the ROA estimation.

One of the common causes of large-signal instability is the improper design of control parameters [

16]. Reasonable and accurate stability analysis can guide the controller design to enlarge the stability margin. In the case of droop-based DC power system, droop control generates the voltage reference, which is fed into voltage-current dual-loop control to generate the duty cycle signal. In practice, voltage-current dual-loop control parameter tuning is conducted in frequency domain, aiming to achieve bandwidth separation between the outer loop and inner loop, thus ensuring overall system dynamic response [8]. Existing control parameter optimization studies are based on eigenvalues analysis, e.g., [17], which aim to optimize control parameters by ensuring that the dominant eigenvalues of the system are in the left half of the s-plane and are far away from the imaginary axis as possible. Essentially, this is a small-signal stability analysis and is not applicable to the large-signal stability analysis when facing scenarios with large disturbances. To further enhance the system stability, the concept of virtual impedance is introduced. Studies in [18] propose a parallel/series virtual impedance control strategy to reshape the input/output impedance of the converter to dampen the instability of the system. Subsequent studies explore the use of nonlinear control methods as alternatives to proportional-integral (PI) control to enhance the large-signal stability of the system. A high-gain observer is designed in [19] to estimate the electrical coupling among subsystems, which is offset by composite controller based on feedback linearization. Then, the global large-signal stability of the system is realized through stabilizing the internal states of subsystems. Similar ideas as “estimation-compensation” can also be found in [20]. Compared with the PI controller, the above-mentioned improved controllers require complex design procedure, heavy computational burden, and high implementation costs. Thus, the generalization of the above-mentioned methods proves to be challenging. In contrast, the design of a PI controller with a focus on large-signal stability is of great practical significance, which is one of the focuses of this study.

The main contributions of this paper are as follows.

1) For large-signal stability analysis, considering the complete circuit-control topology, an ESM model of the droop-based DC power systems is proposed. The proposed ESM model maintains consistent dynamic characteristics with the original system, but significantly reduces the system dimension and the complexity of the large-signal stability analysis.

2) Considering the whole control components, the entire dynamics of CPL, and saturation constraints, a ROA estimation algorithm based on SOS programming is proposed. Compared with other existing methods, the proposed ROA estimation algorithm yields estimation results with lower conservatism.

3) With the objective of expanding the ROA, a control parameter optimization algorithm based on SOS programming is proposed. The droop-based DC system with parameter optimization can withstand greater perturbations, leading to further improvement in large-signal stability.

The structure of this paper is organized as follows. Section II introduces the multi-converter equivalent modeling. In Section III, the ROA estimation algorithm based on SOS programming is presented. Section IV shows the control parameter optimization algorithm based on SOS programming. Section V discusses case studies and numerical validation. Finally, Section VI provides the conclusions of this study.

II. Multi-converter Equivalent Modeling

In this section, considering both the droop control and voltage-current dual-loop control, an ESM model of the droop-based DC power system is developed to maintain consistent dynamic characteristics with the system. Thus, the large-signal stability analysis can be conducted based on the proposed ESM model with lower computational complexity.

For a typical droop-based DC power system, it consists of n DC-DC converters, which are placed in an input-parallel output-parallel (IPOP) form. All the loads can be regarded as a ZIP load, where Z denotes the constant impedance; I denotes the constant current; and P denotes the constant power. The control structure of the droop-based DC power system and the proposed ESM model is shown in Fig. 1.

Fig. 1  Control structure. (a) Droop-based DC power system. (b) Proposed ESM model.

This study makes reasonable assumptions for the droop-based DC power system. Furthermore, we reasonably simplifies the findings in [

21] to derive a more concise ESM model. The equivalent conversion process, including control elements, is briefly introduced as follows. Before deriving the more concise ESM model, the necessary assumptions are presented as:

Assumption 1: the line resistance is far less than the virtual resistance introduced by the droop control. Thus, the total output current of the droop-based DC power system is shared by each converter in proportion to their reciprocals of the droop coefficients [

22].

Assumption 2: the value of Ci is designed under the guidelines in [

23] according to the rated output current Io,i, output voltage ripple Δv, and switch frequency fi.

Ci=Io,ifiΔv (1)

Assumption 3: the control parameters of voltage-current dual-loop control are designed as the widely-used guideline in [

24].

Kcp,iLi=Kcp,i+1Li+1 (2)

For the IPOP form, the equivalent filter inductance and capacitance can be conducted as (3) and (4), respectively.

L=1i=1n1Li (3)
C=i=1nCi (4)

It is worth mentioning that the variables of the proposed ESM model correspond to those of the EMS system. Considering the control structure of current control loop of the droop-based DC power system and the proposed ESM model, combined with Assumptions 1 and 2, and taking into account the current distribution effect achieved by the droop control, the following equation holds:

Kcp(iLref-iL)L=i=1nKcp,ipi(iLref-iL)Li (5)

Notice that iL,i=piiL can be got approximately and i=1npi=1 is satisfied apparently.

The proportional coefficient of the current control loop in the proposed ESM model can be calculated as:

Kcp=Li=1nKcp,ipiLi (6)

Considering the voltage control loop of the ith converter, the inductance current reference can be generated as:

iL,iref=Kp,i(vo,iref-vo,i)+Ki,i(vo,iref-vo,i)dt (7)

The current reference of the proposed ESM model can be generated as:

iLref=Kp(voref-vo)+Ki(voref-vo)dt (8)

By replacing the inductance current reference in (5) with (7) and (8), and then eliminating the inductance current terms by substituting (6), it can be found that there is a one-to-one correspondence between the proportional part and integral part of the droop-based DC power system and the proposed ESM model.

KcpKp(voref-vo)L=i=1nKcp,iKp,i(vo,iref-vo,i)Li (9)
KcpKi(voref-vo)dtL=i=1nKcp,iKi,i(vo,iref-vo,i)dtLi (10)

For the IPOP form, the input voltage, the output voltage of the proposed ESM model, and the output voltage of the ith converter are identical, respectively. Thus, the equivalent control parameters of the voltage control loop for the proposed ESM model can be derived as:

Kp=i=1nKcp,iKp,iLii=1nKcp,ipiLi (11)
Ki=i=1nKcp,iKi,iLii=1nKcp,ipiLi (12)

Since the design of control parameters follows the instruction mentioned in Assumption 2, the relationship between proportional coefficients of different converters, which is shown in Assumption 3 as (2), is used to simplify (11) and (12).

Kp=i=1nKp,i (13)
Ki=i=1nKi,i (14)

Finally, when considering the droop control of the ith converter, under Assumptions 1, the output power of each converter is inversely proportional to the droop coefficient ri.

ri=rpi (15)

There are n DC-DC converters working together to supply the various loads, which can be regarded as one common ZIP load [

25].

io=i=0nio,i=voR+IC+Pvo (16)

Hereinbefore, the proposed ESM model of the droop-based DC power system with n DC-DC converters is given as follows.

For the circuit part, voltage source is thought to be ideal and disturbance can be neglected, which means that vs=Vs.

LdiLdt=dVs-voCdvodt=iL-voR-Ic-Pvo (17)

For the control part, combining (6) with (13)-(15), we can obtain that:

d=KcpKpVref-rvoR+IC+Pvo-vo+KiVref-rvoR+IC+Pvo-vodt-iL (18)

Specially, the integral term in (18) is listed individually as a new state variable ξ with dynamics as:

dξdt=Vref-rvoR+IC+Pvo-vo (19)

Furthermore, the state variables in (16)-(19) are described in a form of x=X+x^, where x is the real value; X is the constant term in the stable state; and x^ is the disturbance term. Rewrite (18) by using the new state variable ξ and substitute it into (17). Then, eliminate the identity terms on both sides of the equation. Finally, the large-signal model of the ESM system is shown in (20).

di^Ldt=-v^o+Vssat(d^)dv^odt=PCVo-PC(Vo+v^o)+i^LC-v^oCRdζ^dt=rPVo-rPVo+v^o-1+rRv^o (20)

In the practical modulation of the duty cycle d, there exists a saturation, which limits d to be in the range of [0,1], thereby bringing in saturated nonlinearity. After the coordinate transformation, the control input can be given as:

d^=-KcpKp1+rRv^o-rPVo-rPVo+v^o+KcpKiξ^-Kcpi^L (21)

The saturation in the large-signal model is:

sat(d^)=-VoVs            d^<-VoVsd^                    -VoVsd^1-VoVs1-VoVs        d^>1-VoVs (22)

It should be pointed out that the value of Vo is derived by solving the expression of (17) in the stable state.

Vo=(Vref-rIC)+(rIC-Vref)2-4rP1+rR21+rR (23)

Based on (20), we can observe that the CPL and duty cycle saturation mainly introduce the nonlinear term into the droop-based DC power system. This term is discussed in Section III.

III. Proposed ROA Estimation Algorithm Based on SOS Programming

The level set of the Lyapunov function serves as an inner estimation of the ROA for the corresponding ESM system [

26]. To obtain the largest estimated ROA of the ESM system, a proper Lyapunov function must be found, and its upper level set must be optimized to satisfy the constraints of the Lyapunov theorem. This section focuses on ESM systems in the form of (20) and presents the proposed ROA estimation algorithm based on SOS programming [27].

Referring to (20), we can consider an ESM system in a general form, which can be given as:

x˙=f0(x)+g(x)·sat (u0(x)) (24)

where x is a vector composed of the all elements x, xRn; u0(x)Rn; g(x) is a multivariable polynomial function; and f0(x) and u0(x) are the polynomial rational functions in n variables x with coefficients in Rn.

Specially, the system in (20) has the same polynomial dominator, which is indicated by h(x). The numerator part of f0(x) and u0(x) are polynomial and expressed as f(x) and u(x), respectively, where u(x)=[u1(x),0,0]. Then, (24) can be transformed as:

x˙=f(x)h(x)+g(x)·satu(x)h(x) (25)
satu(x)h(x)=ulow          u1(x)h(x)>|ulow|u(x)h(x)        |ulow|u1(x)h(x)|uup|uup           u1(x)h(x)>|uup| (26)

In this paper, we choose the expanding interior algorithm to conduct the proposed ROA estimation algorithm based on SOS programming. The detailed introduction can be found in [

28]. In this section, it is improved to be applicable to ESM systems in the form of (25).

The expanding interior algorithm contains a series of set inclusion relationships. Instead of optimizing both the level set and Lyapunov function, the level set is fixed at 1 in this algorithm. We define Ω as the estimated ROA by using an unknown Lyapunov function V(x)Rn with the fixed level set to be 1, where Rn denotes the set of polynomial functions consisting of n variables.

Ω={xRn|V(x)1} (27)

The unknown V(x) is updated during each iteration until the optimization variable converges.

The interior of ROA to be estimated is described by a polynomial shape function p(x) and level set β, which is defined as:

P={xRn|p(x)β} (28)

The proposed ROA estimation algorithm can be summarized as the solution of the following optimization problem.

max βs.t.  PΩ       Ω{xRn|V˙<0} (29)

where V˙=Vxf(x)h(x)+g(x)u(x)h(x).

However, the expanding interior algorithm is designed for pure polynomial system, which is not suitable for the system we studied with saturation and rational parts. Besides, the shape function is fixed, resulting in a single optimization direction [

29]. Thus, an extended form of conditions in (29) is given to cope with the special system form in this paper. Then, a procedure for updating shape function is added to reduce the conservatism problem introduced by the fixed shape function.

Firstly, we carry out the form transformation. Lyapunov stability theorem requires that Lyapunov function is positive definite in its domain. For simplicity, we choose the V(x), which is positive definite for all xRn except x=0. For saturation part, the upper and lower bounds of u0(x) are separately presented as the inclusion of sets. Then, a feasible V(x) for the proposed ROA estimation algorithm of the system in (25) requires that the following constraints hold:

V(x)=0       xRn\0V(0)=0 (30a)
PΩ (30b)
Ω\0{xRn|V˙<0} (30c)
Ω\0xRnu1(x)h(x)|uup| (30d)
Ω\0xRnu1(x)h(x)|ulow| (30e)

However, the system in (25) is not a polynomial system and cannot be directly addressed by SOS programming. Thus, for the rational part, we multiply h2(x) on both sides of (30c)-(30e) to perform an equivalent transformation on the original constraints, as shown in (31).

D\0{xRn|V˙h2(x)<0}D\0{xRn|u1(x)h(x)-|uup|h2(x)0}D\0{xRn|-u1(x)h(x)+|ulow|h2(x)0} (31)

According to the Positivstellensatz theorem [

30], with a feasible set of auxiliary SOS polynomials si, the set constraints (30a)-(30b) and (31) are converted to constraints (32).

V-l1Σn-s1(β-p)-(V-1)Σn-s2(1-V)-s3V˙h2-l2Σn-s4(1-V)-(-|uup|h2+u1h)Σn-s5(1-V)-(|ulow|h2-u1h)Σn (32)

During this conversion, we drop the term x in the expression for simplicity. The constraints of the set D\{0} are not semi-algebraic and are replaced with l10 and l20, where l1,l2Σn.

The SOS optimization problem for the largest estimated ROA is expressed as:

maxVRn,s1,s2,s3,s4,s5Σn βs.t.  (32) (33)

Then, we give the process on design of the improved expanding interior algorithm. Notice that in (32), the constraints are bilinear in si and V, which cannot be solved directly. Drawn on the general idea of coordinate-wise method [

31], the SOS optimization problem is separated into two main subproblems. The SOS programming in the first subproblem is modified for ease of solving, as shown in (34). A detailed explanation of this modification will be given later.

-s2(α-V)-s3V˙h2-l2Σn-s4(α-V)-(-|uup|h2+u1h)Σn-s5(α-V)-(|ulow|h2-u1h)Σn (34)

Denote z as the number of iterations. With the number of iterations indicated in the superscript, the flowchart of the improved expanding interior algorithm is shown in Fig. 2.

Fig. 2  Flow chart of improved expanding interior algorithm.

The main steps of the improved expanding interior algorithm are given as follows.

Step 1:   initialize Lyapunov function V(0) from the linearized system through using Lyapunov’s indirect method and initialize shape function p(0)=V(0). Given the auxiliary polynomials l1 and l2(z), set the degrees of V and si to be solved, and set the tolerance εs. Start the improved expanding interior algorithm with β(0)=0 at z=0.

Step 2:   fix V=V(z) and p=p(z). Do line search on α(z) for feasible s2(z)-s5(z).

Step 3:   record α(z)=α and update auxiliary SOS polynomials s2(z)=s2, s2(z)=s2, s4(z)=s4, s5(z)=s5 in (35).

maxs2,s3,s4,s5Σn αs.t.  34 (35)

Step 4:   update l2(z+1)=l2(z)/α(z), and rescale s4(z), s5(z) with α(z)s4(z) and α(z)s5(z).

Step 5:   update pz with V(z)/α(z).

Step 6:   fix p=pz, s2(z)=s2, s2(z)=s2, s4(z)=s4, s5(z)=s5, and l2=l2(z) and do line search on β for feasible s1 and V. In addition, record β(z+1)=β and update Lyapunov function V(z+1)=V.

maxVRn,s1Σn βs.t.  32 (36)

Step 7:   if |β(z+1)-β(z)|εs, then end the iteration; otherwise, z=z+1 and go to Step 2.

Each step has a specific task as follows. Step 1 is the initialization for the whole improved expanding interior algorithm. Step 2 generates partial auxiliary SOS polynomials for next steps. Since the fixed level set 1 may results in the same solution, an intermediate optimization parameter α is introduced. Step 3 rescales the polynomials from Step 2 to accommodate with the constraints in next steps with level set 1. Step 4 updates shape function with the current Lyapunov function to avoid conservatism problems caused by improper selection of the shape function [

29]. Besides, it also ensures the estimated ROA continues to expand. Step 5 searches for the Lyapunov function candidate based on the results from Steps 2-4. The optimization subproblems in Step 2 and Step 5 are both quasi-convex problems. In this paper, bisection method is used for the line search. It makes sure the optimum can be computed [32], thereby guaranteeing the effectiveness and convergence of the proposed ROA estimation algorithm.

IV. Control Parameter Optimization Algorithm Based on SOS Programming

For the studied system in (20), with the introduction of new state variable ξ, control input u0 has a form as polynomial rational function. The detailed expression for u1 corresponding to the system in (20) in the form of (25) is given as:

u1=-KcpVoi^L+rKcpKpPVo-KcpKpVo-rKcpKpVoRv^o+KcpKiVoξ^-Kcpi^Lv^o+-KcpKp-rKcpKpRv^o2+KcpKiv^oξ^ (37)

In (37), i^L, v^o, ξ^, i^Lv^o, v^o2, and v^oξ^ are the specific monomials. This form of polynomial can be searched during the optimization and enlarge the set of points that are attracted to the equilibrium. By one-to-one correspondence between the polynomial coefficients obtained after optimization and the original expression, the optimized control parameters can be obtained.

For control parameter optimization, note that the form of the original optimization problem remains unchanged, except that u1 is added as a decision variable.

maxV,u1Rn,s1,s2,s3,s4,s5Σn βs.t.  (32) (38)

A tri-linear term exists in s3V˙h2(x) in (32), which expands as s3Vx(fh+guh). To deal with the tri-linear term, a control parameter optimization algorithm is designed. The main structure of the control parameter optimization algorithm is an internal controller synthesis iteration nested in the external Lyapunov function synthesis iteration. This makes the unknown polynomials all appear affine linear in the constraints, so the SOS programming can be done to find the feasible solution. Denote j and k as the number of external iterations and the internal iterations, respectively. With the number of iterations indicated in the superscript, the main steps of the control parameter optimization algorithm are given as follows.

Step 1:   firstly, initialize Lyapunov function V(0) from linearized system by using Lyapunov indirect method. Then, initialize shape function p(0)=V(0). Then, initialize u1(0) with control parameters from traditional method such as Bode plot or control bandwidth design, making sure the close-loop stability of the linearized system. Then, specify the auxiliary polynomials l1 and l2j. For V and si to be solved, set the degrees. For u to be solved, set its degree to be 2. In addition, set the tolerance εC and εL for the controller synthesis iteration and Lyapunov function synthesis iteration, respectively. Finally, start the control parameter optimization algorithm with α(0)=0 and β(0)=0 at j=0 and k=0.

The controller synthesis iteration is given as:

Step 2-1:   fix V=V(j), p=p(j), u1=u1(k), and l2=l2(j) and do line on α(k+1) for feasible s2(k)-s5(k) in (39). In addition, update the auxiliary SOS polynomial s3(k)=s3.

maxs2,s3,s4,s5Σn αs.t.  (34) (39)

Step 2-2:   fix V=V(j), p=p(j), s3=s3(k), and l2=l2(j) and do line search on α(k+1) for u(k+1) in (40). In addition, record α(k+1)=α and update u1(k+1)=u1.

max u1Rn,s2,s4,s5Σn αs.t.  (34) (40)

Step 2-3:   if |α(k+1)-α(k)|εC, go to the Lyapunov function synthesis procedure; otherwise, set k=k+1 and go to Step 2-1 to continue the controller synthesis iteration. After the controller synthesis iteration converges within the tolerance, the Lyapunov function synthesis iteration is carried out.

Step 3-1:   fix V=V(j), p=p(j), u1=u1(k+1), and l2=l2(j) and do line search on α for s2-s5 in (41). Update s2(j+1)=s2, s3(j+1)=s3, s4(j)=s4, and s5(j)=s5. Then, set rs(j)=α.

maxs2,s3,s4,s5Σn αs.t.  (34) (41)

Step 3-2:   rescale the auxiliary SOS polynomials: s4(j+1)=rs(j)s4(j), s5(j+1)=rs(j)s5(j), and l2(j+1)=l2(j)/rs(j).

Step 3-3:   fix p=p(j), u1=u1(k+1), s2=s2(j+1) and l2=l2(j+1) and do line search on β for V(j+1) in (42). Then, record β(j+1)=β and update V(j+1)=V.

maxVRn,s1Σn βs.t.  (32) (42)

Step 3-4:   update p(j+1)=V(j+1);

Step 3-5:   if the variation of β converges and satisfies the tolerance setting, then end the iteration; otherwise, set j=j+1 and go to controller synthesis iteration again.

For the controller synthesis iteration, α is introduced. Step 1 generates the auxiliary SOS polynomials for the next step. Step 2-2 searches for a feasible u to enlarge the level set of fixed Lyapunov function. When the iteration converges and the level set is expanded to εC, a new Lyapunov function is needed. Then, go to the Lyapunov function synthesis iteration. In this iteration, u1 is fixed. The auxiliary polynomials is re-evaluated in Step 3-1. Then, it is rescaled in Step 3-2 to for next steps with fixed level set 1. In Step 3-3, a new Lyapunov function is constructed for estimating ROA with optimized u1. The flowchart of control parameter optimization algorithm is shown in Fig. 3.

Fig. 3  Flowchart of control parameter optimization algorithm.

Similarly, the whole optimization problem consists of several quasi-convex optimization subproblems. Its optimum can be solved by the control parameter optimization algorithm proposed in this paper effectively. The validity and convergence proof of this algorithm are the same as the improved expanding interior algorithm.

With the composition of u1 specified in advance, the result of optimized u1,opt after running the control parameter optimization algorithm is guaranteed in a form as:

u1,opt=k1i^L+k2v^o+k3i^Lv^o+k4v^o2+k5ξ^+k6v^oξ^ (43)

By checking (43) against (37), the optimized control parameters (adding opt in the subscript to distinguish new ones from the original variables) are calculated reversely as:

Kcp,opt=-k3Kp,opt=-(k2-k4Vo)Vork3PKi,opt=-k6k3 (44)

Up to this point, the entire control parameter optimization process is completed.

V. Case Studies and Numerical Validation

In this section, a droop-based DC power system with 2 paralleled DC-DC converters is considered as a case to validate the effectiveness of the proposed ESM model, the proposed ROA estimation algorithm, and controller parameter optimization algorithm.

The parameters of the original system are listed in Table II. All the parameters are designed referring to the guidelines for DC-DC converters mentioned in [

8] and [23]. Thus, the assumptions of the proposed ESM model are satisfied automatically. Parameters of the ESM system are derived and listed in Table III.

TABLE II  Parameters of Origianal System
Parts of studied systemParameterValue
Parts of studied system Vs 800 V
Vref 400 V
Load R 40 Ω
P 5 kW
IC 10 A
DC-DC converter 1 L1 6 mH
C1 0.45 mF
r1 0.0533 V/A
Kvp,1 0.1414 A/V
Kvi,1 4.4414 A/(V·s)
Kcp,1 0.0079 V/A
DC-DC converter 2 L2 2 mH
C2 0.15 mF
r2 0.1600 V/A
Kvp,2 0.0471 A/V
Kvi,2 1.4805 A/(V·s)
Kcp,2 0.0236 V/A
TABLE III  Parameters of ESM System
ParameterValue
Vs 800 V
Vref 400 V
L 1.5 mH
C 0.6 mF
r 0.0400 V/A
Kvp 0.1885 A/V
Kvi 5.9218 A/(V·s)
Kcp 0.0118 V/A

The computation is carried out on an IntelR CoreTM i5-1135G7, 2.40 GHz, 4-cores, 16 GB RAM, MATLAB/Simulink 2019a environment.

Figure 4 shows the comparison of transient process waveforms. Regardless of the transient process, steady-state waveform or instability waveform, the waveforms of the original system and the ESM system exhibit highly consistent dynamic characteristics. This indicates that the original system and its ESM system have the same transient stability. This result verifies the effectiveness of the proposed ESM model and indicates that subsequent large-signal stability analysis can be carried out based on the ESM system.

Fig. 4  Comparison of transient process waveforms. (a) 2 kW resistive load increase. (b) 2.5 kW CPL increase. (c) 20 V voltage sag. (d) 100 V voltage sag.

Secondly, to demonstrate the advancement of the proposed ESM model, this paper uses the ROA obtained by the T-S fuzzy method for comparison. This is because the T-S fuzzy method can effectively handle the rational components and saturation parts in the original system, while other methods or models are not good at these nonlinear parts. At the same time, to verify the accuracy of the proposed ESM model, the real ROA is obtained through Monte Carlo method. Comparison of ROAs obtained by various methods are shown in Fig. 5.

Fig. 5  Comparison of ROAs obtained by various methods. (a) 3D view. (b) i^L-v^o plane at ξ^=0. (c) i^L-v^o plane at ξ^=2. (d) i^L-v^o plane at ξ^=-2.

As shown in Fig. 5(a), under the premise of satisfying the upper and lower bounds, the results obtained by the proposed ESM model completely encompass the results obtained by the T-S fuzzy method. Since the true ROA is obtained through the Monte Carlo method, it is difficult to present it intuitively in a 3D view. Therefore, the estimated results and the true attraction domain are projected onto different i^L-v^o planes for direct comparison. Obviously, the proposed ESM model obtains a less conservative ROA estimation result. It is close to the real ROA within a certain range and, simultaneously, satisfies the saturation constraints.

Finally, the control parameters are tuned by the control parameter optimization algorithm, as shown in Table IV.

TABLE IV  Optimized control parameters
ModuleParameterValue
ESM system Kvp 0.5262 A/V
Kvi 12.0782 A/(V·s)
Kcp 0.0054 V/A
DC-DC converters 1 and 2 Kvp,1 0.3947 A/V
Kvp,2 0.1316 A/V
Kvi,1 9.0587 A/(V·s)
Kvi,2 3.0196 A/(V·s)
Kcp,1 0.0014 V/A
Kcp,2 0.0041 V/A

The T-S fuzzy method only provides a qualitative assessment of the impact of a single parameter on the large-signal stability of the system, which is not suitable for comprehensive controller synthesis. Thus, we compare the conventional design algorithm with the control parameter optimization algorithm in this paper.

The comparison of the ROA with and without optimization is shown in Fig. 6. It is worth noting that with and without optimization in Fig. 6 mean with and without the control parameter optimization algorithm proposed in this paper, respectively. The ROA of the ESM system has certain expansions in both the i^L-axis and v^o-axis directions and completely encompasses the ROA without optimization. Furthermore, as shown in Fig. 6(b)-(d), the true ROA calculated by the Monte Carlo method also verifies this conclusion, indicating that the control parameter optimization algorithm effectively improves the system ROA, thereby enhancing the large-signal stability of the system.

Fig. 6  ROA comparison with and without optimization. (a) 3D view. (b) i^L-v^o plane at ξ^=0. (c) i^L-v^o plane at ξ^=2. (d) i^L-v^o plane at ξ^=-2.

During the implementation of the control parameter optimization algorithm, SOS programming can only solve problems of moderate scale [

33]. To enhance the efficiency of solution-seeking process, users should simplify or decompose the research object reasonably. Alternatively, users can employ more efficient SOS algorithms such as the dominant sum of squares (DSOSs) and the scaled diagonally dominant sum of squares (SDSOSs) [34].

To further demonstrate the practical effect of control parameter optimization algorithm, a comparison experiment is designed using MATLAB/Simulink. One set of the system control parameters is designed using PI control, as shown in Table II, while the other set uses the optimized control parameters obtained by the control parameter optimization algorithm, as shown in Table IV. The results are depicted in Fig. 7.

Fig. 7  Comparison of voltage waveforms with and without optimization after 25 V voltage sag.

At t=0.05 s, a 25 V voltage sag lasting for 10 ms occurs. During the transient process, both voltage waveforms with and without optimization begin to oscillate. After the voltage sag occurs, the waveform with optimization returns to a steady state. In contrast, the amplitude of the waveform without optimization increases gradually, becomes unstable, and then enters a constant amplitude oscillation state. This comparison result indicates that the control parameter optimization algorithm does effectively expand the system ROA and enhance the large-signal stability of the system when facing external disturbances.

VI. Conclusion

A novel SOS programming for analyzing and enhancing the large-signal stability of the droop-based DC power system with saturation constraints is proposed in this paper. The proposed ESM model including the inner-control loops is developed. Considering the nonlinear dynamics of CPL and saturation constraints, an ROA estimation algorithm is designed. The proposed ROA estimation algorithm aims to find an appropriate Lyapunov function and maximize its level set to estimate the ROA based on SOS programming. The results demonstrate that the proposed ROA estimation algorithm provides a more accurate estimation of the ROA compared with other existing methods. Furthermore, a control parameter optimization algorithm is proposed based on SOS programming, which expands the ROA, thereby enhancing the large-signal stability of the droop-based DC power system. The experiment results show the effectiveness of the control parameter optimization algorithm.

Given that the large-signal stability analysis and control parameter optimization algorithm proposed in this paper retain the commonly used voltage-current dual-loop control structure of droop-based DC power systems, they are practical and thus have certain real-world significance and engineering application value. Furthermore, the proposed ESM model can be widely applied to various types of converters such as boost converters and single-phase AC-DC converters. It can also be applied to different converter topologies such as cascaded and parallel configurations. The SOS programming can also be further combined with nonlinear theory to develop the prescribed performance or robust nonlinear controllers, further enhancing the disturbance rejection capability of the system.

References

1

K. Jithin, P. P. Haridev, N. Mayadevi et al., “A review on challenges in DC microgrid planning and implementation,” Journal of Modern Power Systems and Clean Energy, vol. 11, no. 5, pp. 1375-1395, Sept. 2023. [Baidu Scholar] 

2

J. M. Guerrero, J. C. Vasquez, J. Matas et al., “Hierarchical control of droop-controlled AC and DC microgrids: a general approach toward standardization,” IEEE Transactions on Industrial Electronics, vol. 58, no. 1, pp. 158-172, Jan. 2011. [Baidu Scholar] 

3

Y. Xia, W. Wei, Y. Peng et al., “Decentralized coordination control for parallel bidirectional power converters in a grid-connected DC microgrid,” IEEE Transactions on Smart Grid, vol. 9, no. 6, pp. 6850-6861, Nov. 2018. [Baidu Scholar] 

4

H. Tu, H. Yu, and S. Lukic, “Impact of virtual inertia on DC grid stability with constant power loads,” IEEE Transactions on Power Electronics, vol. 38, no. 5, pp. 5693-5699, May 2023. [Baidu Scholar] 

5

F. Göthner, J. Roldán-Pérez, R. E. Torres-Olguin et al., “Reduced-order model of distributed generators with internal loops and virtual impedance,” IEEE Transactions on Smart Grid, vol. 13, no. 1, pp. 119-128, Jan. 2022. [Baidu Scholar] 

6

Q. Fu, W. Du, H. Wang et al., “DC voltage oscillation stability analysis of DC-voltage-droop-controlled multiterminal DC distribution system using reduced-order modal calculation,” IEEE Transactions on Smart Grid, vol. 13, no. 6, pp. 4327-4339, Nov. 2022. [Baidu Scholar] 

7

P. Kundur, N. J. Balu, and M. G. Lauby, Power System Stability and Control. New York: Wiley IEEE Press, 2003. [Baidu Scholar] 

8

S. Bacha, I. Munteanu, and A. I. Bratcu, Power Electronic Converters Modeling and Control: with Case Studies. New York: Springer Science & Business Media Press, 2013. [Baidu Scholar] 

9

Z. Wu, H. Han, Z. Liu et al., “A novel method for estimating the region of attraction for DC microgrids via brayton-moser’s mixed potential theory,” IEEE Transactions on Smart Grid, vol. 14, no. 4, pp. 3313-3316, Jul. 2023. [Baidu Scholar] 

10

Z. Liu, X. Ge, M. Su et al., “Complete large-signal stability analysis of DC distribution network via brayton-moser’s mixed potential theory,” IEEE Transactions on Smart Grid, vol. 14, no. 2, pp. 866-877, Mar. 2023. [Baidu Scholar] 

11

K. G. Murty and S. N. Kabadi, “Some NP-complete problems in quadratic and nonlinear programming,” Mathematical Programming, vol. 39, no. 2, pp. 117-129, Jun. 1987. [Baidu Scholar] 

12

P. A. Parrilo. (2000, May). Structured semidefinite programs and semialgebraic geometry methods in robustness and optimization. [Online]. Available: https://resolver.caltech.edu/CaltechETD:etd-05062004-05 5516 [Baidu Scholar] 

13

Z. Zhang, R. Schuerhuber, L. Fickert et al., “Domain of attraction’s estimation for grid connected converters with phase-locked loop,” IEEE Transactions on Power Systems, vol. 37, no. 2, pp. 1351-1362, Mar. 2022. [Baidu Scholar] 

14

C. Mishra, A. Pal, J. S. Thorp et al., “Transient stability assessment of prone-to-trip renewable generation rich power systems using Lyapunov’s direct method,” IEEE Transactions on Sustainable Energy, vol. 10, no. 3, pp. 1523-1533, Jul. 2019. [Baidu Scholar] 

15

B. Severino and K. Strunz, “Enhancing transient stability of DC microgrid by enlarging the region of attraction through nonlinear polynomial droop control,” IEEE Transactions on Circuits and Systems I: Regular Papers, vol. 66, no. 11, pp. 4388-4401, Nov. 2019. [Baidu Scholar] 

16

L. Xu, J. M. Guerrero, A. Lashab et al., “A review of DC shipboard microgrids: part II: control architectures, stability analysis, and protection schemes,” IEEE Transactions on Power Electronics, vol. 37, no. 4, pp. 4105-4120, Apr. 2022. [Baidu Scholar] 

17

A. Maulik and D. Das, “Stability constrained economic operation of islanded droop-controlled DC microgrids,” IEEE Transactions on Sustainable Energy, vol. 10, no. 2, pp. 569-578, Apr. 2019. [Baidu Scholar] 

18

Y. Yang, C. Li, J. Xu et al., “Virtual inertia control strategy for improving damping performance of DC microgrid with negative feedback effect,” IEEE Journal of Emerging and Selected Topics in Power Electronics, vol. 9, no. 2, pp. 1241-1257, Apr. 2021. [Baidu Scholar] 

19

P. Lin, C. Zhang, P. Wang et al., “A decentralized composite controller for unified voltage control with global system large-signal stability in DC microgrids,” IEEE Transactions on Smart Grid, vol. 10, no. 5, pp. 5075-5091, Sept. 2019. [Baidu Scholar] 

20

X. Li, M. Wang, W. Jiang et al., “Toward large-signal stabilization of interleaved floating multilevel boost converter-enabled high-power DC microgrids supplying constant power loads,” IEEE Transactions on Industrial Electronics, vol. 71, no. 1, pp. 857-869, Jan. 2024. [Baidu Scholar] 

21

X. Zhao, L. Guo, L. Zhu et al., “Power feasible region ensuring transient stability of droop-based multiconverters DC system,” IEEE Transactions on Power Electronics, vol. 38, no. 4, pp. 5442-5455, Apr. 2023. [Baidu Scholar] 

22

H. Wang, M. Han, R. Han et al., “A decentralized current-sharing controller endows fast transient response to parallel DC-DC converters,” IEEE Transactions on Power Electronics, vol. 33, no. 5, pp. 4362-4372, May 2018. [Baidu Scholar] 

23

R. W. Erickson and D. Maksimovic, Fundamentals of Power Electronics. New York: Springer Science & Business Media Press, 2001. [Baidu Scholar] 

24

Q. Xu, Y. Yan, C. Zhang et al., “An offset-free composite model predictive control strategy for DC/DC buck converter feeding constant power loads,” IEEE Transactions on Power Electronics, vol. 35, no. 5, pp. 5331-5342, May 2020. [Baidu Scholar] 

25

J. Liu, W. Zhang, and G. Rizzoni, “Robust stability analysis of DC microgrids with constant power loads,” IEEE Transactions on Power Systems, vol. 33, no. 1, pp. 851-860, Jan. 2018. [Baidu Scholar] 

26

H. Khalil, Nonlinear Systems. New York: Pearson Education Group Press, 2001. [Baidu Scholar] 

27

G. Chesi, Domain of Attraction: Analysis and Control via SOS Programming. New York: Springer Press, 2011. [Baidu Scholar] 

28

Z. W. Jarvis-Wloszek. (2003, Jun.). Lyapunov based analysis and controller synthesis for polynomial systems using sum-of-squares optimization. [Online]. Available: https://xueshu.baidu.com/usercenter/paper/show?paperid=02d53d6dbd290fa167502a6420a60c31&site=xueshu_se [Baidu Scholar] 

29

M. Anghel, F. Milano, and A. Papachristodoulou, “Algorithmic construction of Lyapunov functions for power system stability analysis,” IEEE Transactions on Circuits and Systems I: Regular Papers, vol. 60, no. 9, pp. 2533-2546, Sept. 2013. [Baidu Scholar] 

30

G. Stengle, “A nullstellensatz and a positivstellensatz in semialgebraic geometry,” Mathematische Annalen, vol. 207, pp. 87-97, Jun. 1974. [Baidu Scholar] 

31

D. Bertsekas.(2016, Jan.). Nonlinear programming. [Online]. Available: http://athenasc.com/nonlinbook.html [Baidu Scholar] 

32

P. Seiler and G. J. Balas, “Quasiconvex sum-of-squares programming,” in Proceedings of 49th IEEE Conference on Decision and Control, Atlanta, USA, Dec. 2010, pp. 3337-3342. [Baidu Scholar] 

33

J. Nie and L. Wang, “Regularization methods for SDP relaxations in large-scale polynomial optimization,” SIAM Journal on Optimization, vol. 22, no. 2, pp. 408-428, Jan. 2012. [Baidu Scholar] 

34

A. A. Ahmadi and A. Majumdar, “DSOS and SDSOS optimization: more tractable alternatives to sum of squares and semidefinite optimization,” SIAM Journal on Applied Algebra and Geometry, vol. 3, no. 2, pp. 193-230, Jan. 2019. [Baidu Scholar]