Estimation of composite load model with aggregate induction motor dynamic load for an isolated hybrid power system

Abstract

It is well recognized that the voltage stability of a power system is affected by the load model and hence, to effectively analyze the reactive power compensation of an isolated hybrid wind-diesel based power system, the loads need to be considered along with the generators in a transient analysis. This paper gives a detailed mathematical modeling to compute the reactive power response with small voltage perturbation for composite load. The composite load is a combination of the static and dynamic load model. To develop this composite load model, the exponential load is used as a static load model and induction motors (IMs) are used as a dynamic load model. To analyze the dynamics of IM load, the fifth, third and first order model of IM are formulated and compared using differential equations solver in Matlab coding. Since the decentralized areas have many small consumers which may consist large numbers of IMs of small rating, it is not realistic to model either a single large rating unit or all small rating IMs together that are placed in the system. In place of using a single large rating IM, a group of motors are considered and then the aggregate model of IM is developed using the law of energy conservation. This aggregate model is used as a dynamic load model. For different simulation studies, especially in the area of voltage stability with reactive power compensation of an isolated hybrid power system, the transfer function \mathrm{\Delta }Q/\mathrm{\Delta }V of the composite load is required. The transfer function of the composite load is derived in this paper by successive derivation for the exponential model of static load and for the fifth and third order IM dynamic load model using state space model.

Graphical abstract

Keywords

isolated hybrid power system (IHPS) / composite load model / static load / dynamic load / induction motor load model / aggregate load

1 Introduction

An isolated far located remote area where the availability of grid connected power supply is almost impossible and uneconomical; a wind-diesel based isolated hybrid power system (IHPS) is most promising to provide continuous, efficient, economical and reliable electrical energy [1]. In such a system, wind and diesel based electrical sources may be operated by induction generators and synchronous generators respectively. The advantages of such a hybrid system have been explored [2,3]. Apart from balancing the electric power supply and demand, the ancillary services are also important issues for the system operations. The reactive power issue is one of them [4]. When the load and/or input reactive power demand fluctuates, the voltage will also change. Without any compensation, this voltage variation may go beyond the voltage permissible range and therefore such power would not be acceptable for end users [5], and therefore static and dynamic compensation techniques are required. Electrical energy generation and consumption is a combined process for any power system. In previous works based on the IHPS, the dynamics of wind and diesel based generators have been depicted in detail. However, the load dynamics have not been considered in detail and even ignored by considering the static load model in the system [6].

Loads, however, have a significant impact on a system. It has been concluded that the voltage stability of a power system and choice of compensation techniques significantly depends on the selection of the load model and its parameters [7]. To effectively analyze the dynamics of the IHPS, the loads need to be considered along with the generators in a transient analysis. Power system planners and operators attempt to accurately model loads in order to analyze their systems. However, it is very difficult to exactly describe the loads in a mathematical model. A load consists of several components that have very different dynamic characteristics [8]. The information/knowledge about load model parameters, that properly depict load behavior during electric power system disturbances, makes it possible to properly plan the power system, reliably predict the prospective operating scenarios and select the adequate control actions in order to prevent undesired system behavior and ultimately system instability. Accurate load modeling is important to correctly predict the response of the system to disturbances. With a poor load model, the system needs to be operated with a higher safety margin [9]. Load models, which quantify real and reactive power responses to voltage and frequency disturbances, are generally divided in two groups—static load models (SLM) and dynamic load models (DLM). The SLM and DLM are classified according to the effect of the voltage on the load. If the load variation depends only on the instantaneous voltage input and is unrelated to the preceding voltage inputs, the static load model is used. However, if the load characteristics are affected by all of the voltage inputs over time, the DLM needs to be used [8]. The SLM are generally used for the calculation of steady-state conditions and in steady-state simulations of the power system. Since they do not take into account load dynamics, SLMs alone are not suitable for voltage stability analysis of power systems. The DLM are, therefore, necessary for dynamics studies, i.e., the analysis of power system behavior following small or large disturbances [7]. Particularly for an industrial demand, most of the load is dynamic but in the commercial and the residential demand, the percentage of static load is higher than that of dynamic load and might be in a percent ratio of 4 to 1 [9]. Many studies have been conducted of aggregate or composite load modeling which includes the composite load as a combination of the SLM and DLM. In Ref. [10], three composite load models are depicted, including the ZIP-induction motor model, the exponential-induction motor model and the Z-induction motor model. For the dynamic induction motor model, the third order induction motor model is studied [8].

Since the decentralized areas have many consumers which may consist large numbers of induction motors, it is not realistic to model every induction motor that is in the system. However, it is impractical to accurately represent each individual load due to the intense computation process involved. Hence, aggregate models or single-unit models with a minimum order of induction motor are needed to represent a group of motors [11]. The appropriate dynamic load model aggregation reduces the computation time and provides a faster and efficient model derivation and parameter identification. It is found that the small-scale aggregation model gives acceptably accurate results than the large-scale aggregation model and is good for power system stability analysis [12]. Hence, in place of using a single large rating induction motor, a group of motors are considered and then the aggregate model of the induction motor is developed using the law of energy conservation [11].

Hence, this paper presents the detailed mathematics for composite load which includes the static and dynamic load. The exponential load is modeled as a static load while the aggregate induction motor is modeled as a dynamic load. Induction motor characteristics are compared for its fifth, third and first order model. The transfer function of voltage variations with reactive power change is derived using state space equations for induction motor and then the same is derived for the composite load.

2 Load interaction in IHPS

