Power System Model Identification via the Thermodynamic Formalism

Proceedings of the 18th World Congress The International Federation of Automatic Control Milano (Italy) August 28 - September 2, 2011

Power System Model Identification via the Thermodynamic Formalism Richard M. Kolacinski Charles Stark Draper Laboratory, Cambridge, MA 02134 USA (e-mail: [email protected]). Abstract: This paper focuses upon the development of a methodology for data–driven construction of mesoscopic models of the T&D system for use in real-time monitoring and control. The system dynamics are lifted to a discrete covering space which provides an encoding of the system dynamics within symbol strings. These symbol strings are treated as Bernoulli shifts and are characterized, via the machinery of information theory and formal language theory, as probabilistic automata. As these automata are fundamentally pattern recognizers, they provide a fundamental basis for event/anomaly detection and thus a basis for critical grid monitoring functions such as security state identification. Keywords: Identification, Power systems security, Power systems stability, Statistical inference 1. INTRODUCTION

tion of the system equilibria. The present work extends these notions, interpreting the jump processes used to model the effect of disturbances as modeling transitions in the system’s trajectory relative to the various domains of attraction.

The effective maintenance and reliable operation of T&D systems is tremendously challenging due to their scale, nonlinearity, stochasticity, structural complexity, and coupling between continuous (e.g., voltage, angle) and discrete (i.e., network topology) states. Furthermore, real time operation presupposes the ability to evaluate the dynamic security of the system and to develop timely control strategies in response to evolving conditions, where by security is meant an instantaneous time varying measure of the robustness of the system relative to imminent disturbances, see Fink and Carlson (1978). The concept of security was introduced first in DyLiacco (1968) wherein steady-state security is defined within the context of a multilevel decomposition of the power system and in terms of satisfying a set of equality and inequality constraints over the “next contingency set” of potential disturbances.

A critical component of grid operation is the proper identification of its security state (i.e., normal, alert, emergency, in-extremis, and restorative) as the selection of state estimation filters, contingency sets, and control objectives and policies are dependent upon security state. The current work characterizes security states by the behaviors of the system trajectories. The crucial concept motivating this approach is that system trajectories will possess significant qualitative differences based operant equilibria and, hence, can be used to identify the basin of attraction within which the system trajectory lies. This interpretation of security state as behavior relative to equilibria strongly suggests the application of a thermodynamic formalism where changes in operant equilibria can be viewed as analogs of phase transitions. This paper develops the methodology for characterizing the qualitative behaviors associated with trajectories in the disparate regions of the state space within a thermodynamic framework and, by doing so, provides the basis for identifying the security state of a power system. Specifically, techniques based upon ε–machine reconstruction introduced in Crutchfield and Young (1989) are used to characterize the behavior of disparate trajectories using probabilistic computational automata. The efficacy of the method is then demonstrated via simulation of a sample power system.

This formulation is a pointwise approach; security is assessed for a given set of operating points and thus requires the specification of operating points of interest as well as their associated contingency set that must, in practice, be restricted to a tractable number of “credible” contingencies. Furthermore, this method relates only to steady-state analysis in that only initial disturbances are examined (i.e., not cascading events) and that it assumes that the steady-state equilibrium (operating) point is reachable in the post-fault system. The challenges posed by the dynamics and stochasticity of the T&D system and these new technologies, however, require a more robust framework. DiLiacco’s framework and its extensions provide the basis for Loparo and Abdel-Malek (1990) which presents a concept for a dynamic security region for all possible system topologies, modeling disturbances as finite state jump processes, and introduces a probabilistic security measure that evaluates the probability that the system is in the dynamic security region. Wu and Tsai (1983) characterize dynamic security regions in terms of the domains of attrac978-3-902661-93-7/11/$20.00 © 2011 IFAC

2. PRELIMINARIES The focus of this paper is the characterization of complex behaviors associated with power systems. In particular, an approach for constructing mesoscopic descriptions and models using the thermodynamics formalism is presented. The thermodynamic formalism refers to the application 507

10.3182/20110828-6-IT-1002.03168

