Modeling and Semi-Analytic Stability Analysis for Dynamics of AC Machines

: In this article, a semi-analytical technique is proposed to predict stable sustained periodic responses of AC electrical machines. Based on such desired outputs, the proper selections of machine variables are captured, such as the perturbation parameter arisen from the relative movement between the stationary and rotating parts. Compared to the experimental results, the derived analytical results are relatively well-ﬁtted with the studied practical cases.


Introduction
In electrical engineering, the development of electrical systems, such as AC electrical machines, some reliable mathematical models are needed that take into account their theory. Such a model should cover the different views on an interested system as well as the different degrees of detailing that leads from a hierarchy of model to a hierarchy of design system parameters. The specific views of such a system represent the relevant effects of specific interests, such as dynamics, stability, and control, which are used as interfaces for the design engineers during their works. Thus, a validation of the mathematical model is needed, then, a solution response with the experimental system results over the entire range of system interests or for the considered range of switching system parameters have to be compared, cf. [1][2][3].
The solving process for validation can be managed mathematically by various analytical or numerical approaches. In fact, numerical simulations are often necessary in order to understand the attitude response and control characteristics of complex systems; however, the typical problems in numerical simulations are the short time visualization and error analysis, since it may be prohibitively expensive and time consuming, especially when large numbers of parameters are involved. In contrast, analytical descriptions can be of great help in obtaining a qualitative understanding of the complex dynamical behavior, even simple heuristic analytical results may provide a fast and relatively accurate form for maneuver analysis, cf. [4][5][6][7]. For instance, in most engineering control applications, it appears that the trend for autonomous engineering control schemes dictate the use of compact, simple, and analytical expression responses for evolutionary systems, cf. [8][9][10].
Typically, semi-analytical approaches and perturbation methods are used to solve such kinds of models. Recently, some new and emerging approaches were proposed to handle both nonlinear and linear systems with varying coefficients or to consider large parameters and domains of definition as well. Among these approaches that are employed for solving complex systems are the perturbation-iteration method (PIM) , variational iteration method (VIM), homotopy analysis method (HAM), and the homotopy perturbation method (HPM), cf. [11][12][13][14][15]. Moreover, there are some improved homotopic approaches, such as the optimal homotopy perturbation method (OHPM) and the optimal homotopy asymptotic method (OHAM); it is worth mentioning the modified differential transform method (MDTM), which is particularly settled in solving boundary values problems, cf. [16][17][18][19]. Most of such approaches have been successfully applied in solving numerous linear and nonlinear problems in science and engineering.
In industrial electrical induction/synchronous machines, mathematical formulations have been carried out in the simulation field, including some internal generated effects such as the rotor temperature, the relative movement or the magnetic saturation with a view to obtaining highly accurate results of output or control, cf. [20][21][22][23]. As often as to simplify the electrical construction of these machines, electrical circuits are rigorous simple models used to exhibit the behavior of their outputs. Indeed, it can be considered as very efficient alternative tools to express massive laboratory experiments for their practical purposes, cf. [24][25][26][27][28].
In this work, the simplification of electrical construction of induction/synchronous machines taking into account the relative movement between the stationary and rotating parts using an RLC circuit modeling is considered. In this case, the electrical circuit system has allowed us simply to derive the governing equation containing the effect of relative movement by using varying inductance expression. At that time, using the ease of handling by circuit modeling, semi-analytical expression for the machine response and the configuration stability are captured.
The use of analytical expression via one of the mentioned approaches is limited to isolate some cases in which the natural response of the governing equation can be expressed explicitly in terms of functions in time, or, at least, it can be reduced to cases of quadratures. The variational iteration technique was found a very efficient tool to form the governing equation to its quadrature form, and the solution might be expressed by terms of iterated integrations, cf. [29][30][31]. The obtained solution by this method captures all values of circuit parameters with large domains of definitions. As far as in specific domains of definitions for some circuit parameters, some special approximate expressions of the solution in terms of time can be initiated by using special methods such as Frobenius method and Poincare-Linstedt approach, which are constructed on specific domains of some parameters, cf. [32][33][34][35].
The benefit of the isolated possible cases, obtained from (semi)analytical solutions shared the property of having them as meromorphic functions of time, might be a very efficient indicator for the stability of the modeled machines. Hence, it turned out that necessary conditions for the sustained periodic motions of the electrical machines ought to be explicitly recognized.
In the present paper, we report on the analytical investigation of the AC machine response using a prototype modeling via modified variational iteration methods (VIM) in terms of iterated meromorphic functions of time. Our contribution is concentrated on the mathematical aspects of that method to obtain the stability domains as well as the transition curves of the governing equation for a wide range of values for the related system parameters. More concretely, our attention is focused on the existence of sustained periodic oscillations represented by the boundary curves to describe stable periodic natural responses using the energy-rate function.
The paper is divided into three sections. In Section 2, the governing system is obtained using the circuital RLC series modeling. In Section 3, the construction of semi-analytical solutions are given via the variational iteration method. In Section 4, stability charts using the concluded resulted from the previous sections are drawn via energy-rate approach. In the last section, the conclusion is given.