The voltage control of a wind-diesel based IHPS is discussed in many previous papers [13−16]. The voltage may be controlled by using either the dynamic reactive power compensation alone or the combination of the static and dynamic reactive power compensation. The compensator participation schemes depend on the voltage stability and compensator installation cost. Due to the disturbance in load (\mathrm{\Delta }{Q}_L) and/or input (\mathrm{\Delta }{P}_{\mathrm{in}}), the reactive power demand will change in the IHPS, and the transients will appear in the load voltage of the system. To compensate this voltage, an incremental change in the reactive power of other components will activate and try to stabilize this voltage. The overall incremental change \mathrm{\Delta }Q in reactive power is caused by the incremental change in reactive power generated by the synchronous generator (\mathrm{\Delta }{Q}_{\mathrm{SG}}) (by the compensator (\mathrm{\Delta }{Q}_{\mathrm{com}}), and absorbed by the induction generator \mathrm{\Delta }{Q}_{\mathrm{IG}}.

Since, the net incremental change in reactive power affects the load voltage disturbances (\mathrm{\Delta }V) and electro-magnetic energy absorption in the induction generator (E_m). Equation (1) depicts the reactive power balance equation for an isolated wind-diesel based hybrid power system [5].

(1)\mathrm{\Delta }{Q}_{\mathrm{SG}}+\mathrm{\Delta }{Q}_{\mathrm{com}}-\mathrm{\Delta }{Q}_{\mathrm{IG}}-\mathrm{\Delta }{Q}_L={\displaystyle \frac{d{E}_m}{dt}}+D_v\mathrm{\Delta }V

where D_v is the load transfer function of reactive power change to voltage change. In previous papers exponential form of static and dynamic load model is considered [5,17]. The load pattern in isolated remote area is generally residential and commercial. Reference [18] proposed a composite load model which consists of 80% of the SLM and 20% of the DLM. The composite load structure is illustrated in Fig. 1. The residential, commercial, and small scale industrial consumers may have large number of induction motors of small rating, so it is not realistic to model a single large rating unit. The dynamic load is a combination of several induction motors of small rating and therefore, an aggregate dynamic load model is required. This composite load is the combination of a static load and an aggregate unit for several induction motors of small ratings. The term D_v, for the composite load model (CLM) shown in Fig. 1 can be evaluated by adding the transfer functions of the SLM and the DLM and is given by Eq. (2).

(2){(D_v)}_{\mathrm{CLM}}={(D_v)}_{\mathrm{SLM}}+{(D_v)}_{\mathrm{DLM}}.

3 Composite load model structure

As the composite load model is composed of the static exponential load model in parallel with the dynamic induction motor load model as shown in Fig. 1. The expression of D_{\mathrm{v}} for the CLM can be calculated according to Eq. (2) if the transfer function for the SLM and the DLM are known separately. A detailed methodology for estimating the transfer function \mathrm{\Delta }Q/\mathrm{\Delta }V for the SLM and an induction motor model (DLM) is illustrated in Sub-sections 3.1 and 3.2, respectively. The detailed methodology for getting the aggregate model of the induction motor is illustrated in Sub-section 3.3.

For this study, a total composite load of 250 kW is assumed to be connected at the IHPS. The ratio of static load to dynamic load in the system is taken as 4:1 [18] and hence the SLM and the DLM load powers are considered to be 200 kW and 50 kW, respectively. The other details about load data are provided later. The base value of power, voltage and frequency are set to 250 kW, 326.6 (\sqrt{\mathrm{2}}\times {\displaystyle \frac{\mathrm{400}}{\sqrt{\mathrm{3}}}}) V and 50 Hz.

3.1 SLM analysis

The static load is represented by an exponential function of bus voltage magnitude as given in Eq. (3). The reactive power expression for the SLM is

(3)Q_{\mathrm{SLM}}=Q_{\mathrm{0},\mathrm{SLM}}{\left({\displaystyle \frac{V}{V_{\mathrm{0}}}}\right)}^{n_q}.

The expression for the transfer function of reactive power change to voltage change D_{\mathrm{v}} is given in Eq. (4).

(4){\left(D_{\mathrm{v}}\right)}_{\mathrm{SLM}}={\displaystyle \frac{n_q{Q}_{\mathrm{0}}}{V_{\mathrm{0}}}}.

The exponent coefficient n_q for the reactive power of load may vary from 0 to 7 and 4 to 7 for steady-state and transient conditions, respectively [7,19,20]. The value of load power factor is taken as 0.9 lagging and the value of n_q is set to 3 in this paper [10].

3.2 DLM analysis

Induction motor is the most commonly used dynamic load. The induction motor load may be modeled using the fifth, third and first order model [21]. The load models proposed in Ref. [21] have been verified by developing code in Matlab. The responses of speed, electro-magnetic torque, power, and currents are obtained using differential equation solver ode45 in Matlab. A methodology to achieve the transfer function of reactive power change to voltage change, D_{\mathrm{v}} as depicted in Eq. (2), is explained step by step.

Step 1 Collection of the manufacturer data.

Commonly, the induction motor is defined by the output parameters, given in the catalogs as manufacturer data. These output values can be used for the calculations and to determine the parameters of any proposed equivalent circuits.

Step 2 Identification of equivalent circuit parameters of induction motor.

The equivalent circuit parameters of the induction motor are identified in Fig. 2. The initial reference parameters are used in this paper. The method of getting these initial reference parameters are well documented in Ref. [22]. For squirrel cage, the induction motor rotor is short circuited, and therefore, V_{\mathrm{r}}=\mathrm{0}. R_s, R_r, L_s, L_r and L_m are the circuit parameters of the equivalent circuit of the induction motor referred to the stator side.

Step 3 Estimation of induction motor responses.

The knowledge of induction motor responses is essential for dynamic load modeling [23]. The induction motor dynamic modeling can be achieved by using the first, third or fifth order model. The fifth order model is very close to the real motor while the third and first order model are the simple mathematical version of motor dynamics. To identify the induction motor dynamic responses, the fifth order model is mathematically expressed in Eqs. (5) to (9) [24]. For the balanced operation of squirrel cage induction motor, the most widely used reference frame is synchronous rotating reference frame which is used in this paper [9]. All parameters are in per unit quantities except {\omega }_{\mathrm{s}},{\omega }_{\mathrm{b}},{\omega }_{\mathrm{r}} which are in radian per seconds.

(5){\displaystyle \frac{\mathrm{d}{\phi }_{q\mathrm{s}}}{\mathrm{d}t}}={\omega }_{\mathrm{b}}\left(V_{q\mathrm{s}}-{\displaystyle \frac{{\omega }_{\mathrm{s}}}{{\omega }_{\mathrm{b}}}}{\phi }_{d\mathrm{s}}-R_{\mathrm{s}}{I}_{q\mathrm{s}}\right),

(6){\displaystyle \frac{\mathrm{d}{\phi }_{d\mathrm{s}}}{\mathrm{d}t}}={\omega }_{\mathrm{b}}\left(V_{d\mathrm{s}}+{\displaystyle \frac{{\omega }_{\mathrm{s}}}{{\omega }_{\mathrm{b}}}}{\phi }_{q\mathrm{s}}-R_{\mathrm{s}}{I}_{d\mathrm{s}}\right),

(7){\displaystyle \frac{\mathrm{d}{\phi }_{q\mathrm{r}}}{\mathrm{d}t}}={\omega }_{\mathrm{b}}\left(V_{q\mathrm{r}}-{\displaystyle \frac{{\omega }_{\mathrm{s}}-{\omega }_{\mathrm{r}}}{{\omega }_{\mathrm{s}}}}{\phi }_{d\mathrm{r}}-R_{\mathrm{r}}{I}_{q\mathrm{r}}\right),

(8){\displaystyle \frac{\mathrm{d}{\phi }_{d\mathrm{r}}}{\mathrm{d}t}}={\omega }_{\mathrm{b}}\left(V_{d\mathrm{r}}+{\displaystyle \frac{{\omega }_{\mathrm{s}}-{\omega }_{\mathrm{r}}}{{\omega }_{\mathrm{s}}}}{\phi }_{q\mathrm{r}}-R_{\mathrm{r}}{I}_{d\mathrm{r}}\right),

(9){\displaystyle \frac{\mathrm{d}{\omega }_{\mathrm{r}}}{\mathrm{d}t}}={\displaystyle \frac{{\omega }_{\mathrm{b}}}{\mathrm{2}H}}\left(T_{\mathrm{e}}-T_{\mathrm{L}}\right).

In Eqs. (5) to (9), the direct and quadrature axis voltages are the independent variables and fluxes are the dependent variables. The conversion of supply three phase voltages into d-q axis can be done using the mathematical expression given in Eqs. (10) to (16) [25].

For squirrel cage induction motor, the rotor side is short circuited, so the direct and quadrature axis values of the rotor are zero,

(10)V_{d\mathrm{r}}=V_{q\mathrm{r}}=\mathrm{0.}

For a balanced three phase stator voltages,

(11)V_{\mathrm{s}a}=\sqrt{\mathrm{2}}{V}_{\mathrm{0}}\cos {\omega }_{\mathrm{b}}t,

(12)V_{\mathrm{s}b}=\sqrt{\mathrm{2}}{V}_{\mathrm{0}}\cos \left({\omega }_{\mathrm{b}}t-{\displaystyle \frac{\mathrm{2}\mathrm{\pi }}{\mathrm{3}}}\right),

(13)V_{\mathrm{s}c}=\sqrt{\mathrm{2}}{V}_{\mathrm{0}}\cos \left({\omega }_{\mathrm{b}}t+{\displaystyle \frac{\mathrm{2}\mathrm{\pi }}{\mathrm{3}}}\right).

For \mathrm{constant}={\mathrm{e}}^{-\mathrm{j}\frac{\mathrm{2}\mathrm{\pi }}{\mathrm{3}}} , the stator voltage is given by

(14)V_{\mathrm{s}}={\displaystyle \frac{\mathrm{2}}{\mathrm{3}}}\left(V_{\mathrm{s}a}+a{V}_{\mathrm{s}b}+a^{\mathrm{2}}{V}_{\mathrm{s}c}\right).

This stator voltage is still in stator coordinates. The conversion of this stator voltage into synchronous rotating reference frame is

(15)V=V_{\mathrm{s}}{\mathrm{e}}^{-\mathrm{j}{\omega }_{\mathrm{s}}t}.

The complex variables may be decomposed in the plane along two orthogonal d and q axes rotating at a speed of {\omega }_{\mathrm{b}} to obtain the d-q (Park) model,

(16)V=V_{d\mathrm{s}}+\mathrm{j}{V}_{q\mathrm{s}}.

The dependent variable can be estimated using the input supply voltage. The expressions for fluxes in terms of currents, per unit electro-magnetic torque and per unit load torque are depicted in Refs. [24–25], which are given in Eqs. (17) to (24).

(17)L_{\mathrm{ss}}=L_{\mathrm{s}}+L_{\mathrm{m}},

(18)L_{\mathrm{rr}}=L_{\mathrm{r}}+L_{\mathrm{m}},

(19){\phi }_{q\mathrm{s}}=L_{\mathrm{ss}}{I}_{q\mathrm{s}}+L_{\mathrm{m}}{I}_{q\mathrm{r}},

(20){\phi }_{d\mathrm{s}}=L_{\mathrm{ss}}{I}_{d\mathrm{s}}+L_{\mathrm{m}}{I}_{d\mathrm{r}},

(21){\phi }_{q\mathrm{r}}=L_{\mathrm{rr}}{I}_{q\mathrm{r}}+L_{\mathrm{m}}{I}_{q\mathrm{s}},

(22){\phi }_{d\mathrm{r}}=L_{\mathrm{rr}}{I}_{d\mathrm{r}}+L_{\mathrm{m}}{I}_{d\mathrm{s}},

(23)T_{\mathrm{e}}={\phi }_{q\mathrm{r}}{L}_{d\mathrm{r}}-{\phi }_{d\mathrm{r}}{L}_{q\mathrm{r}},

(24)T_{\mathrm{L}}=B{\displaystyle \frac{{\omega }_{\mathrm{r}}}{{\omega }_{\mathrm{b}}}}.

Among the parameters of the induction motor, the dynamics are largely characterized by the machine inertia H, moment of inertia J and torque-damping factor B. Many researchers have estimated these parameters. A nonlinear least square approach to find these parameters is depicted in Ref. [26]. The B and H values are estimated by analyzing the rotor speed and slip responses with two constraints:

1)‚Slow varying speed of lower overshoots;