18th IFAC World Congress (IFAC'11) Milano (Italy) August 28 - September 2, 2011

of techniques analogous to those used in conventional thermostatistics to the study of nonlinear systems. Its principle components are ensembles, which play a role analogous to sample spaces, and free energies, typically characterized via R´enyi entropies. To this end, this section introduces the mathematical underpinnings of the probabilistic models used to characterize the system behavior and the fundamental theory necessary for their inference. 2.1 Markov Partitions and Symbolic Dynamics In the following development, we make the following assumptions: (A1): Observable phenomena (e.g., phase measurements) can be modeled with sufficient accuracy as an ordinary differential equation (ODE) of the form, x˙ = f (x) (1) where x ∈ Rn is the state vector and f : Rn 7→ Rn . (A2): The system (1) evolves on a compact manifold.

Fig. 1. Refinement of Partitions of B. Letting si ∈ Σ be the index of the domain B ∈ B visited at time i, the itinerary O is identified by the symbolic sequence S = {s0 , s1 , . . . , sk }, where xi ∈ Bsi ∀i = 0, 1, . . . , n. Note that the specification of a sequence implies a restriction T Tof the initial T to a subdomain cell defined by BS ≡ BS0 F−1 (BS1 ) · · · F−k (BSk ), as depicted in Fig. 1. Note that this intersection is not necessarily smaller than B0 , but if the intersection is smaller it provides a refinement.

The vector field f generates a flow Φt (x, u) : X 7→ X which represents the image x(t) of the initial condition x = x(0) ∈ X at time t. The set X is invariant for a flow Φt if Φt (x) ∈ X for any x ∈ X and for all t ∈ R. The pair (f , X) constitutes a dynamical system where the f is the vector field specified in equation (1) and set X is invariant for a flow Φt . The techniques of symbolic dynamics are used to lift system dynamics to a discrete covering space. Within this discrete space, system trajectories are described by a string of symbols wherein the system dynamics are described by the shift dynamics of the string. Substrings provide the basis for ensembles and, hence, the encoding of the dynamics into a symbol string provides entry into the thermodynamic formalism and permits the machinery of information theory and formal language/automata theory to be brought to bear. Symbolic encoding requires a discretization in time that is accomplished via Poincar´e maps. Definition 1. Let Ξ ∈ Rn be an (n − 1)-dimensional surface defined such that Φt (x) ∈ Ξ for some t > 0 and any initial condition x ∈ X, and such that for any x ∈ X, f (x) is not tangent to Ξ. A Poincar´e map F(x) : X 7→ X is given by xk+1 = F(xk ) (2) where F(xk ) = Φtk (xk ) and tk is the elapsed time between successive intersections of the system trajectory xk and xk+1 with Ξ. If Φt is a smooth function, F is a diffeomorphism. Fk denotes the k th composition of F with itself. Negative kvalues represent iterates of the inverse map. The symbolic encoding is accomplished by introducing a partition of the Poincar´e surface and associating a label with each cell of the partition and the set Σ = {0, . . . , b − 1} of labels is known as the alphabet. Definition 2. A topological partition of a metric space X is a finite collection B = {B0 , . . . , Bb−1 } of disjoint T ¯j open sets, i.e., Bi Bj = ∅ ∀i 6= j, whose closures B Sb−1 ¯ together cover X in the sense that j=0 Bj = X.

The first refinement B1 of B under F consists of the T subsets Bsj F−1 (Bsk ) for all sj , sk ∈ Σ where the intersection is nonempty. In general k _ F−i B = B ∨ F−1 B ∨ · · · ∨ F−k B (3) Bn = i=0

−i −i −i where T F B = {F (B0 ), . . . , F (Bb−1 )} and B ∨ C = {Bj Ck : 0 ≤ j ≤ b − 1, 0 ≤ k ≤ b − 1}. The closures ¯ n ∈ Bn of these sets are compact and set inclusion is a B partial ordering, hence the sequence of refinements T∞ ¯ n Bn is a lattice, called the partition lattice and n=0 B 6= ∅. Definition 3. Let (F, X) be an invertible discrete time dynamical system. A topological partition B of X gives a symbolic representation T∞ ¯ n S of (F, X) if for every x ∈ S consists of exactly one point. the intersection n=0 B Such a partition is called generating.

In other words, the symbol sequence associated with a generating partition is a faithful (i.e., unique) encoding of the real trajectory O. This implies an equivalent sequence space. Let Σb = ΣZ indicate the (power) set of all twosided infinite sequences S = {. . . , s−1 , s0 , s1 , . . .} over the alphabet Σ. Definition 4. The (left) shift map σ ˆ : Σb 7→ Σb is defined πi (ˆ σ (S )) = πi+1 (S ) (4) where πi : Σb 7→ Σ, such that si = πi (S ) is the ith projection of the power space Σb onto the base space A. σ ˆ provides the formal mechanism for relating the iterates of the map described by equation (2) to S . The pair (Σb , σ ˆ ) constitute a dynamical system called the full shift (or Bernoulli shift) on b symbols. Shifting S is equivalent to iterating the map F. 2.2 Formal Languages and Information

Under the action of the dynamics, the system (1) describes an orbit O = {x0 , x1 , . . . , xk } which visit the elements

In general, a physical system or a mathematical model does not allow all concatenations of symbols to occur, 508

18th IFAC World Congress (IFAC'11) Milano (Italy) August 28 - September 2, 2011

recurrence of substrings in regular languages. By the Pumping Lemma (for RL), for every regular language there is a string length p such that all strings of this length or greater must represent a path on the associated finite state machine containing a closed loop. This leads directly to the notion of inspecting substrings of increasing lengths for the existence of cycles. Length D subsequences sD = {si . . . si+D−1 : sj = (s)j } of the data stream s are referred to as D-cylinders and the model inference process thus begins with construction of histograms of observed D-cylinders.

hence shift dynamical systems can be classified according to the properties of the sequences that belong to their invariant sets. This can be done in a systematic way by examining all finite subsequences occurring in S . Formal language theory provides the machinery for this. Definition 5. A language L is a set of words (concatenations of symbols) formed from an alphabet Σ. Every formal language can be associated with a discrete automaton that is able to accept (or recognize) all legal words of L when acting on a signal S . In the present work, the focus is restricted to regular languages (RL)/finite state automata (FSA). Definition 6. A finite state automaton (FSA) is a 5–tulple M = (Q, Σ, δ, q0 , F ), where Q is a finite set of states, Σ is a finite input alphabet, q0 ∈ Q is the initial state, F ⊆ Q is the set of final states, and δ : Q × Σ 7→ Q is the transition function.

The natural structure for representing histograms of Dcylinders of increasing length up to some length L is a prefix tree of depth L, denoted T L . Probabilistic structure is added to the tree by recording for each node the number of occurrences N (ωn ) of the associated D-cylinder ωn relative to the total number N (D) of D-cylinders of length l = |ωn |. The estimate of the node probabilities are provided by the frequency estimator pTn = NN(ω(l)n ) . In this context, Shannon entropy measures on the observed sequences of length L, H(X1 , X2 , . . . , XL ), is equivalent to a Shannon entropy measure on the leaf nodes of T L , denoted H(T L ).

There are several extensions to the finite state automata given above. One generalization that is particularly appropriate for the construction of a generative model extends the output from a binary “accept/don’t accept” to an arbitrary alphabet. Definition 7. (Mealy machine). A Mealy machine is a 6– tulple M = (Q, Σ, ∆, δ, λ, q0 ), where Q, Σ, δ, and q0 are as defined for the finite state automata above. ∆ is the output alphabet and λ : Q × Σ 7→ ∆.

If the string has a simple structure (e.g., cyclic), the tree representation of the string will converge at some finite depth L. That is, increases in tree depth will not provide any additional information on the strings content and hence H(XL+1 |X1 , X2 , . . . , XL ) = 0. Applying the chain rule for entropy, we obtain

The descriptive power of language, hence the processing power of the associated automata, are functions of how much storage is required to capture a long sequence. In a probabilistic setting, this measure of processing power can be interpreted as the amount of information contained in a long sequence. The classic measure of self-information of a random variable is given by Shannon entropy. In the context of languages, it is natural to ask how the entropy of a sequence grows with its length. This is given by the entropy rate, the Kolmogorov-Sinai entropy which is the R´enyi entropy of order 1. Definition 8. The entropy rate of a stochastic process {Xi } is defined by 1 h(Σ) = lim H(X1 , X2 , . . . , Xn ) (5) n→∞ n when the limit exists.

H(XL+1 |X1 , X2 , . . . , XL ) = H(X1 , X2 , . . . , XL+1 ) −

L X

H(Xm |X1 , X2 , . . . , Xm−1 )

m=1

(6)

= H(X1 , X2 , . . . , XL+1 ) − H(X1 , X2 , . . . , XL ) = H(T L+1 ) − H(T L ) Thus, if the tree converges at a finite depth, the entropy difference over successive tree depths will vanish. Clearly, equation (6) makes the relationship between prefix tree convergence and the entropy rate on the tree explicit: Definition 9. Let T L be a prefix tree of depth L associated with a symbolic sequence over a fixed alphabet Σ. The entropy growth rate of T is given by  1 H(T L ) − H(T L−1 ) (7) h(T ) = lim L→∞ L We say that the tree representation is convergent in the limit if h(t) = A < ∞.

3. MAIN RESULTS Just as equations of motion and constraints restrict the set of admissible trajectories to lie on a manifold, the “legal” symbol strings are restricted to a subset of all possible strings. This subset is a language and, hence, the modeling problem can be cast as inferring the membership rules to this language or its grammar. In this section, ε– machine inference techniques presented in Crutchfield and Young (1989) are formalized and basic results concerning convergence are given.

This suggests a metric for monitoring tree convergence for finite sequences, namely the approximation   ˆ L ) = 1 H(T L ) − H(T L−1 ) . h(T (8) L Lemma 10. Given an infinite symbol sequence S over a fixed alphabet Σ, the associated tree representation is convergent in the limit.

The language inference techniques described here are predicated upon a priori knowledge of the language’s structure. That is to say, that the different languages possess different characteristics which may be exploited to deduce their grammars, or more accurately, the automata which recognize them. Of particular importance is the

Proof. Let S = σ1 σ2 · · · · with σi ∈ Σ, i ∈ Z+ and let S L = X1 X2 · · · XL with Xi ∈ Σ, i ∈ {1, 2, . . . , L} be a finite substring (D–cylinder) extracted from S with T L the associated L-depth prefix tree. By the chain rule for entropy, equation (5) may be expressed 509

18th IFAC World Congress (IFAC'11) Milano (Italy) August 28 - September 2, 2011 n 1 X H(Xm |X1 , X2 , . . . , Xm−1 ). n→∞ n m=1

Corollary 11 suggests a slightly simpler metric for monitoring tree convergence, namely the source entropy

H(Σ) = lim

∆H(T L ) = H(T L ) − H(T L−1 ).

Since H(·|·) ≥ 0 by definition and m ≤ n

If the information in the tree structure grows without bound, that is, if the entropy growth rate does not vanish, the reconstruction process must proceed to the next level in representations and explicitly construct FSA associated with the input string. The template of machinery necessary to construct this representation is provided by the Myhill–Nerode Theorem which asserts the existence of a finite partition of the set of legal strings into equivalences classes and associates each class with a state of the automaton. Moreover, as the automaton’s states and the equivalence classes are isomorphic, the automaton’s state can be seen as indexing the equivalence class of the part of the string seen by the automaton.

n 1 X H(Xm |X1 , X2 , . . . , Xm−1 ) n m=1

= ≤

n X

m=1 n X

n−1 H(Xm |X1 , X2 , . . . , Xm−1 ) (9) m

−1

H(Xm |X1 , X2 , . . . , Xm−1 )

m=1

Thus, if the limit exists, n 1 X H(Xm |X1 , X2 , . . . , Xm−1 ) n→∞ n m=1

lim

n X

(18)

(11)

In order to infer probabilistic automata, the machinery of the Myhill–Nerode theorem must be extended to operate on sets of strings and their associated histograms, or, equivalently, on subtrees of the prefix tree. Definition 12. (Subtree similarity). A pair of D-level subtrees, TnD1 and TnD2 , are said to be similar, denoted TnD1 ∼ TnD1 if their branch topologies are identical and the associated branching probabilities are equal. Lemma 13. Subtree similarity is an equivalence relation on subtrees.

Since the series (11) dominates the entropy rate, the entropy rate is defined if the series (11) is convergent. By d’Alembert’s Ratio test, the series is absolutely convergent if n − 1 H(Xn |X1 , X2 , . . . , Xn−1 ) < 1 (12) lim sup n→∞ n H(Xn−1 |X1 , X2 , . . . , Xn−2 )

Proof. For each D-level subtree there are k = |Σ|D possible dependent nodes. We associate with each D-level subtree and ordered k-tuple tD ni where the jth element D of tn1 lies in the closed interval [0, 1] and corresponds to the branching probability of the jth node depending from node ni . Subtree similarity can thus be seen as vector equality which has the reflexive, symmetric, and transitive properties. 2

≤ lim

n→∞

m−1 H(Xm |X1 , X2 , . . . , Xm−1 ).

(10)

m=1

The right hand side of equation (10) is the limit of partial sums Sn of an infinite series, thus, lim

n→∞

n X

m−1 H(Xm |X1 , X2 , . . . , Xm−1 )

m=1

= lim Sn = n→∞

∞ X

m−1 H(Xm |X1 , X2 , . . . , Xm−1 )

m=1

Applying Definition 9 to equation (12) yields lim max

L→∞

h(T L ) <1 h(T L−1 )

Denoting the set of all D-level subtrees by T D , the equivalence relation ∼ induces a partition on T D , the cells of which are the elements of the quotient set T D / ∼ = {[TnDi ]|TnDi ∈ T D } where [TnDi ] = {TnDj |TnDi ∼ TnDj } is the equivalence class containing TnDi . We will denote [TnDi ] by Ci and T D / ∼ by QT . The set of D-level subtrees in Ci is isomorphic to the set of cylinders ωni preceding each TnDj ∈ Ci . Thus Ci is also an equivalence class of Dcylinders and we say ωnj ∈ Ci if ωnj precedes TnDj ∈ Ci .

(13)

H(T L ) takes on its maximum value if all L-cylinders are possible and equally likely, L 1 log |Σ|−L = − L log |Σ| (14) Hmax (T L ) = L |Σ| |Σ| Direct substitution of equation (14) into equation (7) and algebraic reduction produces log |Σ| hmax (T L ) = [(L − 1)|Σ| − L] (15) L|Σ|L

Definition 14. We define on QT a probabilistic next state function P (Ci , Cj , sα ) = Pijα where we define Pijα = Pr{ωn |ωn ∈ Ci , ωm ∈ Cj , and ωn = ωm sα }, i.e., the probability of transitioning from a D-level subtree in Ci to a D-level subtree in Cj on a symbol sα .

The ratio test may then be expressed (|Σ| − 1)L2 + (1 − 2|Σ|)L + 1 |Σ| − 1 lim max = <1 L→∞ |Σ|2 L2 + (|Σ| − 2|Σ|2 )L |Σ|2 (16) 2 Corollary 11. Given an infinite symbol sequence S over a fixed alphabet Σ,   (17) lim H(T L ) − H(T L−1 ) = 0.

The next state functions also define an next symbol function as the probability associated with each legal output symbol. In fact, the next state function P (Ci , Cj , sα ) can be thought of as a collection of next state functions for a set of automata, each associated with a particular member of Σ, running simultaneously a la nondeterministic finite automata. This collection of next state functions can be collapsed into a single P function for each equivalence class P (Ci , Cj ) = Pij = α∈|Σ| Pijα .

L→∞

Proof. This is a direct consequence of the nth term test for divergence. 2

Definition 15. The probabilistic state vector is defined

510

18th IFAC World Congress (IFAC'11) Milano (Italy) August 28 - September 2, 2011

of the representation’s suitability can be constructed that is analogous to the use of source entropy to evaluate the suitability of a tree representation. Definition 17. For a Finite State Automata, the source entropy is given by the finitary complexity: X Cµ (L) = − pi log2 pi (24)

(19) pi = Pr{ωn ∈ [TnDi ]} The state transition functions on symbol sα are then defined pt+1 = Tijsα ptj with the associated state transition i matrices defined (20) Tijsα = Pijα , The probabilistic state connection matrix is defined X (21) Tijsα , Tij =

Ci ∈QT

where the probabilities pi are those associated with the equivalence classes Ci ∈ QT .

sα ∈Σ

with the associated total state transition function pt+1 = i Tij ptj .

The finitary complexity is the total Shannon information contained in the subtrees associated with equivalence classes. As in the case of tree machines, these metrics can be used to determine when the FSA representation is inappropriate. Again, the convergence of the total information contained in the machine, i.e., the growth rate of the total information provides this indicator. Definition 18. The growth rate of total information in an automata is given by the graph entropy rate, i h (25) cµ (L) = lim 2Cµ (L) − 2Cµ (L−1) .

The input stream of interest is assumed to have a stochastic component and hence, a “legal” string cannot be assumed to have a determinate starting symbol nor to emanate from a start state of the signal source. Any model of this signal must be capable of addressing this uncertainty. The addition of -moves provides the machinery for addressing this. That is, we assume that the underlying is capable of changing state while emitting the empty string, . The automata models to be inferred are thus assumed to begin in some neutral state that stands in lieu of the actual start state and, with the first symbol generated, transitions to an appropriate state on . To this end, we augment the quotient space with a neutral class C0 = ∅, the alphabet with , i.e., Σ∗ , and the mappings P (Ci , Cj ) and P (Ci , Cj , sα ) are extended to accommodate the augmented quotient space and alphabet. Theorem 16. The minimal probabilistic finite automata that accepts S is the quintuple M = (Q∗T , Q∗T , P ∗ (Ci , Cj ), C0 , Q∗T ), (22) ∗ ∗ ∗ ∗ where QT = QT ∪ ∅, and P (Ci , Cj ) : QT × QT 7→ Q∗T . Moreover, the minimal Mealy machine capable of generating S is the sextuple M = (Q∗T , Q∗T , Σ, P ∗ (Ci , Cj ), P ∗ (Ci , Cj , sα ), C0 ), (23) where P ∗ (Ci , Cj , sα ) : Q∗T × Q∗T × Σ 7→ Q∗T .

L→∞

Theorem 19. If the observed symbol stream comes from a regular language, the sequence of finite automata M produced with increasing reconstruction length is convergent, and cµ = 0 Proof. By the minimality of the construction process given by Theorem 16, the minimal automata corresponding to T L is inferred at each step in the reconstruction process. By the pumping lemma at some point the reconstruction length L will become sufficiently large to capture all possible cycles in the automaton and subsequent increases in the reconstruction length will produce no additional states and the sequence of finite automata will converge. A straight-forward restatement of Theorem 10, mutatis mutandis, will show that if the sequence of finite automata is convergent, then cµ = 0. 2 Thus as representations built on longer and longer sequences subsequences fail to converge, the inappropriateness of the representations can be determined and an increase in representation complexity can be effected.

Proof. By definition, Q∗T is closed under P ∗ (Ci , Cj ). Furthermore, by Definition 14 and Definition 15, it can be seen that P ∗ (Ci , Cj ) is a linear transformation with Tij , defined in equation (21), as the matrix of the transformation. Thus successive applications of P ∗ (Ci , Cj ) are equivalent to a matrix multiplication and hence P ∗ (Ci , Cj ) is associative and (Q∗T , P ∗ (Ci , Cj )) is a semigroup.

4. APPLICATION This section introduces a simplified model of power systems that captures behaviors of interest and, using this model, demonstrates the efficacy of the model identification techniques developed in the prequel.

With every semigroup A is associated an automaton, M = M (A) known as the automaton of the semigroup A. Thus, with Q∗T we define the group automaton M = (Q∗T , Q∗T , P ∗ (Ci , Cj )) where both the set of states and the alphabet are Q∗T and the group operation P ∗ (Ci , Cj ) provides the state transition function. The alphabet of the automaton of the semigroup is its set of states, that is, the alphabet indexes the state reached by the input string thus far. By the Myhill-Nerode, this automata is minimal. By the closure of Q∗T under P ∗ (Ci , Cj ), the set of output states is Q∗T . Augmenting this group automaton with the start state C0 completes the construction of the finite state automata. Equipping this automata with an output alphabet Σ and an associated output map produces the desired Mealy machine. 2

4.1 Power System Model Kopell and Washburn (1982) studied a special case of the three interconnected synchronous generators where one of the machines has significantly lower inertia and is weakly coupled to the larger machines. The dynamic (swing) equations can be simplified to δ (t) = δ¯1 + γ cos(Ωt) + O(γ 2 ) (26) p1 where Ω = β12 [1 + µ1 ] and δ˙3 = ω3 ω˙ 3 = {α3 − β31 sin(δ3 − δ¯1 ) − β32 sin(µ1 δ¯1 + δ3 )} γ{β31 cos(δ3 − δ¯1 ) − µ1 β32 sin(µ1 δ¯1 + δ3 )} (27)

By extending the concept of Shannon entropy to characterize the information storage of digraphs and its rate of increase with increasing tree depth, a characterization

· cos(Ωt) + O(γ 2 ). 511

18th IFAC World Congress (IFAC'11) Milano (Italy) August 28 - September 2, 2011

Poincare Section HΩ3 =0L for Γ=0.525

System Trajectory with Γ=0.525

Poincare Section HΩ3 =0L for Γ=1.175

System Trajectory with Γ=1.175 Ω3

Ω3

1.5 1.0

1.5

2

1.0 1.0

0.5

1

0.5

-1.0

0.5

-0.5

1.0

∆3

∆3

∆3

0.5

0.0 -1.5

-1.0

0.5

-0.5

1.0

1.5

∆3

0.0

-0.5

-0.5

-0.5

-1

-1.0

-1.0

-1.0

-2

-1.5

-1.5

0.0

0.2

0.4

0.6

0.8

1.0

1.2

0.0

0.2

∆1

0.4

0.6

0.8

1.0

1.2

∆1

Fig. 2. (a) (δ3 –ω3 ) phase plane trajectory of stable (γ = 0.525) three machine swing equation and (b) δ1 –δ3 projection of Poincar´e map

Fig. 3. (a) (δ3 –ω3 ) phase plane trajectory of unstable (γ = 1.175) three machine swing equation and (b) δ1 –δ3 projection of Poincar´e map

where δi and ωi are defined as above, γ is a perturbation of the small machine’s equilibrium angle δ¯i , αi are normalized power inputs, βij are normalized impedances, and µi are normalized inertias. This system exhibits chaotic motion and has a solution manifold that is isomorphic to T2 for sufficiently small γ. 4.2 Simulation Results An implementation of language inference method presented in the prequel is used to examine the behaviors of trajectories generated by the system (26)–(27) for several initial conditions that produce markedly different behaviors. The system (27) has an equilibrium point at (δ¯3 , ω¯3 ) = (0.716738, 0). In this section, the trajectories of equations (26)–(27) are examined for several initial conditions of the form (δ¯3 + γ, ω¯3 ) and their behaviors characterized via probabilistic automata.

Fig. 4. Finite State Automata (FSA) representations of symbolic dynamics for (a) stable power system, and (b) unstable power system (two equivalent automata) REFERENCES Crutchfield, J.P. and Young, K. (1989). Inferring statistical complexity. Physical Review Letters, 63(2), 105–108. DyLiacco, T. (1968). Control of power systems via the multi-level concept. Technical Report SRC-68-19, Case Western Reserve University, Systems Research Center. Fink, L. and Carlson, K. (1978). Operating under stress and strain. IEEE Spectrum, 15, 48–53. Kopell, N. and Washburn, R. (1982). Chaotic motions in the two degree-of-freedom swing equations. IEEE Trans. Circuits Syst., CAS-29(11), 738–746. Loparo, K. and Abdel-Malek, F. (1990). A probabilistic approach to dynamic power system security. IEEE Trans. Circuits Syst., 37(6). Wu, F. and Tsai, Y. (1983). Probabilistic dynamic security assessment of power systems: Part I – Basic model. IEEE Trans. Circuits Syst., CAS-30(3).

The (δ3 –ω3 ) phase plane trajectory and Poincar´e map (strobed at ω3 = 0) for system trajectories with γ = 0.525 are presented in Fig. 2. The behavior seen in both the phase and Poincar´e plane reveals a cyclic behavior that is completely captured via a deterministic Finite State Automata consisting of a single closed loop, shown in Subfigure 4(a). The phase plane trajectory and associated Poincar´e map is shown for the case γ = 1.175 in Fig. 3. This trajectory displays more complicated behavior than the stable case and this is reflected in the associated FSA, depicted in Subfigure 4(b). Note that two equivalent automata are recovered here reflecting a temporal shift in the observation process. Either FSA may be used for pattern recognition. 5. CONCLUSIONS The machinery of formal language/automata theory and information theory have been used to develop a methodology for infering probabilistic models of dynamical systems. Moreover, this technique has been applied to a nonlinear power system model and its ability to distinguish between system behaviors in different regions of the state space have been demonstrated. The probabilistic models inferred from the data are light weight and hence promising for use in real-time monitoring applications such as identification of the security state of a power system. 512