The Governing Dynamic Model
Typically, the basic operation of many types of electrical induction/synchronous machines relies on periodic variation of inductance in time. Self-excitation of induction generators can be achieved when connecting machine terminals to sufficient amount of capacitance that depends on the driving rotor speed and load conditions, cf. [36]. The three-phase machines have been successfully analyzed using well-known approaches that depend on d-q transformation, cf. [37]. However, single-phase operation is difficult to analyze using conventional methods. A single-phase parametric generator is modeled as an RLC circuit with periodically time varying inductance, cf. [38]. This type of generator is expected to be suitable for energy harvesting applications, e.g., electrical generation from wave energy. Practically, the inductance variation can be attained through variation of magnetic coupling between two series-connected coils; one is fixed while the other is rotating as shown in Figure 1. From the basic laws of electrical circuits and Figure 2, the voltage equation of the RLC-series circuit is given by the following equation, where t, i(t), V s (t) are the time, the electric current, and the supply, respectively. The coefficients L t (i, t), R t (i, t) and C t (i, t) are the inductance, the resistance and capacitance functions in time, and the electric current. By considering that the resistance and the capacitance are constant and the inductance varies periodically with time due to the existence of relative movement between the stator and the rotor of such AC machine, some loop inductances of AC machines depend upon the rotor position. Hence, the governing equation takes the following form regardless of the power measurement and control accessories. For simplicity, let us define the following new variables in terms of the angle x = ωt in the natural response of the system, where q and ω are the charge and the angular frequency respectively. Let us introduce the time varying inductance as, where L 0 is the reference value of the inductance and F(x) is considered as predictable function for the studied AC machine.
Including the new variables in the homogeneous part, we obtain Integrating both sides once, taking into the account the natural response only, hence we get Hence, dividing both sides by L 0 , it yields By considering that If the predictable function of inductance is considered periodically time dependent as follows where h is the perturbation parameter arisen from the dynamic perturbation on the inductance due the relative movement between stator and rotor of the machine. Now, the governing equation of the modeling circuit described by a linear model as follows (1 + h cos 2x)y + Q α y + 1 Hence, the governing equation of the natural response reads where c 1 and c 2 are real constants representing the values of variants y and its derivative at the initial angular position of the rotor oscillatory movement. Indeed, the solution (y(x)) of Equation (11) represents the net positive charge that moved past an instantaneous point in a specified direction in the shown circuit under the specified initial conditions and the absence of the voltage supply (V s ), cf. [5].