2)‚Slip should reach to its steady-state value maintaining positive value of slip at all instants of response.

These parameters H and B are evaluated by developing codes on MATLAB 7.10. To solve these fifth order differential equations, the ode45 syntax is used in Matlab coding. To reduce it into the third order model, the stator direct and quadrature axis flux is set to zero in Eqs. (5) and (6). To further reduce it into the first order model, the rotor direct and quadrature axis flux is set to zero in Eqs. (7) and (8) [9].

Step 4 State space model of induction motor.

The transfer function of the induction motor load discussed in Eq. (2) is developed by using the linearizing techniques of the state space model. Figure 3 represents a generalized block diagram for state space model of any order (fifth, third or first) induction motor model. The control and disturbance vectors for model study are ΔV and ΔQ, respectively, and the relation is shown in Eq. (25) [10]. Since the value of V_{q\mathrm{s}} is very close to zero, ΔV is assumed to be equal to ΔV_ds and therefore V_{q\mathrm{s}} is neglected in the expression of reactive power (Q).

(25)Q=V_{q\mathrm{s}}{I}_{d\mathrm{s}}-V_{d\mathrm{s}}{I}_{q\mathrm{s}}.

The fifth, third and first order induction motor model has 5, 3 and 1 state vectors, respectively. The corresponding state variables are represented in their state space equation.

