Beyond the Spectral Theorem: Spectrally Decomposing Arbitrary Functions of Nondiagonalizable Operators

Nonlinearities in finite dimensions can be linearized by projecting them into infinite dimensions. Unfortunately, often the linear operator techniques that one would then use simply fail since the operators cannot be diagonalized. This curse is well known. It also occurs for finite-dimensional linear operators. We circumvent it by developing a meromorphic functional calculus that can decompose arbitrary functions of nondiagonalizable linear operators in terms of their eigenvalues and projection operators. It extends the spectral theorem of normal operators to a much wider class, including circumstances in which poles and zeros of the function coincide with the operator spectrum. By allowing the direct manipulation of individual eigenspaces of nonnormal and nondiagonalizable operators, the new theory avoids spurious divergences. As such, it yields novel insights and closed-form expressions across several areas of physics in which nondiagonalizable dynamics are relevant, including memoryful stochastic processes, open non unitary quantum systems, and far-from-equilibrium thermodynamics. The technical contributions include the first full treatment of arbitrary powers of an operator. In particular, we show that the Drazin inverse, previously only defined axiomatically, can be derived as the negative-one power of singular operators within the meromorphic functional calculus and we give a general method to construct it. We provide new formulae for constructing projection operators and delineate the relations between projection operators, eigenvectors, and generalized eigenvectors. By way of illustrating its application, we explore several, rather distinct examples.


I. INTRODUCTION
Decomposing a complicated system into its constituent parts-reductionism-is one of science's most power- * pmriechers@ucdavis.edu † chaos@ucdavis.edu ful strategies for analysis and understanding. Largescale systems with linearly coupled components give one paradigm of this success. Each can be decomposed into an equivalent system of independent elements using a similarity transformation calculated by the linear algebra of the system's eigenvalues and eigenvectors. The physics of linear wave phenomena, whether of classical light or quantum mechanical amplitudes, sets the standard of complete reduction rather high. The dynamics is captured by an "operator" whose allowed or exhibited "modes" are the elementary behaviors out of which composite behaviors are constructed by simply weighing each mode's contribution and adding them up. However, one should not reduce a composite system more than is necessary nor, as is increasingly appreciated these days, more than one, in fact, can. Indeed, we live in a complex, nonlinear world whose constituents are strongly interacting. Often their key structures and memoryful behaviors emerge only over space and time. These are the complex systems. Yet, perhaps surprisingly, many complex systems with nonlinear dynamics correspond to linear operators in abstract highdimensional spaces [2][3][4]. And so, there is a sense in which even these complex systems can be reduced to the study of independent nonlocal collective modes. Reductionism, however, faces its own challenges even within its paradigmatic setting of linear systems: linear operators may have interdependent modes with irreducibly entwined behaviors. These irreducible components correspond to so-called nondiagonalizable subspaces. No similarity transformation can reduce them.
In this view, reductionism can only ever be a guide. The actual goal is to achieve a happy medium, as Einstein reminds us, of decomposing a system only to that level at which the parts are irreducible. To proceed, though, begs the original question, What happens when reductionism fails? To answer this requires revisiting one of its more successful implementations, spectral decomposition of completely reducible operators.

A. Spectral Decomposition
Spectral decomposition-splitting a linear operator into independent modes of simple behavior-has greatly accelerated progress in the physical sciences. The impact stems from the fact that spectral decomposition is not only a powerful mathematical tool for expressing the organization of large-scale systems, but also yields predictive theories with directly observable physical consequences [5]. Quantum mechanics and statistical mechanics identify the energy eigenvalues of Hamiltonians as the basic objects in thermodynamics: transitions among the energy eigenstates yield heat and work. The spectrum of eigenvalues reveals itself most directly in other kinds of spectra, such as the frequency spectra of light emitted by the gases that permeate the galactic filaments of our universe [6]. Quantized transitions, an initially mystifying feature of atomic-scale systems, correspond to distinct eigenvectors and discrete spacing between eigenvalues. The corresponding theory of spectral decomposition established the quantitative foundation of quantum mechanics.
The applications and discoveries enabled by spectral decomposition and the corresponding spectral theory fill a long list. In application, direct-bandgap semiconducting materials can be turned into light-emitting diodes (LEDs) or lasers by engineering the spatiallyinhomogeneous distribution of energy eigenvalues and the occupation of their corresponding states [7]. Before their experimental discovery, anti-particles were anticipated as the nonoccupancy of negative-energy eigenstates of the Dirac Hamiltonian [8].
The spectral theory, though, extends far beyond physical science disciplines. In large measure, this arises since the evolution of any object corresponds to a linear dynamic in a sufficiently high-dimensional state space. Even nominally nonlinear dynamics over several variables, the canonical mechanism of deterministic chaos, appear as linear dynamics in appropriate infinite-dimensional shift-spaces [4]. A nondynamic version of rendering nonlinearities into linearities in a higherdimensional feature space is exploited with much success today in machine learning by support vector machines, for example [9]. Spectral decomposition often allows a problem to be simplified by approximations that use only the dominant contributing modes. Indeed, human-face recognition can be efficiently accomplished using a small basis of "eigenfaces" [10].
Certainly, there are many applications that highlight the importance of decomposition and the spectral theory of operators. However, a brief reflection on the mathematical history will give better context to its precise results, associated assumptions, and, more to the point, the generalizations we develop here in hopes of advancing the analysis and understanding of complex systems.
Following on early developments of operator theory by Hilbert and co-workers [11], the spectral theorem for normal operators reached maturity under von Neumann by the early 1930s [12,13]. It became the mathematical backbone of much progress in physics since then, from classical partial differential equations to quantum physics. Normal operators, by definition, commute with their Hermitian conjugate: A † A = AA † . Examples include symmetric and orthogonal matrices in classical mechanics and Hermitian, skew-Hermitian, and unitary operators in quantum mechanics.
The spectral theorem itself is often identified as a collection of related results about normal operators; see, e.g., Ref. [14]. In the case of finite-dimensional vector spaces [15], the spectral theorem asserts that normal operators are diagonalizable and can always be diagonalized by a unitary transformation; that left and right eigenvectors (or eigenfunctions) are simply related by complexconjugate transpose; that these eigenvectors form a complete basis; and that functions of a normal operator reduce to the action of the function on each eigenvalue. Most of these qualities survive with only moderate provisos in the infinite-dimensional case. In short, the spectral theorem makes physics governed by normal operators tractable.
The spectral theorem, though, appears powerless when faced with nonnormal and nondiagonalizable operators. What then are we to do when confronted, say, by complex interconnected systems with nonunitary time evolution, by open systems, by structures that emerge on space and time scales different from the equations of motion, or by other frontiers of physics governed by nonnormal and notnecessarily-diagonalizable operators? Where is the comparably constructive framework for calculations beyond the standard spectral theorem? Fortunately, portions of the necessary generalization have been made within pure mathematics [16], some finding applications in engineering and control [17,18]. However, what is available is incomplete. And, even that which is available is often not in a form adapted to perform calculations that lead to quantitative predictions.

B. Synopsis
Here, we build on previous work in functional analysis and operator theory to provide both a rigorous and constructive foundation for physically relevant calculations involving not-necessarily-diagonalizable operators. In effect, we extend the spectral theorem for normal operators to a broader setting, allowing generalized "modes" of nondiagonalizable systems to be identified and manipulated. The meromorphic functional calculus we develop extends Taylor series expansion and standard holomorphic functional calculus to analyze arbitrary functions of not-necessarily-diagonalizable operators. It readily handles singularities arising when poles (or zeros) of the function coincide with poles of the operator's resolvent-poles that appear precisely at the eigenvalues of the operator. Pole-pole and pole-zero interactions substantially modify the complex-analytic residues within the functional calculus. A key result is that the negative-one power of a singular operator exists in the meromorphic functional calculus. It is the Drazin inverse, a powerful tool that is receiving increased attention in stochastic thermodynamics.
Taken altogether, the functional calculus, Drazin inverse, and methods to manipulate particular eigenspaces, are key to a thorough-going analysis of many complex systems, many now accessible for the first time. Indeed, the framework has already been fruitfully employed by the authors in several specific applications, including closed-form expressions for signal processing and information measures of hidden Markov processes [19][20][21] and for compressing stochastic processes over a quantum channel [22]. However, the techniques are sufficiently general they will be much more widely useful. We envi-sion new opportunities for similar detailed analyses, ranging from biophysics to quantum field theory, wherever restrictions to normal operators and diagonalizability have been roadblocks.
With this broad scope in mind, we develop the mathematical theory first without reference to specific applications and disciplinary terminology. We later give pedagogical (yet, we hope, interesting) examples, exploring several niche, but important applications to finite hidden Markov processes, basic stochastic process theory, nonequilibrium thermodynamics, signal processing, and nonlinear dynamical systems. At a minimum, the examples and their breadth serve to better acquaint readers with the basic methods required to employ the theory.
We introduce the meromorphic functional calculus in §III through §IV, after necessary preparation in §II. §V A further explores eigenprojectors, which we refer to here simply as projection operators. §V B makes explicit their relationship with eigenvectors and generalized eigenvectors. §V B 4 then discusses simplifications of the functional calculus for special cases, while §VI A takes up the spectral properties of transition operators. The examples are discussed at length in §VI before we close in §VII with suggestions on future applications and research directions.