Construction of Semi-Analytical Solutions via VIM
In this section, the explicit algebraic solutions as meromorphic functions of time and the rough analytical analysis of stability are given. Particularly, in accordance with the existence of large domains of perturbation parameter (h), approximate forms of the solution are needed to construct. Consequently, the curves containing a stable periodic characterization might be drawn.
Variational iteration method is needed as a general method to capture the behavior of solutions for a large scale of h. This method is an efficient technique to construct general solutions and not influenced by inherent variation of problem parameters, cf. [39]. Consequently, the yielded solutions describe the behavior of the governing equation in neighborhoods of boundary curves and upon them. Consider the following differential equation (DE) where L and N are linear and nonlinear operators, respectively, and f (x) is the source of inhomogeneous term. The method admits the use of a correction functional for Equation (12) in the form where λ is a general Lagrange's multiplier, which can be identified optimally via variational theory,ỹ(t) is a restricted variational such that δỹ(t) = 0. Lagrange's multiplier is very critical in this method and it can be constant or function according to equation and conditions of it, the zeroth approximation y 0 can be any function, using the initial values y(0), y (0) and y (0) are preferable used for the selective zeroth approximation y(0). The final solution is given by Now the general second order differential equation can be written as Following [14,15,31,40,41], then VIM admits the use of the correction functional for this equation by y n+1 (x) = y n (x) + x 0 λ(t)(y n (t) + a(t)ỹ n (t), +b(t)ỹ n (t))dt, where y(0) = c 1 , y (0) = c 2 .
Taking the variation of both sides of Equation (16) with respect to independent variable y n (x), we find δy n+1 δy n = 1 + δ δy n x 0 λ(t)(y n (t) + a(t)ỹ n (t) + b(t)ỹ n (t))dt, (17) or equivalently Using δỹ n (t) = 0, δỹ n (t) = 0 that gives Integrating Equation (19) by parts twice gives This is equivalent to δy n+1 (x) = 0. This means that the left side of Equation (21) is 0, and the right side of Equation (21) should be 0 as well. Thus, this yields the following stationary conditions This in turn gives Substituting this values of λ into Equation (16) gives the iteration formula such that y(0) = c 1 , y (0) = c 2 and We find the solution of the system Equation (15) by using the iteration formula Equation (26) y n+1 (x) = y n (x) + x 0 (t − x)(y n (t) + Q α(1 + h cos 2t) y n (t) where y 0 (x) = y 0 where Li s is the polylogarithm function and is defined by a power series in z, which is also a Dirichlet series in s: Thus, to overcome the difficulty of integration term, we can use the transformation then From the first approximation (Equation (33)), it seems very predictable for the aperiodic solutions, but for the periodic ones, the treatment should involve them. In the case of |h| < 1, we can use the transformation then y 1 (x) =y 0 − y 0 36α 2 (9(h 2 + 2)x 2 + 9h(2h 2 − h + 2) sin 2 (x) − 3h 2 (−4h + 3) sin 4 (x) − 8h 3 sin 6 (x)). (35) For further simplification, it reads y 1 (x) =y 0 + y 0 288α 2 (h 2 (56h + 9) + 72h − 72(h 2 + 2)x 2 − 18h(3h 2 + 4) cos 2x + 9h 2 cos 4x − 2h 3 cos 6x) + .... (36) The second approximation (n = 1) can be obtained as where Regarding the case of |h| < 1, the modified VIM method (MVIM) using the linear function of Lagrange multiplier can be modified to obtain a general periodic characterization. We should search for Lagrange' multiplier to maintain this property. To overcome this difficulty, we can use the following transformation represented by Then, the governing equation is converted to the following y + 1 Following [29,40,41], the modified VIM admits the use of the correction functional by where y(0) = c 1 , y (0) = c 2 Taking the variation of both sides of Equation (42) with respect to independent variable y n (x), we find δy n+1 δy n = 1 + δ δy n ( x 0 λ(t)(y n (t) + 1 α 2 y n (t) + a(t)ỹ n (t) + b(t)ỹ n (t)dt), (43) or equivalently Using δỹ n (t) = 0, δỹ n (t) = 0 gives Integrating Equation (45) by parts twice to give This is equivalent to δy n+1 (x) = 0, this mean that the left side of Equation (47) is 0 and so the right side of the equation should be 0 as well. This yields the following stationary conditions λ + 1 This in turn gives Substituting this value of λ into Equation (42) to give the iteration formula.
such that y(0) = c 1 , y (0) = c 2 and y 0 ( where y 0 (x) = y 0 . Hence, The verification of the obtained semi-analytical solutions with the corresponding numerical solutions are vividly shown in Figure 3. The asymptotic characterization of the response of Equation (11) is strongly verified numerically using the 4th order Rung-Kutta algorithm with the shown initial values. As it is clearly shown in the comparison, the semi-analytical solution is comparable with the numerical one at different values of the perturbation parameter. All results exhibit mostly that the stability characterization achieved merely in the domain of definition: |h| < 1, which is coinciding with the results in [5,38].

Construction of Stability Domains
Several techniques are used to predict the stability of systems, cf. [42][43][44][45][46][47]. Based on the constructed semi-analytical solutions, it is readily to apply the energy-rate method to obtain three dimensional illustration of stability surface showing a sense of relative stability of parameter plane, cf. [48][49][50]. Based on the shown figure, Figure 3, the stability regions lie on the domain |h| < 1. So that the intersection of a zero plane with stability surface coming from the solution in |h| < 1 determines the boundary curves of the system. The energy-rate function for our system reads where E is the mechanical energy of the system, and represents the instantaneous rate of generated or absorbed energy by the applied force Then the energy-rate function can be calculated by The energy-rate function is zero for free vibration of a conservative system, as well as the steady state T-periodic response of the system. Equation (57) includes three parameters, h, Q, and α; therefore, the parameter space is two-dimensional if Q is fixed. Determination of the stability domains is carried out by finding the pairs of (α, h) located on the boundary of stable domains for every fixed value of Q. Plugging the semi-analytical solution represented by Equation (53) and Equation (54) for the case of |h| < 1 in Equation (57) and by using the iterative symbolic computations of Maple program, then we get the value of the energy-rate (Γ) as function of h, α and Q(omitted here for brevity due to its long form).
Analogous to the mechanical system if Γ > 0, then the energy is being inserted to the system and the point (α, h) belongs to an unstable domain. In contrast, if Γ < 0, then the energy is being extracted from the system and the point (α, h) belongs to a stable domain. In the case of Γ = 0, common boundaries of these two domains are appeared called transition curves. So that, in the case of zero energy-rate, we look for the transition curves by the relation between α and h at fixed values of Q with the period T = π or T = 2π according to Floquet's theory, cf. [5]. The resulting figure at Γ in the domain |h| < 1 is shown in Figure 4 starting with the value of (0, 1). These results have relative agreement with the experimental work on AC synchronous machines shown in [26,38,51].

Conclusions
In this work, the modeling of AC machines has been dealt with, taking into account the existence of the effect of relative movements between the stators and rotors, such as the RLC electric circuits with periodically time varying inductance. Semi-analytic solutions are obtained for the circuit responses via the modified variational iteration method. In that case of study, the resulting solutions are likely to be valid for large domains of definitions of parameters in the governing system. Within the domain |h| < 1, the stability domains are recognized using the zero surface of energy-rate function to represent the transition curves between the stable and unstable regions. This result was shown on a wide variety of practical and experimental literature of the studied problem. Moreover, the relationships among the circuit parameters for stable sustained periodic outputs are configured. Funding: This research did not receive specific funding, but was performed as part of the employment of the authors: University of Tanta and University of Kafrelsheikh (Faculty of Engineering), Egypt.