The state space equations for the fifth order model are

(26){\left[\mathrm{\Delta }{\dot{\phi }}_{q\mathrm{s}}\mathrm{\Delta }{\dot{\phi }}_{d\mathrm{s}}\mathrm{\Delta }{\dot{\phi }}_{q\mathrm{r}}\mathrm{\Delta }{\dot{\phi }}_{d\mathrm{r}}\mathrm{\Delta }{\dot{\omega }}_{\mathrm{r}}\right]}^{\mathrm{T}}=A{q}_{\mathrm{5}}{\left[\mathrm{\Delta }{\phi }_{q\mathrm{s}}\mathrm{\Delta }{\phi }_{d\mathrm{s}}\mathrm{\Delta }{\phi }_{q\mathrm{r}}\mathrm{\Delta }{\phi }_{d\mathrm{r}}\mathrm{\Delta }{\omega }_{\mathrm{r}}\right]}^{\mathrm{T}}+B{q}_{\mathrm{5}}\mathrm{\Delta }V,

(27)\mathrm{\Delta }Q=C{q}_{\mathrm{5}}{\left[\mathrm{\Delta }{\phi }_{q\mathrm{s}}\mathrm{\Delta }{\phi }_{d\mathrm{s}}\mathrm{\Delta }{\phi }_{q\mathrm{r}}\mathrm{\Delta }{\phi }_{d\mathrm{r}}\mathrm{\Delta }{\omega }_{\mathrm{r}}\right]}^{\mathrm{T}}+D{q}_{\mathrm{5}}\mathrm{\Delta }V.

The state space equations for the third order model are

(28){\left[\mathrm{\Delta }{\dot{\phi }}_{q\mathrm{r}}\mathrm{\Delta }{\dot{\phi }}_{d\mathrm{r}}\mathrm{\Delta }{\dot{\omega }}_{\mathrm{r}}\right]}^{\mathrm{T}}=A{q}_{\mathrm{3}}{\left[\mathrm{\Delta }{\phi }_{q\mathrm{r}}\mathrm{\Delta }{\phi }_{d\mathrm{r}}\mathrm{\Delta }{\omega }_{\mathrm{r}}\right]}^{\mathrm{T}}+B{q}_{\mathrm{3}}\mathrm{\Delta }V,