II. SPECTRAL PRIMER
The following is relatively self-contained, assuming basic familiarity with linear algebra at the level of Refs. [15,17]-including eigen-decomposition and knowledge of the Jordan canonical form, partial fraction expansion (see Ref. [23]), and series expansion-and basic knowledge of complex analysis-including the residue theorem and calculation of residues at the level of Ref. [24]. For those lacking a working facility with these concepts, a quick review of §VI's applications may motivate reviewing them. In this section, we introduce our notation and, in doing so, remind the reader of certain basic concepts in linear algebra and complex analysis that will be used extensively in the following.
To begin, we restrict attention to operators with finite representations and only sometimes do we take the limit of dimension going to infinity. That is, we do not consider infinite-rank operators outright. While this runs counter to previous presentations in mathematical physics that consider only infinite-dimensional operators, the upshot is that they-as limiting operators-can be fully treated with a countable point spectrum. We present examples of this later on. Accordingly, we restrict our attention to operators with at most a countably infinite spectrum. Such operators share many features with finite-dimensional square matrices, and so we recall several elementary but essential facts from matrix theory used extensively in the main development.
If A is a finite-dimensional square matrix, then its spectrum is simply the set Λ A of its eigenvalues: where det(·) is the determinant of its argument and I is the identity matrix. The algebraic multiplicity a λ of eigenvalue λ is the power of the term (z − λ) in the characteristic polynomial det(zI − A). In contrast, the geometric multiplicity g λ is the dimension of the kernel of the transformation A − λI or, equivalently, the number of linearly independent eigenvectors associated with the eigenvalue. The algebraic and geometric multiplicities are all equal when the matrix is diagonalizable.
Since there can be multiple subspaces associated with a single eigenvalue, corresponding to different Jordan blocks in the Jordan canonical form, it is structurally important to distinguish the index of the eigenvalue associated with the largest of these subspaces [25]. Definition 1. Eigenvalue λ's index ν λ is the size of the largest Jordan block associated with λ.
If z / ∈ Λ A , then ν z = 0. Note that the index of the operator A itself is sometimes discussed [26]. In such contexts, the index of A is ν 0 . Hence, ν λ corresponds to the index of A − λI.
The index of an eigenvalue gives information beyond what the algebraic and geometric multiplicities themselves yield. Nevertheless, for λ ∈ Λ A , it is always true that ν λ − 1 ≤ a λ − g λ ≤ a λ − 1. In the diagonalizable case, a λ = g λ and ν λ = 1 for all λ ∈ Λ A .
The following employs basic features of complex analysis extensively in conjunction with linear algebra. Let us therefore review several elementary notions in complex analysis. Recall that a holomorphic function is one that is complex differentiable throughout the domain under consideration. A pole of order n at z 0 is a singularity that behaves as h(z)/(z − z 0 ) n as z → z 0 , where h(z) is holomorphic within a neighborhood of z 0 and h(z 0 ) = 0. We say that h(z) has a zero of order m at z 1 if 1/h(z) has a pole of order m at z 1 . A meromorphic function is one that is holomorphic except possibly at a set of isolated poles within the domain under consideration.
Defined over the continuous complex variable z ∈ C, A's resolvent: captures all of A's spectral information through the poles of R(z; A)'s matrix elements. In fact, the resolvent contains more than just A's spectrum: we later show that the order of each pole gives the index ν of the corresponding eigenvalue.
The spectrum Λ A can be expressed in terms of the resolvent. Explicitly, the point spectrum (i.e., the set of eigenvalues) is the set of complex values z at which zI −A is not a one-to-one mapping, with the implication that the inverse of zI − A does not exist: where inv(·) is the inverse of its argument. Later, via our investigation of the Drazin inverse, it should become clear that the resolvent operator can be self-consistently defined at the spectrum, despite the lack of inverse.
For infinite-rank operators, the spectrum becomes more complicated. In that case, the right point spectrum (the point spectrum of A) need not be the same as the left point spectrum (the point spectrum of A's dual A ). Moreover, the spectrum may grow to include noneigenvalues z for which the range of zI − A is not dense in the vector space it transforms or for which zI − A has dense range but the inverse of zI − A is not bounded. These two settings give rise to the so-called residual spectrum and continuous spectrum, respectively [27]. To mitigate confusion, it should be noted that the point spectrum can be continuous, yet never coincides with the continuous spectrum just described. Moreover, understanding only countable point spectra is necessary to follow the developments here.
Each of A's eigenvalues λ has an associated projection operator A λ , which is the residue of the resolvent as z → λ [14]. Explicitly: where Res( · , z → λ) is the element-wise residue of its first argument as z → λ. The projection operators are orthonormal: and sum to the identity: The following discusses in detail and then derives several new properties of projection operators.

III. FUNCTIONAL CALCULI
In the following, we develop an extended functional calculus that makes sense of arbitrary functions f (·) of a linear operator A. Within any functional calculus, one considers how A's eigenvalues map to the eigenvalues of f (A); which we call a spectral mapping. For example, it is known that holomorphic functions of bounded linear operators enjoy an especially simple spectral mapping theorem [28]: To fully appreciate the meromorphic functional calculus, we first state and compare the main features and limitations of alternative functional calculi.

A. Taylor series
Inspired by the Taylor expansion of scalar functions: a calculus for functions of an operator A can be based on the series: where f (n) (ξ) is the n th derivative of f (z) evaluated at z = ξ. This is often used, for example, to express the exponential of A as: This particular series-expansion is convergent for any A since e z is entire, in the sense of complex analysis. Unfortunately, even if it exists there is a limited domain of convergence for most functions. For example, suppose f (z) has poles and choose a Maclaurin series; i.e., ξ = 0 in Eq. (3). Then the series only converges when A's spectral radius is less than the radius of the innermost pole of f (z). Addressing this and related issues leads directly to alternative functional calculi.

B. Holomorphic functional calculus
Holomorphic functions are well behaved, smooth functions that are complex differentiable. Given a function f (·) that is holomorphic within a disk enclosed by a counterclockwise contour C, its Cauchy integral formula is given by: Taking this as inspiration, the holomorphic functional calculus performs a contour integration of the resolvent to extend f (·) to operators: where C Λ A is a closed counterclockwise contour that encompasses Λ A . Assuming that f (z) is holomorphic at z = λ for all λ ∈ Λ A , a nontrivial calculation [25] shows that Eq. (5) is equivalent to the holomorphic calculus defined by: After some necessary development, we will later derive Eq. (6) as a special case of our meromorphic functional calculus, such that Eq. (6) is valid whenever f (z) is holomorphic at z = λ for all λ ∈ Λ A . The holomorphic functional calculus was first proposed in Ref. [25] and is now in wide use; e.g., see Ref. [17, p. 603]. It agrees with the Taylor-series approach whenever the infinite series converges, but gives an functional calculus when the series approach fails. For example, using the principal branch of the complex logarithm, the holomorphic functional calculus admits log(A) for any nonsingular matrix, with the satisfying result that e log(A) = A. Whereas, the Taylor series approach fails to converge for the logarithm of most matrices even if the expansion for, say, log(1 − z) is used.
The major shortcoming of the holomorphic functional calculus is that it assumes f (z) is holomorphic at Λ A . Clearly, if f (z) has a pole at some z ∈ Λ A , then Eq. (6) fails. An example of such a failure is the negative-one power of a singular operator, which we take up later on.
Several efforts have been made to extend the holomorphic functional calculus. For example, Refs. [29] and [30] define a functional calculus that extends the standard holomorphic functional calculus to include a certain class of meromorphic functions that are nevertheless still required to be holomorphic on the point spectrum (i.e., on the eigenvalues) of the operator. However, we are not aware of any previous work that introduces and develops the consequences of a functional calculus for functions that are meromorphic on the point spectrum-which we take up in the next few sections.

C. Meromorphic functional calculus
Meromorphic functions are holomorphic except at a set of isolated poles of the function. The resolvent of a finite-dimensional operator is meromorphic, since it is holomorphic everywhere except for poles at the eigenvalues of the operator. We will now also allow our function f (z) to be meromorphic with possible poles that coincide with the poles of the resolvent.
Inspired again by the Cauchy integral formula of Eq. (4), but removing the restriction to holomorphic functions, our meromorphic functional calculus instead employs a partitioned contour integration of the resolvent: where C λ is a small counterclockwise contour around the eigenvalue λ. This and a spectral decomposition of the resolvent (to be derived later) extends the holomorphic calculus to a much wider domain, defining: The contour is integrated using knowledge of f (z) since meromorphic f (z) can introduce poles and zeros at Λ A that interact with the resolvent's poles.
The meromorphic functional calculus agrees with the Taylor-series approach whenever the series converges and agrees with the holomorphic functional calculus whenever f (z) is holomorphic at Λ A . However, when both the previous functional calculi fail, the meromorphic calculus extends the domain of f (A) to yield surprising, yet sensible answers. For example, we show that within it, the negative-one power of a singular operator is the Drazin inverse-an operator that effectively inverts everything that is invertible.
The major assumption of our meromorphic functional calculus is that the domain of operators must have a spectrum that is at most countably infinite-e.g., A can be any compact operator. A related limitation is that singularities of f (z) that coincide with Λ A must be isolated singularities. Nevertheless, we expect that these restrictions can be lifted with proper treatment, as discussed in fuller context later.

IV. MEROMORPHIC SPECTRAL DECOMPOSITION
The preceding gave an overview of the relationship between alternative functional calculi and their trade-offs, highlighting the advantages of the meromorphic functional calculus. This section leverages these advantages and employs a partial fraction expansion of the resolvent to give a general spectral decomposition of almost any function of any operator. Then, since it plays a key role in applications, we apply the functional calculus to investigate the negative-one power of singular operators-thus deriving, what is otherwise an operator defined axiomatically, the Drazin inverse from first principles.

A. Partial fraction expansion of the resolvent
The elements of A's resolvent are proper rational functions that contain all of A's spectral information. (Recall that a proper rational function r(z) is a ratio of polynomials in z whose numerator has degree strictly less than the degree of the denominator.) In particular, the resolvent's poles coincide with A's eigenvalues since, for z / ∈ Λ A : where a λ is the algebraic multiplicity of eigenvalue λ and C is the matrix of cofactors of zI − A. That is, C's transpose C is the adjugate of zI − A: , whose elements will be polynomial functions of z of degree less than λ∈Λ A a λ . Recall that the partial fraction expansion of a proper rational function r(z) with poles in Λ allows a unique decomposition into a sum of constant numerators divided by monomials in z − λ up to degree a λ , when a λ is the order of the pole of r(z) at λ ∈ Λ [23]. Equation (8) thus makes it clear that the resolvent has the unique partial fraction expansion: where {A λ,m } is the set of matrices with constant entries (not functions of z) uniquely determined elementwise by the partial fraction expansion. However, R(z; A)'s poles are not necessarily of the same order as the algebraic multiplicity of the corresponding eigenvalues since the entries of C, and thus of C , may have zeros at A's eigenvalues. This has the potential to render A λ,m equal to the zero matrix 0. The Cauchy integral formula indicates that the constant matrices {A λ,m } of Eq. (9) can be obtained by the residues: where the residues are calculated elementwise. The projection operators A λ associated with each eigenvalue λ were already referenced in §II, but can now be properly introduced as the A λ,0 matrices: Since R(z; A)'s elements are rational functions, as we just showed, it is analytic except at a finite number of isolated singularities-at A's eigenvalues. In light of the residue theorem, this motivates the Cauchy-integral-like formula that serves as the starting point for the meromorphic functional calculus: Let's now consider several immediate consequences.

B. Decomposing the identity
Even the simplest applications of Eq. (13) yield insight. Consider the identity as the operator function f (A) = A 0 = I that corresponds to the scalar function f (z) = z 0 = 1. Then, Eq. (13) implies: This shows that the projection operators are, in fact, a decomposition of the identity, as anticipated in Eq. (2).

C. Dunford decomposition, decomposed
For f (A) = A, Eqs. (13) and (10) imply that: We denote the important set of nilpotent matrices A λ,1 that project onto the generalized eigenspaces by relabeling them: Equation (14) is the unique Dunford decomposition [16]: The special case where A is diagonalizable implies that N = 0. And so, Eq. (14) simplifies to:

D. The resolvent, resolved
As shown in Ref. [14] and can be derived from Eqs. (12) and (16): Due to these, our spectral decomposition of the Dunford decomposition implies that: Moreover: It turns out that for m > 0: A λ,m = N m λ . (See also Ref. [14, p. 483].) This leads to a generalization of the projection operator orthonormality relations of Eq. (1). Most generally, the operators of {A λ,m } are mutually related by: Finally, if we recall that the index ν λ is the dimension of the largest associated subspace, we find that the index of λ characterizes the nilpotency of N λ : N m λ = 0 for m ≥ ν λ . That is: Returning to Eq. (9), we see that all A λ,m with m ≥ ν λ are zero-matrices and so do not contribute to the sum. Thus, we can rewrite Eq. (9) as: or: for z / ∈ Λ A . The following sections sometimes use A λ,m in place of A λ A − λI m . This is helpful both for conciseness and when applying Eq. (19). Nonetheless, the equality in Eq. (18) is a useful one to keep in mind.

E. Meromorphic functional calculus
In light of Eq. (13), Eq. (21) together with Eq. (18) allow us to express any function of an operator simply and solely in terms of its spectrum (i.e., its eigenvalues for the finite dimensional case), its projection operators, and itself: In obtaining Eq. (23) we finally derived Eq. (7), as promised earlier in § III C. Effectively, by modulating the modes associated with the resolvent's singularities, the scalar function f (·) is mapped to the operator domain, where its action is expressed in each of A's independent subspaces.

F. Evaluating the residues
Interpretation aside, how does one use this result? Equation (23) says that the spectral decomposition of f (A) reduces to the evaluation of several residues, where: So, to make progress with Eq. (23), we must evaluate function-dependent residues of the form: If f (z) were holomorphic at each λ, then the order of the pole would simply be the power of the denominator. We could then use Cauchy's differential formula for holomorphic functions: for f (z) holomorphic at a. And, the meromorphic calculus would reduce to the holomorphic calculus. Often, f (z) will be holomorphic at least at some of A's eigenvalues. And so, Eq. (24) is still locally a useful simplification in those special cases. In general, though, f (z) introduces poles and zeros at λ ∈ Λ A that change their orders. This is exactly the impetus for the generalized functional calculus. The residue of a complex-valued function g(z) around its isolated pole λ of order n + 1 can be calculated from: (23) says that we can explicitly derive the spectral decomposition of powers of the operator A. Of course, we already did this for the special cases of A 0 and A 1 . The goal, though, is to do this in general.
For f (A) = A L → f (z) = z L , z = 0 can be either a zero or a pole of f (z), depending on the value of L. In either case, an eigenvalue of λ = 0 will distinguish itself in the residue calculation of A L via its unique ability to change the order of the pole (or zero) at z = 0. For example, at this special value of λ and for integer L > 0, λ = 0 induces poles that cancel with the zeros of f (z) = z L , since z L has a zero at z = 0 of order L. For integer L < 0, an eigenvalue of λ = 0 increases the order of the z = 0 pole of f (z) = z L . For all other eigenvalues, the residues will be as expected. Hence, from Eq. (23) and inserting f (z) = z L , for any L ∈ C: where L m is the generalized binomial coefficient: and [0 ∈ Λ A ] is the Iverson bracket which takes on value 1 if zero is an eigenvalue of A and 0 if not. A λ,m was replaced by A λ (A − λI) m to suggest the more explicit calculations involved with evaluating any A L . Equation (25) applies to any linear operator with only isolated singularities in its resolvent. If L is a nonnegative integer such that L ≥ ν λ − 1 for all λ ∈ Λ A , then: where L m is now reduced to the traditional binomial coefficient L!/(m!(L − m)!).

H. Drazin inverse
If L is any negative integer, then −|L| m can be written as a traditional binomial coefficient (−1) m |L|+m−1 m , yielding: Thus, negative powers of an operator can be consistently defined even for noninvertible operators. In light of Eqs. (25) and (28), it appears that the zero eigenvalue does not even contribute to the function. It is well known, in contrast, that it wreaks havoc on the naive, oft-quoted definition of a matrix's negative power: since this would imply dividing by zero. If we can accept large positive powers of singular matrices-for which the zero eigenvalue does not contribute-it seems fair to also accept negative powers that likewise involve no contribution from the zero eigenvalue. Editorializing aside, we note that extending the definition of A −1 to the domain including singular operators via Eqs. (25) and (28) implies that: which is a very sensible and desirable condition. Moreover, we find that AA −1 = I − A 0 . Specifically, the negative-one power of any square matrix is in general not the same as the matrix inverse since inv(A) need not exist. However, it is consistently defined via Eq. (28) to be: This is the Drazin inverse A D of A. Note that it is not the same as the Moore-Penrose pseudo-inverse [31,32]. Although the Drazin inverse is usually defined axiomatically to satisfy certain criteria [33], it is naturally derived as the negative one power of a singular operator in the meromorphic functional calculus. We can check that it indeed satisfies the axiomatic criteria for the Drazin inverse, enumerated according to historical precedent: which gives rise to the Drazin inverse's moniker as the {1 ν0 , 2, 5}-inverse [33]. While A −1 always exists, the resolvent is nonanalytic at z = 0 for a singular matrix. Effectively, the meromorphic functional calculus removes the nonanalyticity of the resolvent in evaluating A −1 . As a result, as we can see from Eq. (29), the Drazin inverse inverts what is invertible; the remainder is zeroed out.
Of course, whenever A is invertible, A −1 is equal to inv(A). However, we should not confuse this coincidence with equivalence. Moreover, despite historic notation there is no reason that the negative-one power should in general be equivalent to the inverse. Especially, if an operator is not invertible! To avoid confusing A −1 with inv(A), we use the notation A D for the Drazin inverse of A. Still, A D = inv(A), whenever 0 / ∈ Λ A . Amusingly, this extension of previous calculi lets us resolve an elementary but fundamental question: What is 0 −1 ? It is certainly not infinity. Indeed, it is just as close to negative infinity! Rather: Although Eq. (29) is a constructive way to build the Drazin inverse, it imposes more work than is actually necessary. Using the meromorphic functional calculus, we can derive a new, simple construction of the Drazin inverse that requires only the original operator and the eigenvalue-0 projector.
First, assume that λ is an isolated singularity of R(z; A) with finite separation at least distance from the nearest neighboring singularity. And, consider the operator-valued function f λ defined via the RHS of: with λ+ e iφ defining an -radius circular contour around λ. Then we see that: where [z ∈ C : |z − λ| < ] is the Iverson bracket that takes on value 1 if z is within -distance of λ and 0 if not.
Second, we use this to find that, for any c ∈ C \ {0}: where we asserted that the contour C 0 exists within the finite -ball about the origin. Third, we note that A + cA 0 is invertible for all c = 0; this can be proven by multiplying each side of Eq.
Finally, multiplying each side of Eq. (31) by I − A 0 , and recalling that A 0,0 A 0,m = A 0,m , we find a useful expression for calculating the Drazin inverse of any linear operator A, given only A and A 0 . Specifically: which is valid for any c ∈ C \ {0}. Eq. (32) generalizes the result found specifically for c = −1 in Ref. [34].
For the special case of c = −1, it is worthwhile to also consider the alternative construction of the Drazin inverse implied by Eq. (31): By a spectral mapping (λ → 1 − λ, for λ ∈ Λ T ), the Perron-Frobenius theorem and Eq. (31) yield an important consequence for any stochastic matrix T . The Perron-Frobenius theorem guarantees that T 's eigenvalues along the unit circle are associated with a diagonalizable subspace. In particular, ν 1 = 1. Spectral mapping of this result means that T 's eigenvalue 1 maps to the eigenvalue 0 of I − T and T 1 = (I − T ) 0 . Moreover: since ν 0 = 1. This corollary of Eq. (31) (with c = 1) corresponds to a number of important and well known results in the theory of Markov processes. Indeed, Z ≡ (I − T + T 1 ) −1 is called the fundamental matrix in that setting [35].

I. Consequences and generalizations
For an infinite-rank operator A with a continuous spectrum, the meromorphic functional calculus has the natural generalization: where the contour C Λ A encloses the (possibly continuous) spectrum of A without including any unbounded contributions from f (z) outside of C Λ A . The function f (z) is expected to be meromorphic within C Λ A . This again deviates from the holomorphic approach, since the holomorphic functional calculus requires that f (z) is analytic in a neighborhood around the spectrum; see § VII of Ref. [36]. Moreover, Eq. (34) allows an extension of the functional calculus of Refs. [29,30,37], since the function can be meromorphic at the point spectrum in addition being meromorphic on the residual and continuous spectra.
In either the finite-or infinite-rank case, whenever f (z) is analytic in a neighborhood around the spectrum, the meromorphic functional calculus agrees with the holomorphic. Whenever f (z) is not analytic in a neighborhood around the spectrum, the function is undefined in the holomorphic approach. In contrast, the meromorphic approach extends the function to the operator-valued domain and does so with novel consequences.
In particular, when f (z) is not analytic in a neighborhood around the spectrum-say f (z) is nonanalytic within A's spectrum at Ξ f ⊂ Λ A -then we expect to lose both homomorphism and spectral mapping properties: A simple example of both losses arises with the Drazin inverse, above. There, f 1 (z) = z −1 . Taking this and f 2 (z) = z combined with singular operator A leads to the loss of homomorphism: A D A = I. As for the second property, the spectral mapping can be altered for the candidate spectra at Ξ f via pole-pole or pole-zero interactions in the complex contour integral. For f (A) = A −1 , how does A's eigenvalue of 0 get mapped into the new spectrum of A D ? A naive application of the spectral mapping theorem might seem to yield an undefined quantity. But, using the meromorphic functional calculus selfconsistently maps the eigenvalue as 0 −1 = 0. It remains to be explored whether the full spectral mapping is preserved for any function f (A) under the meromorphic interpretation of f (λ).
It should now be apparent that extending functions via the meromorphic functional calculus allows one to express novel mathematical properties, some likely capable of describing new physical phenomena. At the same time, extra care is necessary. The situation is reminiscent of the loss of commutativity in non-Abelian operator algebra: not all of the old rules apply, but the gain in nuance allows for mathematical description of important phenomena.
We chose to focus primarily on the finite-rank case here since it is sufficient to demonstrate the utility of the general projection-operator formalism. Indeed, there are ample nontrivial applications in the finite-rank setting that deserve attention. To appreciate these, we now turn to address the construction and properties of general eigenprojectors.

V. CONSTRUCTING DECOMPOSITIONS
At this point, we see that projection operators are fundamental to functions of an operator. This prompts the practical question of how to actually calculate them. The next several sections address this by deriving expressions with both theoretical and applied use. We first address the projection operators associated with index-one eigenvalues. We then explicate the relationship between eigenvectors, generalized eigenvectors, and projection operators for normal, diagonalizable, and general matrices. Finally, we address how the general results specialize in several common cases of interest. After these, we turn to examples and applications.

A. Projection operators of index-one eigenvalues
To obtain the projection operators associated with each index-one eigenvalue λ ∈ {ζ ∈ Λ A : ν ζ = 1}, we apply the meromorphic calculus to an appropriately chosen function of A, finding: Therefore, if ν λ = 1: As convenience dictates in our computations, we let ν ζ → a ζ − g ζ + 1 or even ν ζ → a ζ in Eq. (35), since multiplying A λ by (A − ζI)/(λ − ζ) has no effect for ζ ∈ Λ A \ {λ} if ν λ = 1. Equation (35) generalizes a well known result that applies when the index of all eigenvalues is one. That is, when the operator is diagonalizable, we have: To the best of our knowledge, Eq. (35) is original. Since eigenvalues can have index larger than one, not all projection operators of a nondiagonalizable operator can be found directly from Eq. (35). Even so, it serves three useful purposes. First, it gives a practical reduction of the eigen-analysis by finding all projection operators of index-one eigenvalues. Second, if there is only one eigenvalue that has index larger than one-what we call the almost diagonalizable case-then Eq. (35), together with the fact that the projection operators must sum to the identity, does give a full solution to the set of projection operators. Third, Eq. (35) is a powerful theoretical tool that we can use directly to spectrally decompose functions, for example, of a stochastic matrix whose eigenvalues on the unit circle are guaranteed to be index-one by the Perron-Frobenius theorem.
Although index-one expressions have some utility, we need a more general procedure to obtain all projection operators of any linear operator. Recall that, with full generality, projection operators can also be calculated directly via residues, as in Eq. (12).
An alternative procedure-one that extends a method familiar at least in quantum mechanics-is to obtain the projection operators via eigenvectors. However, quantum mechanics always concerns itself with a subset of diagonalizable operators. What is the necessary generalization? For one, left and right eigenvectors are no longer simply conjugate transposes of each other. More severely, a full set of spanning eigenvectors is no longer guaranteed and we must resort to generalized eigenvectors. Since the relationships among eigenvectors, generalized eigenvectors, and projection operators are critical to the practical calculation of many physical observables of complex systems, we collect these results in the next section.

B. Eigenvectors, generalized eigenvectors, and projection operators
Two common questions regarding projection operators are: Why not just use eigenvectors? And, why not use the Jordan canonical form? First, the eigenvectors of a defective matrix do not form a complete basis with which to expand an arbitrary vector. One needs generalized eigenvectors for this. Second, some functions of an operator require removing, or otherwise altering, the contribution from select eigenspaces. This is most adroitly handled with the projection operator formalism where different eigenspaces (correlates of Jordan blocks) can effectively be treated separately. Moreover, even for simple cases where eigenvectors suffice, the projection operator formalism simply can be more calculationally or mathematically convenient.
That said, it is useful to understand the relationship between projection operators and generalized eigenvectors. For example, it is often useful to create projection operators from generalized eigenvectors. This section clarifies their connection using the language of matrices. In the most general case, we show that the projection operator formalism is usefully concise. Since an eigenvalue λ ∈ Λ A has algebraic multiplicity a λ ≥ 1, λ appears a λ times in the ordered tuple.
Assuming A is normal, each projection operator A λ can be constructed as the sum of all ket-bra pairs of righteigenvectors corresponding to λ composed with their conjugate transpose. We later introduce bras and kets more generally via generalized eigenvectors of the operator A and its dual A . However, since the complex-conjugate transposition rule between dual spaces is only applicable to a ket basis derived from a normal operator, we put off using the bra-ket notation for now so as not to confuse the more familiar "normal" case with the general case.
To explicitly demonstrate this relationship between projection operators, eigenvectors, and their Hermitian adjoints in the case of normality, observe that: Evidently, for normal matrices A: And, since u † i u j = δ i,j , we have an orthogonal set {A λ } λ∈Λ A with the property that: Moreover: and so on. All of the expected properties of projection operators can be established again in this restricted setting.
The rows of U −1 = U † are A's left-eigenvectors. In this case, they are simply the conjugate transpose of the righteigenvector. Note that conjugate transposition is the familiar transformation rule between ket and bra spaces in quantum mechanics (see e.g., Ref. [38])-a consequence of the restriction to normal operators, as we will show. Importantly, a more general formulation of quantum mechanics would not have this same restricted correspondence between the dual ket and bra spaces.
To elaborate on this point, recall that vector spaces admit dual spaces and dual bases. However, there is no sense of a dual correspondence of a single ket or bra without reference to a full basis [15]. Implicitly in quantum mechanics, the basis is taken to be the basis of eigenstates of any Hermitian operator, nominally since observables are self-adjoint.
To allude to an alternative, we note that u † j u j is not only the Hermitian form of inner product u j , u j (where ·, · denotes the inner product) of the right eigenvector u j with itself, but importantly also the simple dotproduct of the left eigenvector u † j and the right eigenvector u j , where u † j acts as a linear functional on u j . Contrary to the substantial effort devoted to the innerproduct-centric theory of Hilbert spaces, this latter interpretation of u † j u j -in terms of linear functionals and a left-eigenvector basis for linear functionals-is what generalizes to a consistent and constructive framework for the spectral theory beyond normal operators, as we will see shortly.

Diagonalizable matrices
By definition, diagonalizable matrices can be diagonalized, but not necessarily via a unitary transformation. All diagonalizable matrices can nevertheless be diagonalized via the transformation: A = P ΛP −1 , where the columns of the square matrix P are the not-necessarilyorthogonal right eigenvectors of A corresponding to the eigenvalues ordered along the diagonal matrix Λ and where the rows of P −1 are A's left eigenvectors. Importantly, the left eigenvectors need not be the Hermitian adjoint of the right eigenvectors. As a particular example, this more general setting is required for almost any transition dynamic of a Markov chain. In other words, the transition dynamic of any interesting complex network with irreversible processes serves as an example of a nonnormal operator.
With these definitions in place, the projection operators for a diagonalizable matrix can be written: Then: So, we see that the projection operators introduced earlier in a coordinate-free manner have a concrete representation in terms of left and right eigenvectors when the operator is diagonalizable.

Any matrix
Not all matrices can be diagonalized, but all square matrices can be put into Jordan canonical form via the transformation: A = Y JY −1 [17]. Here, the columns of the square matrix Y are the linearly independent right eigenvectors and generalized right eigenvectors corresponding to the Jordan blocks ordered along the diagonal of the block-diagonal matrix J. And, the rows of Y −1 are the corresponding left eigenvectors and generalized left eigenvectors, but reverse-ordered within each block, as we will show.
Let there be n Jordan blocks forming the n-tuple (J 1 , J 2 , . . . , J n ), with 1 ≤ n ≤ M . The k th Jordan block J k has dimension m k -by-m k : Note that eigenvalue λ ∈ Λ A corresponds to g λ different Jordan blocks, where g λ is the geometric multiplicity of the eigenvalue λ. Indeed: Moreover, the index ν λ of the eigenvalue λ is defined as the size of the largest Jordan block corresponding to λ. So, we write this in the current notation as: If the index of any eigenvalue is greater than one, then the conventional eigenvectors do not span the Mdimensional vector space. However, the set of M generalized eigenvectors does form a basis for the vector space [39].
Given the n-tuple of possibly-degenerate eigenvalues (Λ A ) = (λ 1 , λ 2 , . . . , λ n ), there is a corresponding ntuple of m k -tuples of linearly-independent generalized right-eigenvectors: and a corresponding n-tuple of m k -tuples of linearlyindependent generalized left-eigenvectors: and: for 0 ≤ m ≤ m k − 1, where |λ for m, n, ∈ {0, 1, . . . , m k } and ≥ 0. At the same time, it is then easy to show that: where m, n ∈ {0, 1, . . . , m k }. Imposing appropriate normalization, we find that: Hence, we see that the left eigenvectors and generalized eigenvectors are a dual basis to the right eigenvectors and generalized eigenvectors. Interestingly though, within each Jordan subspace, the most generalized left eigenvectors are dual to the least generalized right eigenvectors, and vice versa.
(To be clear, in this terminology "least generalized" eigenvectors are the standard eigenvectors. For example, the λ   With these details worked out, we find that the projection operators for a nondiagonalizable matrix can be written as: And, we see that a projection operator includes all of its left and right eigenvectors and all of its left and right generalized eigenvectors. This implies that the identity operator must also have a decomposition in terms of both eigenvectors and generalized eigenvectors: Let |λ (m) k m k m=1 denote the column vector: and let λ Then, using the above results, and the fact that Eq. (37) implies that λ k |, we derive the explicit generalized-eigenvector decomposition of the nondiagonalizable operator A: where, defining Y as: we are forced by Eq. (38) to recognize that: since then Y −1 Y = I, and we recall that the inverse is guaranteed to be unique. The above demonstrates an explicit construction for the Jordan canonical form. One advantage we learn from this explicit decomposition is that the complete set of left eigenvectors and left generalized eigenvectors (encapsulated in Y −1 ) can be obtained from the inverse of the matrix of the complete set of right eigenvectors and generalized right eigenvectors (encoded in Y ) and vice versa. One unexpected lesson, though, is that the generalized left eigenvectors appear in reverse order within each Jordan block.
Using Eqs. (39) and (18) with Eq. (37), we see that the nilpotent operators A λ,m with m > 0 further link the various generalized eigenvectors within each subspace k. Said more suggestively, generalized modes of a nondiagonalizable subspace are necessarily cooperative.
It is worth noting that the left eigenvectors and generalized left eigenvectors form a basis for all linear functionals of the vector space spanned by the right eigenvectors and generalized right eigenvectors. Moreover, the left eigenvectors and generalized left eigenvectors are exactly the dual basis to the right eigenvectors and generalized right eigenvectors by their orthonormality properties. However, neither the left nor right eigen-basis is a priori more fundamental to the operator. Sympathetically, the right eigenvectors and generalized eigenvectors form a (dual) basis for all linear functionals of the vector space spanned by the left eigenvectors and generalized eigenvectors.

Simplified calculi for special cases
In special cases, the meromorphic functional calculus reduces the general expressions above to markedly simpler forms. And, this can greatly expedite practical calculations and provide physical intuition. Here, we show which reductions can be used under which assumptions.
For functions of operators with a countable spectrum, recall that the general form of the meromorphic functional calculus is: Equations (18) and (39) gave the method to calculate A λ,m in terms of eigenvectors and generalized eigenvectors.
When the operator is diagonalizable (not necessarily normal), this reduces to: where A λ can now be constructed from conventional right and left eigenvectors, although λ j | is not necessarily the conjugate transpose of |λ j . When the function is analytic on the spectrum of the (not necessarily diagonalizable) operator, then our functional calculus reduces to the holomorphic functional calculus: When the function is analytic on the spectrum of a diagonalizable (not necessarily normal) operator this reduces yet again to: When the function is analytic on the spectrum of a diagonalizable (not necessarily normal) operator with no degeneracy this reduces even further to: Finally, recall that an operator is normal when it commutes with its conjugate transpose. If the function is analytic on the spectrum of a normal operator, then we recover the simple form enabled by the spectral theorem of normal operators familiar in physics. That is, Eq. (43) is applicable, but now we have the extra simplification that λ j | is simply the conjugate transpose of |λ j : λ j | = |λ j † .

VI. EXAMPLES AND APPLICATIONS
To illustrate the use and power of the meromorphic functional calculus, we now adapt it to analyze a suite of applications from quite distinct domains. First, we point to a set of example calculations for finite-dimensional operators of stochastic processes. Second, we show that the familiar Poisson process is intrinsically nondiagonalizable, and hint that nondiagonalizability may be common more generally in semi-Markov processes. Third, we illustrate how commonly the Drazin inverse arises in nonequilibrium thermodynamics, giving a roadmap to developing closed-from expressions for a number of key observables. Fourth, we turn to signal analysis and comment on power spectra of processes generated by nondiagonalizable operators. Finally, we round out the applications with a general discussion of Ruelle-Frobenius-Perron and Koopman operators for nonlinear dynamical systems.

A. Spectra of stochastic transition operators
The preceding employed the notation that A represents a general linear operator. In the following examples, we reserve the symbol T for the operator of a stochastic transition dynamic. If the state-space is finite and has a stationary distribution, then T has a representation that is a nonnegative row-stochastic-all rows sum to unitytransition matrix.
The transition matrix's nonnegativity guarantees that for each λ ∈ Λ T its complex conjugate λ is also in Λ T . Moreover, the projection operator associated with the complex conjugate of λ is the complex conjugate of T λ : If the dynamic induced by T has a stationary distribution over the state space, then the spectral radius of T is unity and all of T 's eigenvalues lie on or within the unit circle in the complex plane. The maximal eigenvalues have unity magnitude and 1 ∈ Λ T . Moreover, an extension of the Perron-Frobenius theorem guarantees that eigenvalues on the unit circle have algebraic multiplicity equal to their geometric multiplicity. And, so, ν ζ = 1 for all ζ ∈ {λ ∈ Λ T : |λ| = 1}.
T 's index-one eigenvalue of λ = 1 is associated with stationarity of the associated Markov process. T 's other eigenvalues on the unit circle are roots of unity and correspond to deterministic periodicities within the process.
All of these results carry over from discrete to continuous time. In continuous time, where e tG = T t0→t0+t , T 's stationary eigenvalue of unity maps to G's stationary eigenvalue of zero. If the dynamic has a stationary distribution over the state space, then the rate matrix G is row-sum zero rather than row-stochastic. T 's eigenvalues, on or within the unit circle, map to G's eigenvalues with nonpositive real part in the left-hand side of the complex plane.
To reduce ambiguity in the presence of multiple operators, functions of operators, and spectral mapping, we occasionally denote eigenvectors with subscripted operators on the eigenvalues within the bra or ket. For example, |0 G = |1 T = |0 G = |1 T = |0 T disambiguates the identification of |0 when we have operators G, T , G, and T with T = e τ G , T = e τ G , and 0 ∈ Λ G , Λ G , Λ T .

B. Randomness and memory in correlated processes
The generalized spectral theory developed here has recently been applied to give the first closed-form expressions for many measures of complexity for stochastic processes that can be generated by probabilistic finite automata [19]. Rather than belabor the Kolmogorov-Chaitin notion of complexity which is inherently uncomputable [40], the new analytic framework infuses computational mechanics [41] with a means to compute very practical answers about an observed system's organization and to address the challenges of prediction.
For example, we can now answer the obvious questions regarding prediction: How random is a process? How much information is shared between the past and the future? How far into the past must we look to predict what is predictable about the future? How much about the observed history must be remembered to predict what is predictable about the future? And so on. The Supplementary Materials of Ref. [19] exploit the meromorphic functional calculus to answer these (and more) questions for the symbolic dynamics of a chaotic map, the spacetime domain for an elementary cellular automata, and the chaotic crystallographic structure of a close-packed polytypic material as determined from experimental Xray diffractograms.
In the context of the current exposition, the most notable feature of the analyses across these many domains is that our imposed questions, which entail tracking an observer's state of knowledge about a process, necessarily induce a nondiagonalizable metadynamic that becomes the central object of analysis in each case. (This metadynamic is the so-called mixed-state presentation of Refs. [42,43].) This theme, and the inherent nondiagonalizability of prediction, is explored in greater depth elsewhere [44]. We also found that another nondiagonalizable dynamic is induced even in the context of quantum communication when determining how much memory reduction can be achieved if we generate a classical stochastic process using quantum mechanics [22].
We mention the above nondiagonalizable metadynamics primarily as a pointer to concrete worked-out examples where the meromorphic functional calculus has been employed to analyze finitary hidden Markov processes via explicitly calculated, generalized eigenvectors and projection operators. We now return to a more self-contained discussion, where we show that nondiagonalizability can be induced by the simple act of counting. Moreover, the theory developed is then applied to deliver quick and powerful results.

C. Poisson point processes
The meromorphic functional calculus leads naturally to a novel perspective on the familiar Poisson counting process-a familiar stochastic process class used widely across physics and other quantitative sciences to describe "completely random" event durations that occur over a continuous domain [45][46][47][48]. The calculus shows that the basic Poisson distribution arises as the signature of a simple nondiagonalizable dynamic. More to the point, we derive the Poisson distribution directly, without requiring the limit of the discrete-time binomial distribution, as conventionally done [24].
Consider all possible counts, up to some arbitrarily large integer N . The dynamics among these first N + 1 counter states constitute what can be called the truncated Poisson dynamic. We recover the full Poisson distribution as N → ∞. A Markov chain for the truncated Poisson dynamic is shown in Fig. 1. The corresponding rate matrix G, for any arbitrarily large truncation N of the possible count, is: where G ij is the rate of transitioning to state (count) j given that the system is in state (count) i. Elements not on either the main diagonal or first superdiagonal are zero. This can be rewritten succinctly as: where I is the identity operator in N -dimensions and D 1 is the upshift-by-1 matrix in N -dimensions, with zeros everywhere, except 1s along the first superdiagonal. Let us also define the upshift-by-m matrix D m with zeros everywhere except 1s along the m th superdiagonal, such that D m = D m 1 and D n m = D m·n , with D 0 = I. Operationally, if δ | is the probability distribution over counter states that is peaked solely at state , then δ | D m = δ +m |.
For any arbitrarily large N , G's eigenvalues are given by det(G − λI) = (−r − λ) N +1 = 0, from which we see that its spectrum is the singleton: Λ G = {−r}. Moreover, since it has algebraic multiplicity a −r = N + 1 and geometric multiplicity g −r = 1, the index of the −r eigenvalue is ν −r = N + 1. Since −r is the only eigenvalue, and all projection operators must sum to the identity, we must have the eigenprojection: G −r = I. The lesson is that the Poisson point process is highly nondiagonalizable.

Homogeneous Poisson processes
When the transition rate r between counter states is constant in time, the net counter state-to-state transition operator from initial time 0 to later time t is given simply by: The functional calculus allows us to directly evaluate e tG for the Poisson nondiagonalizable transition-rate operator G; we find: Consider the orthonormality relation δ i |δ j = δ i,j between counter states, where |δ j is represented by 0s everywhere except for a 1 at counter-state j. It effectively measures the occupation probability of counter-state j. Employing the result for T (t), we find the simple consequence that: That is, the probability that the counter is incremented by n in a time interval t is independent of the initial count and given by: (rt) n e −rt /n!.
Let us emphasize that these steps derived the Poisson distribution directly, rather than as the typical limit of the binomial distribution. Our derivation depended critically on spectral manipulations of a highly nondiagonalizable operator. Moreover, our result for the transition dynamic T (t) allows a direct analysis of how distributions over counts evolve in time, as would be necessary, say, in a Bayesian setting with unknown prior count. This type of calculus can immediately be applied to the analysis of more sophisticated processes, for which we can generally expect nondiagonalizability to play an important functional role.

Inhomogeneous Poisson processes
Let us now generalize to time-inhomogeneous Poisson processes, where the transition rate r between count events is instantaneously uniform, but varies in time as r(t). Conveniently, the associated rate matrices at different times commute with each other. Specifically, with G a = −aI + aD 1 and G b = −bI + bD 1 , we see that: Therefore, the net counter state-to-state transition operator from time t 0 to time t f is given by: where ∆t = t f − t 0 is the time elapsed and: is the average rate during that time. Given Eq. (45), the functional calculus proceeds just as in the timehomogeneous case to give the analogous net transition dynamic: The probability that the count is incremented by n during the time interval ∆t follows directly: With relative ease, our calculus allowed us to derive an important result for stochastic process theory that is nontrivial to derive by other means. Perhaps surprisingly, we see that the probability distribution over final counts induced by any rate trajectory r(t) is the same as if the transition rate were held fixed at mean r throughout the duration. Moreover, we can directly analyze the net evolution of distributions over counts using the derived transition operator T t0,t f .
Note that the nondiagonalizability of the Poisson dynamic is robust in a physical sense. That is, even varying the rate parameter in time in an erratic way, the inherent structure of counting imposes a fundamental nondiagonalizable nature. That nondiagonalizability can be robust in a physical sense is significant, since one might otherwise be tempted to argue that nondiagonalizability is extremely fragile due to numerical perturbations within any matrix representation of the operator. This is simply not the case since such perturbations are physically forbidden. Rather, this simple example challenges us with the fact that some processes, even those familiar and widely used, are intrinsically nondiagonalizable. On the positive side, it appears that spectral methods can now be applied to analyze them. And, this will be particularly important in more complex, memoryful processes [49][50][51], including the hidden semi-Markov processes [45] that are, roughly speaking, the cross-product of hidden finite-state Markov chains and renewal processes.

D. Stochastic thermodynamics
The previous simple examples started to demonstrate the methods of the meromorphic functional calculus. Next, we show a novel application of the meromorphic functional calculus to environmentally driven mesoscopic dynamical systems, selected to give a new set of results within nonequilibrium thermodynamics. In particular, we analyze functions of singular transition-rate operators. Notably, we show that the Drazin inverse arises naturally in the general solution of Green-Kubo relations. We mention that it also arises when analyzing moments of the excess heat produced in the driven transitions atop either equilibrium steady states or nonequilibrium steady states [52].

Dynamics in independent eigenspaces
An important feature of the functional calculus is its ability to address particular eigenspaces independently when necessary. This feature is often taken for granted in the case of normal operators; say, in physical dynamical systems when analyzing stationary distributions or dominant decay modes. Consider a singular operator L that is not necessarily normal and not necessarily diagonalizable and evaluate the simple yet ubiquitous integral τ 0 e tL dt. Via the meromorphic functional calculus we find: where L D is the Drazin inverse of L, discussed earlier.
The pole-pole interaction (z −1 with z −m−1 ) at z = 0 distinguished the 0-eigenspace in the calculations and required the meromorphic functional calculus for direct analysis. The given solution to this integral will be useful in the following.
Next, we consider the case where L is the transitionrate operator among the states of a structured stochastic dynamical system. This leads to several novel consequence within stochastic thermodynamics.

Green-Kubo relations
Let us reconsider the above integral in the case when the singular operator L-let us call it G-is a transitionrate operator that exhibits a single stationary distribution. By the spectral mapping ln Λ e G of the eigenvalue 1 ∈ Λ e G addressed in the Perron-Frobenius theorem, G's zero eigenmode is diagonalizable. And, by assuming a single attracting stationary distribution, the zero eigenvalue has algebraic multiplicity a 0 = 1. Equation (46) then simplifies to: Since G is a transition-rate operator, the above integral corresponds to integrated time evolution. The Drazin inverse G D concentrates on the transient contribution beyond the persistent stationary background. In Eq. (47), the subscript within the left and right eigenvectors explicitly links the eigenvectors to the operator G, to reduce ambiguity. Specifically, the projector |0 G 0 G | maps any distribution to the stationary distribution.
Green-Kubo-type relations [53,54] connect the outof-steady-state transport coefficients to the time integral of steady-state autocorrelation functions. They are thus very useful for understanding out-of-steady-state dissipation due to steady-state fluctuations. (Steady state here refers to either equilibrium or nonequilibrium steady state.) Specifically, the Green-Kubo relation for a transport coefficient, κ say, is typically of the form: where A(0) and A(t) are some observable of the stationary stochastic dynamical system at time 0 and time t, respectively, and the subscript · s.s. emphasizes that the expectation value is to be taken according to the steadystate distribution. Using: the transport coefficient κ can be written more explicitly in terms of the relevant transition-rate operator G for the stochastic dynamics: Thus, we learn that relations of Green-Kubo form are direct signatures of the Drazin inverse of the transitionrate operator for the stochastic dynamic. The result of Eq. (48) holds quite generally. For example, if the steady state has some number of periodic flows, the result of Eq. (48) remains valid. Alternatively, in the case of nonperiodic chaotic flows-where G will be the logarithm of the Ruelle-Frobenius-Perron operator, as described later in § VI F 1-|0 G 0 G | still induces the average over the steady-state trajectories.
In the special case where the transition-rate operator is diagonalizable, − A G D A s.s. is simply the integrated contribution from a weighted sum of decaying exponentials. Transport coefficients then have a solution of the simple form: Note that the minus sign keeps κ positive since Re(λ) < 0 for λ ∈ Λ G \ {0}. Also, recall that G's eigenvalues with nonzero imaginary part occur in complex-conjugate pairs and G λ = G λ . Moreover, if G i,j is the classical transition-rate from state i to state j (to disambiguate from the transposed possibility), then 0 G | is the station-ary distribution. (The latter is sometimes denoted π| in the Markov process literature.) And, |0 G is a column vector of all ones (sometimes denoted |1 ) which acts to integrate contributions throughout the state space. A relationship of the form of Eq. (48), between the Drazin inverse of a classical transition-rate operator and a particular Green-Kubo relation was recently found in Ref. [55] for the friction tensor for smoothly-driven transitions atop nonequilibrium steady states. Subsequently, a truncation of the eigen-expansion of the form of Eq. (49) was recently used in a similar context to bound a universal tradeoff between power, precision, and speed [56]. Equation (48) shows that a fundamental relationship between a physical property and a Drazin inverse is to be expected more generally whenever the property can be related to integrated correlation.
Notably, if a Green-Kubo-like relation integrates a cross-correlation, say between A(t) and B(t) rather than an autocorrelation, then we have only the slight modification: The foregoing analysis bears on both classical and quantum dynamics. G may be a so-called linear superoperator in the quantum regime [57]; for example, the Lindblad superoperator [58,59] that operates on density operators. A Liouville-space representation [60] of the superoperator, though, exposes the superficiality of the distinction between superoperator and operator. At an abstract level, time evolution can be discussed uniformly across subfields and reinterpretations of Eq. (50) will be found in each associated physical theory.
Reference [52] presents additional constructive results that emphasize the ubiquity of integrated correlation and Drazin inverses in the transitions between steady states [61], relevant to the fluctuations within any physical dynamic. Overall, these results support the broader notion that dissipation depends on the structure of correlation.
Frequency-dependent generalizations of integrated correlation have a corresponding general solution. To be slightly less abstract, later on we give novel representative formulae for a particular application: the general solution to power spectra of a process generated by any countable-state hidden Markov chain.

E. Power spectra
A signal's power spectrum quantifies how its power is distributed across frequency [62]. For a discrete-domain process it is: where ω is the angular frequency and X n is the random variable for the observation at time n. For a widesense stationary stochastic process, the power spectrum is also determined from the signal's autocorrelation function γ(τ ): where the autocorrelation function for a wide-sense stationary stochastic process is defined: The windowing function N − |τ | appearing in Eq. (52) is a direct consequence of Eq. (51). It is not imposed externally, as is common practice in signal analysis. This is important to subsequent derivations. The question we address is how to calculate the correlation function and power spectrum given a model of the signal's generator. To this end, we briefly introduce hidden Markov models as signal generators and then use the meromorphic calculus to calculate their autocorrelation and power spectra in closed-form. This leads to several lessons. First, we see that the power spectrum is a direct fingerprint of the resolvent of the generator's time-evolution operator, analyzed along the unit circle. Second, spectrally decomposing the not-necessarilydiagonalizable time evolution operator, we derive the most general qualitative behavior of the autocorrelation function and power spectra. Third, contributions from eigenvalues on the unit circle must be extracted and dealt with separately. Contributions from eigenvalues . 3: Simple 3-state state-emitting HMM that generates a stochastic process according to the state-to-state transition dynamic T and the probability density functions (pdfs) {p(x|s)} s∈S associated with each state. Theorem 1 asserts that its power spectrum will be the same (with only constant offset) as the power spectrum generated from the alternative process where the pdfs in each state are solely concentrated at the Platonic average value x p s (x) of the former pdf associated with the state.
on the unit circle correspond to Dirac delta functionsthe analog of Bragg peaks in diffraction. Whereas, eigen-contributions from inside the unit circle correspond to diffuse peaks, which become sharper for eigenvalues closer to the unit circle. Finally, nondiagonalizable eigenmodes yield qualitatively different line profiles than their diagonalizable counterparts. In short, when applied to signal analysis our generalized spectral decomposition has directly measurable consequences. This has been key to analyzing low-dimensional disordered materials, for example, when adapted to X-ray diffraction spectra [20,21,63]. Let the 4-tuple M = S, A, P, T be some discretetime state-emitting hidden Markov model (HMM) that generates the stationary stochastic process . . . X −2 X −1 X 0 X 1 X 2 . . . according to the following. S is the (finite) set of latent states of the hidden Markov chain and A ⊆ C is the observable alphabet. S t is the random variable for the hidden state at time t that takes on values s ∈ S. X t is the random variable for the observation at time t that takes on values x ∈ A. Given the latent state at time t, the possible observations are distributed according to the conditional probability density functions: P = p(X t = x|S t = s) s∈S . For each s ∈ S, p(X t = x|S t = s) may be abbreviated as p(x|s) since the probability density function in each state is assumed not to change over t. Finally, the latent-state-to-state stochastic transition matrix T has elements T i,j = Pr(S t+1 = s j |S t = s i ), which give the probability of transitioning from latent state s i to s j given that the system is in state s i , where s i , s j ∈ S. It is important for the subsequent derivation that we use Pr(·) to denote a probability in contrast to p(·) which denotes a probability density. The Bayes network diagram of Fig. 2 depicts the structure of conditional independence among the random variables.
1. Continuous-value, discrete-state and -time processes Figure 3 gives a particular HMM with continuous observable alphabet A = R distributed according to the probability density function shown within each latent state. Processes generated as the observation of a function of a Markov chain can be of either finite or infinite Markov order. (They are, in fact, typically infinite Markov order in the space of processes [64].) Directly calculating, one finds that the autocorrelation function, for τ > 0, for any such HMM is: where: holds by definition of conditional probability. The decomposition of: for τ = 0 follows from the conditional independence in the relevant Bayesian network shown in Fig. 2. Moreover, the equality: can be derived by marginalizing over all possible intervening state sequences. Note that |δ s is the column vector of all 0s except for a 1 at the index corresponding to state s and δ s | is simply its transpose. Recall that π| = 1 T | is the stationary distribution induced by T over latent states and |1 = |1 T is a column vector of all ones. Note also that π|δ s = Pr(s) and δ s |1 = 1.
Since the autocorrelation function is symmetric in τ and: we find the full autocorrelation function is given by: where Ω is the |S|-by-|S| matrix defined by: The power spectrum is then calculated via Eq. (52) using the meromorphic calculus. In particular, the power spectrum decomposes naturally into a discrete part and a continuous part. Full details will be given elsewhere, but the derivation is similar to that given in Ref. [20] for the special case of diffraction patterns from HMMs. We note that it is important to treat individual eigenspaces separately, as our generalized calculus naturally accommodates. The end result, for the continuous part of the power spectrum, is: All of the ω-dependence is in the resolvent. Using the spectral expansion of the resolvent given by Eq. (21) allows us to better understand the qualitative possibilities for the shape of the power spectrum: Note that π| Ω T λ,m Ω |1 is a complex-valued scalar and all of the frequency dependence now handily resides in the denominator. The discrete portion (delta functions) of the power spectrum is: where ω λ is related to λ by λ = e iω λ . An extension of the Perron-Frobenius theorem guarantees that the eigenvalues of T on the unit circle have index ν λ = 1.
When plotted as a function of the angular frequency ω around the unit circle, the power spectrum suggestively appears to emanate from the eigenvalues λ ∈ Λ T of the hidden linear dynamic. See Fig. 4 for the analysis of an example parametrized process and the last two panels for this display mode for the power spectra.
Eigenvalues of T on the unit circle yield Dirac delta functions in the power spectrum. Eigenvalues of T within the unit circle yield more diffuse line profiles, increasingly diffuse as the magnitude of the eigenvalues retreats toward the origin. Moreover, the integrated magnitude of each contribution is determined by projecting pairwise observation operators onto the eigenspace emanating the contribution. Finally, we note that nondiagonalizable eigen-modes yield qualitatively different line profiles.
Remarkably, the power spectrum generated by such process is the same as the that generated by a potentially much simpler one-a process that is a function of the same underlying Markov chain but instead emits the state-dependent expectation value of the observable within each state: Theorem 1. Let P = p s (x) s∈S be any set of probability distribution functions over the domain A ⊆ C. Let B = x ps(x) s∈S and let Q = δ(x − x ps(x) ) s∈S . Then, the power spectrum generated by any hidden Markov model M = S, A, P, T differs at most by a constant offset from the power spectrum generated by the hidden Markov model M = S, B, Q, T that has the same latent Markov chain but in any state s ∈ S emits, with probability one, the average value x ps(x) of the state-conditioned probability density function p s (x) ∈ P of M.
This shows that all HMMs that share the same T and x p(x|s) } s∈S have the same power spectrum P (ω) = P c (ω) + P d (ω) besides a constant offset determined by differences in |x| 2 . One immediate consequence is that any hidden Markov chain with any arbitrary set of zero-mean distributions attached to each state, i.e.: P ∈ {p(x|s)} s∈S : x p(x|s) = 0 for all s ∈ S , generates a flat power spectrum with the appearance of white noise. On the one hand, this strongly suggests to data analysts to look beyond power spectra when attempting to extract a process' full architecture. On the other, whenever a process's power spectrum is structured, it is a direct fingerprint of the resolvent of the hidden linear dynamic. In short, the power spectrum is a filtered image of the resolvent along the unit circle.
The power spectrum of a particular stochastic process is shown in Fig. 4 and using coronal spectrograms, introduced in Ref. [20], it illustrates how the observed spectrum can be thought of as emanating from the spectrum of the hidden linear dynamic, as all power spectra must. Figure 4a shows the state-emitting HMM with state-tostate transition probabilities parametrized by b; the mean values x p(x|s) of each state's pdf p(x|s) are indicated as the blue number inside each state. The process generated depends on the actual pdfs and the transition parameter b although, and this is our point, the power spectrum is ignorant to the details of the pdfs. The evolution of the eigenvalues Λ T of the transition dynamic among latent states is shown from thick blue to thin red markers in Fig. 4b, as we sweep the transition parameter b from 1 to 0. A subset of the eigenvalues pass continuously but very quickly through the origin of the complex plane as b passes through 1/2. The continuity of this is not immediately apparent numerically, but can be revealed with a finer increment of b near b ≈ 1/2. Notice the persistent eigenvalue of λ T = 1, which is guaranteed by the Perron-Frobenius theorem.
In Fig. 4c and again, at another parameter setting, in Fig. 4d, we show the continuous part of the power spectrum P c (ω) (plotted around the unit circle in solid blue) and the eigen-spectrum (plotted as red dots on and within the unit circle) of the state-to-state transition matrix for the 11-state hidden Markov chain (leftmost panel) that generates it. There is also a δ-function contribution to the power spectrum at ω = 0 (corresponding to λ T = 1). This is not shown. These coronal spectrograms illustrate how the power spectrum emanates from the HMM's eigen-spectrum, with sharper peaks when the eigenvalues are closer to the unit circle. This observation is fully explained by Eq. (55). The integrated magnitude of each peak depends on π| Ω |λ λ| Ω |1 .
Interestingly, the apparent continuous spectrum component is the shadow of the discrete spectrum of nonunitary dynamics. This suggests that resonances in various physics domains concerned with a continuous spectrum can be modeled as simple consequences of nonunitary dynamics. Indeed, hints of this appear in the literature [65][66][67].

Continuous-time processes
We close this exploration of conventional signal analysis methods using the meromorphic calculus by commenting on continuous-time processes. Analogous formulae can be derived with similar methods for continuous-time hidden Markov jump processes and continuous-time deterministic (possibly chaotic) dynamics in terms of the generator G of time evolution. For example, the continuous part P c (ω) of the power spectrum from a continuoustime deterministic dynamic has the form: Appealing to the resolvent's spectral expansion again allows us to better understand the possible shapes of their power spectrum: Since all of the frequency-dependence has been isolated in the denominator and π| Ω G λ,m Ω |1 is a frequencyindependent complex-valued constant, peaks in P c (ω) can only arise via contributions of the form Re c (iω−λ) n for c ∈ C, ω ∈ R, λ ∈ Λ G , and n ∈ Z + . This provides a rich starting point for application and further theoretical investigation. For example, Eq. (57) helps explain the shapes of power spectra of nonlinear dynamical systems, as have appeared, e.g., in Ref. [68]. Furthermore, it suggests an approach to the inverse problem of inferring the spectrum of the hidden linear dynamic via power spectra.
In the next section, however, we develop a more general proposal for inferring eigenvalues from a time series. Fur-Since trajectories in state-space can be generated independently of each other, any nonlinear dynamic corresponds to a linear operation on an infinite-dimensional vector-space of complex-valued distributions (in the sense of generalized functions) over the original state-space. For example, the well known Lorenz ordinary differential equations [69] are nonlinear over its three given statespace variables-x, y, and z. Nevertheless, the dynamic is linear in the infinite-dimensional vector space D(R 3 ) of distributions over R 3 . Although D(R 3 ) is an unwieldy state-space, the dynamics there might be well approximated by a finite truncation of its modes.

Ruelle-Frobenius-Perron and Koopman operators
The preceding operator formalism applies, in principle at least. The question, of course, is, Is it practical and does it lead to constructive consequences? Let's see. The right eigenvector is either |0 G or |1 T with T = e τ G as the Ruelle-Frobenius-Perron transition operator [70,71]. Equivalently, it is also π, the stationary distribution, with support on attracting subsets of R 3 in the case of the Lorenz dynamic. The corresponding lefteigenvector 1, either 0 G | or 1 T |, is uniform over the space. Other modes of the operator's action, according to the eigenvalues and left and right eigenvectors and generalized eigenvectors, capture the decay of arbitrary distributions on R 3 .
The meromorphic spectral methods developed above give a view of the Koopman operator and Koopman modes of nominally nonlinear dynamical systems [4] that is complementary to the Ruelle-Frobenius-Perron operator. The Koopman operator K is the adjoint-in the sense of vector spaces, not inner product spaces-of the Ruelle-Frobenius-Perron operator T : effectively the transpose K = T . Moreover, it has the same spectrum with only right and left swapping of the eigenvectors and generalized eigenvectors.
The Ruelle-Frobenius-Perron operator T is usually associated with the evolution of probability density, while the Koopman operator K is usually associated with the evolution of linear functionals of probability density. The duality of perspectives is associative in nature: f | T n |ρ 0 corresponds to the Ruelle-Frobenius-Perron perspective with T acting on the density ρ and Both approaches suffer when their operators are defective. Given the meromorphic calculus' ability to work around a wide class of such defects, adapting it the Ruelle-Frobenius-Perron and Koopman operators suggests that it may lift their decades-long restriction to only analyzing highly idealized (e.g., hyperbolic) chaotic systems.

Eigenvalues from a time series
Let's explore an additional benefit of this view of the Ruelle-Frobenius-Perron and Koopman operators, by proposing a novel method to extract the eigenvalues of a nominally nonlinear dynamic. Let O N (f, z) be (z −1 times) the z-transform [72, pp. 257-262] of a length-N sequence of τ -spaced type-f observations of a dynamical system: as N → ∞ for |z| = r > 1. Note that f |T n |ρ 0 is simply the f -observation of the system at time nτ , when the system started in state ρ 0 . We see that this z-transform of observations automatically induces the resolvent of the hidden linear dynamic. If the process is continuous-time, then T = e τ G implies λ T = e τ λ G , so that the eigenvalues should shift along the unit circle if τ changes; but the eigenvalues should be invariant to τ in the appropriate τ -dependent conformal mapping of the inside of the unit circle of the complex plane to the left half complex plane. Specifically, for any experimentally accessible choice of inter-measurement temporal spacing τ , the fundamental set of continuous time eigenvalues Λ G can be obtained from λ G = 1 τ ln λ T , where each λ T ∈ Λ T is extrapolated from c/(re iω − λ T ) n curves fit to O N (f, re iω ) for c ∈ C, large N , and fixed r.
The square magnitude of O N (f, z) is related to the power spectrum generated by f -type observations of the system. Indeed, the power spectrum generated by any type of observation of a nominally nonlinear system is a direct fingerprint of the eigenspectrum and resolvent of the hidden linear dynamic. This suggests many opportunities for inferring eigenvalues and projection operators from frequency-domain transformations of a time series.

VII. CONCLUSION
The original, abstract spectral theory of normal operators rose to central importance when, in the early development of quantum mechanics, the eigenvalues of Hermitian operators were detected experimentally in the optical spectra of energetic transitions of excited electrons. We extended this powerful theory by introducing the meromorphic functional calculus, providing the spectral theory of nonnormal operators. Our straightforward examples suggest that the spectral properties of these general operators should also be experimentally accessible in the behavior of complex-open, strongly interactingsystems. We see a direct parallel with the success of the original spectral theory of normal operators as it made accessible the phenomena of the quantum mechanics of closed systems. This turns on nondiagonalizability and appreciating how ubiquitous it is.
Nondiagonalizability has consequences for settings as simple as counting, as shown in § VI C. Moreover, there we found that nondiagonalizability can be robust. The Drazin inverse, the negative-one power in the meromorphic functional calculus, is quite common in the nonequilibrium thermodynamics of open systems, as we showed in § VI D. Finally, we showed that the spectral character of nonnormal and nondiagonalizable operators manifests itself physically, as illustrated by Figs. 4c and 4d of § VI E.
From the perspective of functional calculus, nonunitary time evolution, open systems, and non-Hermitian generators are closely related concepts since they all rely on the manipulation of nonnormal operators. Moreover, each domain is gaining traction. Nonnormal operators have recently drawn attention, from the nonequilibrium thermodynamics of nanoscale systems [73] to largescale cosmological evolution [74]. In another arena entirely, complex directed networks [75] correspond to nonnormal and not-necessarily-diagonalizable weighted digraphs. There are even hints that nondiagonalizable network structures can be optimal for implementing certain dynamical functionality [76]. The opportunity here should be contrasted with the well established field of spectral graph theory [77] that typically considers consequences of the spectral theorem for normal operators applied to the symmetric (and thus normal) adjacency matrices and Laplacian matrices. It seems that the meromorphic calculus and its generalized spectral theory will enable a spectral weighted digraph theory beyond the purview of current spectral graph theory.
Even if the underlying dynamic is diagonalizable, particular questions or particular choices of observable often induce a nondiagonalizable hidden linear dynamic. The examples already showed this arising from the simple imposition of counting or assuming a Poissonian dynamic. In more sophisticated examples, we recently found nondiagonalizable dynamic structures in quantum memory reduction [22] and classical complexity measures [19].
Our goal has been to develop tractable, exact analytical techniques for nondiagonalizable systems. We did not discuss numerical implementation of algorithms that naturally accompany its practical application. Nevertheless, the theory does suggest new algorithms-for the Drazin inverse, projection operators, power spectra, and more. Guided by the meromorphic calculus, such algorithms can be made robust despite the common knowledge that numerics with nondiagonalizable matrices is sensitive in certain ways.
The meromorphic calculus complements attempts to address nondiagonalizability, e.g., via pseudospectra [78,79]. It also extends and simplifies previously known results, especially as developed by Dunford [16]. Just as the spectral theorem for normal operators enabled much theoretical progress in physics, we hope that our generalized and tractable analytic framework yields rigorous understanding for much broader classes of complex system. Importantly, the analytic framework should enable new theory of complex systems beyond the limited purview of numerical investigations.
While the infinite-dimensional theory is in principle readily adaptable from the present framework, special care must be taken to guarantee a similar level of tractability and generality. Nevertheless, even the finitedimensional theory enables a new level of tractability for analyzing not-necessarily-diagonalizable systems, including nonnormal dynamics. Future work will take full advantage of the operator theory, with more emphasis on infinite-dimensional systems and continuous spectra. Another direction forward is to develop creation and annihilation operators within nondiagonalizable dynamics. In the study of complex stochastic information processing, for example, this would allow analytic study of infinitememory processes generated by, say, stochastic pushdown and counter automata [51,[80][81][82]. In a physical context, such operators may aid in the study of open quantum field theories. One might finally speculate that the Drazin inverse will help to tame the divergences that arise there.