(29)\mathrm{\Delta }Q=C{q}_{\mathrm{3}}{\left[\mathrm{\Delta }{\phi }_{q\mathrm{r}}\mathrm{\Delta }{\phi }_{d\mathrm{r}}\mathrm{\Delta }{\omega }_{\mathrm{r}}\right]}^{\mathrm{T}}+D{q}_{\mathrm{3}}\mathrm{\Delta }V.

The state space equations for the first order model are

(30)\mathrm{\Delta }{\dot{\omega }}_{\mathrm{r}}=A{q}_{\mathrm{1}}\mathrm{\Delta }{\omega }_{\mathrm{r}}+B{q}_{\mathrm{1}}\mathrm{\Delta }V,

(31)\mathrm{\Delta }Q=C{q}_{\mathrm{1}}\mathrm{\Delta }{\omega }_{\mathrm{r}}+D{q}_{\mathrm{1}}\mathrm{\Delta }V.

The linear state space model for the fifth, third and first order model can be derived from Eqs. (5) to (25). The A, B, C and D are the constant matrices of the appropriate dimensions associated with the above control, state and disturbances vectors. The order of each matrix of the state space equation is given in Table 1. The elements of these matrices shown in Table 2 are evaluated by the knowledge of parameters identified through Steps 1, 2 and 3 in Sub-section 3.2.

The linearization of these multistate models is performed by Taylor series expansion of the nonlinear functions and neglecting the high order terms. If ε notation is used for state vectors and k notation is used for denoting the model order, the value of matrices can be formulated as

(32)A{q}_k=\left[\begin{array}{cccc}\hfill {\displaystyle \frac{\partial \dot{{ϵ}_{\mathrm{1}}}}{\partial {ϵ}_{\mathrm{1}}}}\hfill & \hfill {\displaystyle \frac{\partial \dot{{ϵ}_{\mathrm{1}}}}{\partial {ϵ}_{\mathrm{2}}}}\hfill & \hfill \cdots \hfill & \hfill {\displaystyle \frac{\partial \dot{{ϵ}_{\mathrm{1}}}}{\partial {ϵ}_k}}\hfill \\ \hfill ⋮\hfill & \hfill ⋮\hfill & \hfill \cdots \hfill & \hfill ⋮\hfill \\ \hfill {\displaystyle \frac{\partial \dot{{ϵ}_k}}{\partial {ϵ}_{\mathrm{1}}}}\hfill & \hfill {\displaystyle \frac{\partial \dot{{ϵ}_k}}{\partial {ϵ}_{\mathrm{2}}}}\hfill & \hfill \cdots \hfill & \hfill {\displaystyle \frac{\partial \dot{{ϵ}_k}}{\partial {ϵ}_k}}\hfill \end{array}\right],

(33)B{q}_k=\left[\begin{array}{c}\hfill {\displaystyle \frac{\partial \dot{{ϵ}_{\mathrm{1}}}}{\partial V}}\hfill \\ \hfill ⋮\hfill \\ \hfill {\displaystyle \frac{\partial \dot{{ϵ}_k}}{\partial V}}\hfill \end{array}\right],

(34)C{q}_k=\left[\begin{array}{ccc}\hfill {\displaystyle \frac{\partial Q}{\partial {ϵ}_{\mathrm{1}}}}\hfill & \hfill \cdots \hfill & \hfill {\displaystyle \frac{\partial Q}{\partial {ϵ}_k}}\hfill \end{array}\right],

(35)D{q}_k={\displaystyle \frac{\partial Q}{\partial V}}.

Step 5 Configuration of transfer function of reactive power change to voltage change.

These state space models in Step 4 are used to evaluate the transfer function in Matlab. The two syntaxes used in Matlab coding for finding the transfer function of induction motor for order k (where, k=\left[\mathrm{1},\mathrm{3},\mathrm{5}\right]) are given in Eqs. (36) and (37).

(36)\left[\mathrm{num}{q}_k,\mathrm{den}{q}_k\right]=\mathrm{ss}\mathrm{2}\mathrm{tf}\left(A{q}_k,B{q}_k,C{q}_k,D{q}_k\right),

(37){\left(D_{\mathrm{v}}\right)}_{\mathrm{IM}}=\mathrm{tf}\left(\mathrm{num}{q}_k,\mathrm{den}{q}_k\right).

3.3 Aggregate DLM analysis

The step by step method for estimating D_{\mathrm{v}} for an induction motor is documented in Sub-section 3.2. This procedure may be adopted to find {\left(D_{\mathrm{v}}\right)}_{\mathrm{IM}}\mathrm{=tf}\left(\mathrm{num}{q}_k,\mathrm{den}{q}_k\right)D_{\mathrm{v}} for any single unit of the induction motor. The real load has a variety of small induction motors which are connected in any residential or commercial load profile. Hence, the dynamic load may be defined as the combination of several induction motors of small rating in any power system. To analyze any simulation model, an aggregate motor model may be used to demonstrate the set of several induction motors of small rating which are connected at load end. The parameters of the aggregate motor are derived using the energy conservation law. Figure 2 illustrates the equivalent circuit of a single unit of the induction motor. If n units of induction motors of different ratings are connected to any load end as shown in Fig. 1, the following procedure may be used to identify the aggregate model parameters [11,27].

If S_i is the rated kVA of each motor, the kVA power absorbed by the aggregate motor is equal to the sum of the individual motor kVA absorbed. The aggregate kVA is

(38)S^{\mathrm{agg}}=\sum _{k=\mathrm{1}}^n{S}_k.

Similarly, the aggregate stator and rotor current are

(39)\overrightarrow{I_{\mathrm{s}}^{\mathrm{agg}}}=\sum _{k=\mathrm{1}}^n\overrightarrow{I_{\mathrm{s}k}},

(40)\overrightarrow{I_{\mathrm{r}}^{\mathrm{agg}}}=\sum _{k=\mathrm{1}}^n\overrightarrow{I_{\mathrm{r}k}}.

Using the energy conservation law again, the aggregate motor load equivalent circuit parameters are

(41)R_{\mathrm{s}}^{\mathrm{agg}}={\displaystyle \frac{{\displaystyle \sum _{k=\mathrm{1}}^n{\left|\overrightarrow{I_{\mathrm{s}k}}\right|}^{\mathrm{2}}{R}_{\mathrm{s}k}}}{{\left|\overrightarrow{I_{\mathrm{s}}^{\mathrm{agg}}}\right|}^{\mathrm{2}}}},

(42)R_{\mathrm{r}}^{\mathrm{agg}}={\displaystyle \frac{{\displaystyle \sum _{k=\mathrm{1}}^n{\left|\overrightarrow{I_{\mathrm{r}k}}\right|}^{\mathrm{2}}{R}_{\mathrm{r}k}}}{{\left|\overrightarrow{I_{\mathrm{r}}^{\mathrm{agg}}}\right|}^{\mathrm{2}}}},

(43)X_{\mathrm{s}}^{\mathrm{agg}}={\displaystyle \frac{{\displaystyle \sum _{k=\mathrm{1}}^n{\left|\overrightarrow{I_{\mathrm{s}k}}\right|}^{\mathrm{2}}{X}_{\mathrm{s}k}}}{{\left|\overrightarrow{I_{\mathrm{s}}^{\mathrm{agg}}}\right|}^{\mathrm{2}}}},

(44)X_{\mathrm{r}}^{\mathrm{agg}}={\displaystyle \frac{{\displaystyle \sum _{k=\mathrm{1}}^n{\left|\overrightarrow{I_{\mathrm{r}k}}\right|}^{\mathrm{2}}{X}_{\mathrm{r}k}}}{{\left|\overrightarrow{I_{\mathrm{r}}^{\mathrm{agg}}}\right|}^{\mathrm{2}}}},

(45)X_{\mathrm{m}}^{\mathrm{agg}}={\displaystyle \frac{{\displaystyle \sum _{k=\mathrm{1}}^n{\left|\overrightarrow{I_{sk}}-\overrightarrow{I_{rk}}\right|}^{\mathrm{2}}{X}_{\mathrm{m}k}}}{{\left|\overrightarrow{I_{\mathrm{s}}^{\mathrm{agg}}}-\overrightarrow{I_{\mathrm{r}}^{\mathrm{agg}}}\right|}^{\mathrm{2}}}}.

The air gap power of the aggregate motor load is expressed as

(46)P_{\mathrm{airgap}}^{\mathrm{agg}}={\displaystyle \sum }_{k=\mathrm{1}}^n\left(\mathrm{Re}\left(\overrightarrow{V_{\mathrm{s}k}}\overrightarrow{I_{\mathrm{s}k}^{*}}\right)-{\left|\overrightarrow{I_{\mathrm{s}k}}\right|}^{\mathrm{2}}{R}_{\mathrm{s}k}\right).

The slip of the aggregate motor can then be computed by

(47){\mathrm{slip}}^{\mathrm{agg}}={\displaystyle \frac{{\left|\overrightarrow{I_{\mathrm{r}}^{\mathrm{agg}}}\right|}^{\mathrm{2}}{R}_{\mathrm{r}}^{\mathrm{agg}}}{P_{\mathrm{airgap}}^{\mathrm{agg}}}},

(48)H^{\mathrm{agg}}={\displaystyle \frac{{\displaystyle \sum _{k=\mathrm{1}}^n{H}_k{S}_k}}{S^{\mathrm{agg}}}}.

The moment of inertia and inertia constant of the motor has the following relation,

(49)H_k{S}_k={\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}{J}_k{\omega }_{\mathrm{s}k}^{\mathrm{2}}.

Therefore, the aggregate motor model can be estimated for the given set of several motors. This aggregate motor model can be used to identify the transfer function of aggregate DLM as described in Eq. (2).

4 Results and discussion

The Matlab program is developed for estimating, analyzing and identifying the results of this paper. The most realistic and probable nature of the composite load is the combination of the static and dynamic load. In an isolated far located remote area, the load pattern depends on the power requirement of residential and commercial sector load. In such cases, the participation of the static and dynamic load is taken as 80% and 20%, respectively.

Equation (3) shows the dependence of the instantaneous reactive power upon the voltage of the system for the exponential model of static load. Equation (4) gives the transfer function of the static load model {(D_v)}_{\mathrm{SLM}}. Induction motors are considered for the dynamic load model. Six induction motors of different ratings are assumed to be connected at the load end as shown in Fig. 1. The aggregate motor rating is 50 kW, which is equal to the sum of the rating of individual motor. A single unit of a 50 kW induction motor is also listed among six motors for verifying the results of the aggregate motor model by the energy conservation law. The manufactures data for all induction motors are assumed to be equal for simplifying the system study and listed in Table 2.

The equivalent circuit diagram of the squirrel cage induction motor is shown in Fig. 2. The initial reference parameters for the machines are calculated using the Matlab program, and the results are listed in Table 3. The responses of the fifth, third and first order model of the induction motor are found using Eqs. (5) to (24). The command ‘ode45’ is used in coding to solve the differential equations given in Eqs. (5) to (9). All time dependent state variables as mentioned in Eqs. (25) to (30) have been evaluated. The responses of rotor speed, slip and torques are used to estimate the dynamic parameters J/H and B of the induction motor as discussed in Step 3 in Sub-section 3.2. The estimated values of J and B for the individual induction motor are also given in Table 3 along with other calculate parameters of the machines.

Based on the data provided by the manufacturer in Table 2 and calculated parameters of induction motor in Table 3, the fifth, third and first order transient responses for speed, slip, electromagnetic torque, load torque, direct/quadrature axis stator and rotor currents, direct/quadrature axis stator and rotor flux, active power and reactive power are achieved through ode45 solver in the Matlab coding. The rotor speed, electromagnetic torque, active power and reactive power responses of the 50 kW induction motor are shown in Figs. 4 to 7. The responses of motor for the fifth, third and first order model are also compared in Figs. 4 to 7.

Concluding from Figs. 4 to 7, the behavior of the third order model is quite similar to that of the fifth order model while the first order model is out of context. So, it is advisable to use either the fifth or the third order model of the induction motor. As the real induction motor behaves like a fifth order model, the fifth order model should be considered for the exact load modeling of the induction motor. The third order model is achieved by neglecting the stator flux transients. So, the fifth order induction motor model responses have high transients and more settling time compare to the third order model. The step response comparison after evaluating the transfer function for fifth and the third order 50 kW IM model are compared and given in Fig. 8, which shows an almost similar behavior of both the models of the induction motor. The mathematical model of the fifth order model involves 5 differential equations to solve. The state space equations for the fifth order model are more complex as shown in Table 1, and the transfer function of the fifth order load model is more complex to compute than that of the third order model. To reduce this computational complexity, the third order load model of the induction motor may also be preferred and therefore it is also used to find the transfer function of the composite load in this paper.

Equations (26) to (31) are the state space equations for predefined state vectors, control vector and disturbance vector corresponding to the induction motor load model of different orders. Equations (32) to (35) are the Jacobian matrices that help to identify the transfer functions by applying two syntaxes in Matlab coding given in Eqs. (36) and (37). The bode plots for the fifth and the third order model, showing the effects of voltage change on reactive power change for linearized transfer function of 50 kW induction motor are shown here for better understanding of the response of both orders in Figs. 9 and 10 respectively.

The actual dynamic load may be defined as the combination of several induction motors of small rating. The real dynamic load is represented by an aggregate induction motor load. The load parameters of this aggregate induction motor are calculated by clubbing the induction motor of six different ratings (IM1-IM6) as given in Table 2. The dynamic load model parameters of the aggregate motor are identified by the law of energy conservation. Equations (38) to (49) show the mathematics for calculating the parameters of the aggregate model of the induction motor. Their corresponding results are listed in Table 3. The parameters of the 50 kW single unit induction motor and the aggregate motor model of 50 kW are given in Table 3. The circuit parameters of the aggregate model are the same as that of the 50 kW motor except the value of the moment of inertia and torque-damping factor which denotes that the aggregate model will also behave like the induction motor model but at different dynamics.

The bode plot comparison between the single unit and the aggregate unit of 50 kW IM load model also verifies the aggregate procedure in Figs. 11 and 12. The variation in the magnitude and phase for the single and the aggregate induction motor is caused by the different magnitude of the torque-damping factor and the moment of inertia. The transfer function of the reactive power change to the voltage change for this aggregate motor model is evaluated in the same way as discussed before.

The procedure for getting the transfer function of the aggregate motor load model of 50 kW is documented in detail in this paper. The composite load transfer function is achieved using Eq. (2) by adding the static load transfer function with this aggregate dynamic load transfer function. Table 4 illustrates the results of derived transfer functions of the load.

The transfer function of the static load model, the third order IM based dynamic load model, the fifth order IM based dynamic load model and the composite load by adding these static and dynamic load models are given in Table 4. The corresponding step responses for the third and the fifth order model of the aggregate dynamic load, the static load and the composite load model are compared in Figs. 13 and 14. As the composite load is designed with 80% of the static load and 20% of the dynamic load compositions, most of the reactive power demand results from the static load. However, the transients in the reactive power demand are caused by the presence of the dynamic load demand. If the composition of the dynamic load increases, these transients will also increase. The reactive power compensation scheme will become more important in the case of the composite load connected in the power system. Hence, the composite load model has a significant impact on the study of the isolated hybrid power system. It may give more realistic results of the power system voltage stability and reactive power compensation. The transfer function derived for the composite load shown in Table 4 may further be used for reactive power balance schemes, voltage stability and compensation economics study in isolated hybrid power systems.

5 Conclusions

In the past, researchers did not focus on the load model when studying the reactive power compensation in the isolated hybrid power system. It is well documented that the load model is also an important factor for calculating the reactive power requirement and voltage stability of such systems. To analyze the system, a transfer function of load reactive power change to voltage change is required. A step by step method is developed and explained in this paper. A composite load model is considered which includes the combination of the static and the dynamic load model. The aggregate induction motor model is taken as the dynamic load model. This aggregate model of the induction motor is developed by clubbing several induction motors of small rating using the law of energy conservation. The aggregate motor model is compared with the single unit induction motor of equivalent rating. The fifth and the third order induction motor load model is developed using the solver ode45. The responses of the fifth and the third order are compared. The transfer functions for both the orders are evaluated using sate space equations. The composite load model transfer function is finally achieved by adding the transfer function of both the static and the dynamic load model. Both the fifth and the third order transfer functions are derived for the composite load. Though the behavior of both order transfer functions are almost the same, the final choice about the load model depend on the requirements of individuals.

References

[1] Sharma P, Kumar Saxena N, Ramakrishna K S S, Bhatti T S. Reactive power compensation of isolated wind-diesel hybrid power systems with STATCOM and SVC. International Journal on Electrical Engineering and Informatics, 2010, 2(3): 192–203

[2] Hunter R, Elliot G. Wind-diesel systems, a guide to the technology and its implementation. Cambridge: Cambridge University Press, 1994

[3] Cardenas R, Pena R, Perez M, Clare J, Asher G, Vargas F. Vector control of frond end converters for variable speed wind-diesel systems. IEEE Transactions on Industrial Electronics, 2006, 53(4): 1127–1136

[4] Bansal R C, Bhatti T S, Kumar V. Reactive power control of autonomous wind diesel hybrid power systems using ANN. In: Proceedings of the International Power Engineering Conference, Singapore, 2007, 982–987

[5] Bansal R C. Automatic reactive power control of autonomous hybrid power system. Dissertation for the Doctoral Degree. Delhi: Indian Institute of Technology, 2002

[6] Sharma P, Saxena N K, Bhatti T S. Study of autonomous hybrid power system using SVC and STATCOM. In: Proceedings of International Conference on Power Systems. Kharagpur, India, 2009, 27–29

[7] Stojanović D P, Korunović L M, Milanović J V. Dynamic load modelling based on measurements in medium voltage distribution network. Electric Power Systems Research, 2008, 78(2): 228–238

[8] Kim B H, Kim H, Lee B, 0. Kim H, Lee B. Parameter estimation for the composite load model. Journal of International Council on Electrical Engineering, 2012, 2(2): 215–218

[9] Parveen T. Composite load model decomposition: induction motor contribution. Dissertation for the Doctoral Degree. Brisbane: Queensland University of Technology, 2009

[10] Choi B K, Chiang H D, Li Y, Chen Y T, Huang D H, Lauby M G. Development of composite load models of power systems using on-line measurement data. Journal of Electrical Engineering & Technology, 2006, 1(2): 161–169

[11] Aree P. Aggregating method of induction motor group using energy conservation law. ECTI Transactions on Electrical & Computer Engineering, 2014, 12(1): 1–6

[12] Muriuki J K, Muriithi C M. Comparison of aggregation of small and large induction motors for power system stability study. Global Engineers & Technologists Review, 2013, 3(2): 9–13

[13] Sharma P, Bhatti T S. Performance investigation of isolated wind-diesel hybrid power systems with WECS having PMIG. IEEE Transactions on Industrial Electronics, 2013, 60(4): 1630–1637

[14] Bansal R C, Bhatti T S, Kothari D P. A novel mathematical modelling of induction generator for reactive power control of isolated hybrid power systems. International Journal of Modelling and Simulation, 2004, 24(1): 1–7

[15] Sharma P, Sulkowski W, Hoff B. Dynamic stability study of an isolated wind-diesel hybrid power system with wind power generation using IG, PMIG and PMSG: a comparison. International Journal of Electrical Power and Energy Systems, 2013, 53: 857–866

[16] Vachirasricirikul S, Ngamroo I, Kaitwanidvilai S. Coordinated SVC and AVR for robust voltage control in a hybrid wind-diesel system. Energy Conversion and Management, 2010, 51(12): 2383–2393

[17] Saxena N, Kumar A. Load modeling interaction on hybrid power system using STATCOM. In: Proceedings of Annual IEEE India Conference. Kolkata, India, 2010.

[18] Kosterev D, Meklin A. Load modelling in WECC. Power Systems Conference and Exposition (PSCE’06), 2006, 576–581

[19] Fahmy O M, Attia A S, Badr M A L. A novel analytical model for electrical loads comprising static and dynamic components. Electric Power Systems Research, 2007, 77(10): 1249–1256

[20] Hiskens I A, Milanovic J V. Load modelling in studies of power system damping. IEEE Transactions on Power Systems, 1995, 10(4): 1781–1788

[21] Krause P C, Wasynczuk O, Sudhoff S D. Analysis of Electric Machinery and Drive Systems, 2nd Ed. New York: John Wiley & Sons Publication-IEEE Press, 2002

[22] Lehtla T. Parameter identification of an induction motor using fuzzy logic controller. 2014-10

[23] Sandhu K S, Pahwa V. A novel approach to incorporate the main flux saturation effect in a three-phase induction machine during motoring and plugging. International Journal of Computer and Electrical Engineering, 2011, 3(3): 443–448

[24] Kundur P. Power System Stability and Control. India: Tata-Mcgraw-Hill, 2006

[25] Boldea I, Nasar S A. The Induction Machine Handbook. New York, USA: CRC Press LLC, 2001, Chapter 13

[26] Wang K, Chiasson J, Bodson M, Tolbert L M. A nonlinear least-squares approach for identification of the induction motor parameters. In: Proceedings of the 43rd IEEE Conference on Decision and Control. Atlantis, Paradise Island, Bahamas, 2004, 14–17

[27] Pedra J, Sainz L, Corcoles F. Study of aggregate models for squirrel-cage induction motors. IEEE Transactions on Power Systems, 2005, 20(3): 1519–1527

RIGHTS & PERMISSIONS

Higher Education Press and Springer-Verlag Berlin Heidelberg