Beyond the spectral theorem: Spectrally decomposing arbitrary functions of nondiagonalizable operators

Decomposing a complicated system into its constituent parts—reductionism—is one of science’s most powerful strategies for analysis and understanding. Large-scale 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 from 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 weighting 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 high-dimensional spaces.^2–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.

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.⁵ 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 eigenvalue spectrum 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.⁶ 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 spatially-inhomogeneous distribution of energy eigenvalues and the occupation of their corresponding states.⁷ Before their experimental discovery, anti-particles were anticipated as the nonoccupancy of negative-energy eigenstates of the Dirac Hamiltonian.⁸ 

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.⁴ A nondynamic version of rendering nonlinearities into linearities in a higher-dimensional feature space is exploited with much success today in machine learning by support vector machines, for example Ref. 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”.¹⁰ 

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,¹¹ 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,¹⁵ 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 complex-conjugate 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 by, say, 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 novel physics governed by nonnormal and not-necessarily-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,¹⁶ 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.

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 operator’s eigenvalues. 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 and elsewhere. Furthermore, we derive consequences from the more familiar holomorphic functional calculus that readily allow spectral decomposition of nondiagonalizable operators in terms of spectral projections and left and right generalized eigenvectors—decanting the abstract mathematical theory into a more tractable framework for analyzing complex physical systems.

Taken together, 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 in several specific applications, including closed-form expressions for signal processing and information measures of hidden Markov processes,^19–23 compressing stochastic processes over a quantum channel,^24,25 and stochastic thermodynamics.^26,27 However, the techniques are sufficiently general they will be much more widely useful. We envision 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 Sec. III through Sec. IV, after necessary preparation in Sec. II. Section V A further explores and gives a new formula for eigenprojectors, which we refer to here simply as projection operators. Section V B makes explicit their general relationship with eigenvectors and generalized eigenvectors and clarifies the orthonormality relationship among left and right generalized eigenvectors. Section V B 4 then discusses simplifications of the functional calculus for special cases, while Sec. VI A takes up the spectral properties of transition operators.

Section VI then turns to several applications that demonstrate the theoretical results’ broad utility. In particular, the advantage of the meromorphic approach over holomorphic is demonstrated in deriving Eq. (46) of Sec. VI D on stochastic thermodynamics. This enables a general operator approach to analyzing the nonequilibrium thermodynamics of complex systems. Generally, though, each application benefits from the meromorphic approach via the analytical forms derived from it. Specifically, the meromorphic calculus allows one to derive: (i) arbitrary powers of an operator (Eq. (25)), explicitly including the qualitatively distinct contribution of the zero eigenspace; (ii) negative powers (Eq. (28)) and the Drazin inverse (Eq. (29)); and (iii) key formulae for the spectral projection operators, including Eq. (35). These results make the applications tractable. In this way, the applications demonstrate how the meromorphic toolset is used to analyze linear operators in a range of settings. The result is an efficient and intuitive analysis of the structure and randomness in complex systems. Finally, Sec. VII closes with suggestions for future applications and research directions.

The following is relatively self-contained, assuming basic familiarity with linear algebra at the level of Refs. 15 and 17—including eigen-decomposition and knowledge of the Jordan canonical form, partial fraction expansion (see Ref. 28), and series expansion—and basic knowledge of complex analysis—including the residue theorem and calculation of residues at the level of Ref. 29. For those lacking a working facility with these concepts, a quick review of Sec. 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 repeatedly 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.³⁰ 

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:³³ 

To fully appreciate the meromorphic functional calculus, we first state and compare the main features and limitations of alternative functional calculi.

Inspired by the Taylor expansion of scalar functions:

a calculus for functions of an operator A can be based on the series:

where f⁽ⁿ⁾(ξ) 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.

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³⁰ 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. 30 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 a 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. 34 and 35 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.

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.

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.

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 $∑λ∈ΛAaλ$⁠.

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 λ ∈ Λ.²⁸ 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 as the residues:

where the residues are calculated elementwise. The projection operators A_λ associated with each eigenvalue λ were already referenced in Sec. 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.

Even the simplest applications of Eq. (13) yield insight. Consider the identity as the operator function f(A) = A⁰ = I that corresponds to the scalar function f(z) = z⁰ = 1. Then, Eq. (13) implies:

This shows that the projection operators are, in fact, a decomposition of the identity, as anticipated in Eq. (2).

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:¹⁶ A = D + N, where $D≡∑λ∈ΛAλAλ$ is diagonalizable, $N≡∑λ∈ΛANλ$ is nilpotent, and D and N commute: [D, N] = 0. This is also known as the Jordan–Chevalley decomposition.

The special case where A is diagonalizable implies that N = 0. And so, Eq. (14) simplifies to:

As shown in Ref. 14 and can be derived from Eqs. (12) and (16):

and

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.

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 Sec. 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.

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:

Equation (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⁰ and A¹. 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 $Lm$ 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.

The eigen-decomposition of polynomials implied by Eq. (25) makes the contribution of the zero eigenvalue more explicit than previous treatments and enables closed-form expressions, e.g., for correlation functions, where the zero eigenvalue makes a qualitatively distinct contribution.^21,22 Consequentially, this formulation can lead to the recognition of coexistent finite and infinite range physical phenomena of different mechanistic origin.^23,24

If L is a nonnegative integer such that L ≥ ν_λ − 1 for all λ ∈ Λ_A, then:

where $Lm$ is now reduced to the traditional binomial coefficient L!/(m!(L − m)!).

If L is any negative integer, then $−|L|m$ can be written as a traditional binomial coefficient $(−1)mL+m−1m$⁠, yielding:

for $−L∈{−1,−2,−3,…}$⁠.

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⁻¹ 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⁻¹ = I − A₀.

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$AD$ of A. Note that it is not the same as the Moore–Penrose pseudo-inverse.^36,37

Although the Drazin inverse is usually defined axiomatically to satisfy certain criteria,³⁸ 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.³⁸ The analytical form of Eq. (29) has been teased out previously by other means; see, e.g., Ref. 38 and for other settings see Refs. 39 and 40. Nevertheless, due to its utility in application, it is noteworthy and appealing that the Drazin inverse falls out organically in the meromorphic functional calculus, as the negative-one power, in contrast to its otherwise rather esoteric axiomatic origin.

While A⁻¹ 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⁻¹. 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⁻¹ 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⁻¹ with inv(A), we use the notation $AD$ for the Drazin inverse of A. Still, $AD=inv(A)$⁠, whenever 0 ∉ Λ_A.

Amusingly, this extension of previous calculi lets us resolve an elementary but fundamental question: What is 0⁻¹? It is certainly not infinity. Indeed, it is just as close to negative infinity! Rather: 0⁻¹ = 0 ≠ inv(0).

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₀ exists within the finite ϵ-ball about the origin.

Third, we note that A + cA₀ is invertible for all c ≠ 0; this can be proven by multiplying each side of Eq. (31) by A + cA₀. Hence, $(A+cA0)−1=inv(A+cA0)$ for all c ≠ 0.

Finally, multiplying each side of Eq. (31) by I − A₀, and recalling that A_0,0A_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₀. Specifically:

which is valid for any $c∈C\{0}$⁠. Eq. (32) generalizes the result found specifically for c = −1 in Ref. 41.

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. Spectral mapping of this result means that T’s eigenvalue 1 maps to the eigenvalue 0 of I − T and T₁ = (I − T)₀. Moreover:

since ν₀ = 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+T1)−1$ is called the fundamental matrix in that setting.⁴² 

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 Sec. VII of Ref. 43. Moreover, Eq. (34) allows an extension of the functional calculus of Refs. 34, 35, and 44, 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:

• Loss of homomorphism: f₁(A)f₂(A) ≠ ( f₁ · f₂)(A);

• Loss of naive spectral mapping: f(Λ_A∖Ξ_f) ⊂ Λ_f(A).

A simple example of both losses arises with the Drazin inverse, above. There, f₁(z) = z⁻¹. Taking this and f₂(z) = z combined with singular operator A leads to the loss of homomorphism: $ADA≠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⁻¹, how does A’s eigenvalue of 0 get mapped into the new spectrum of $AD$? A naive application of the spectral mapping theorem might seem to yield an undefined quantity. But, using the meromorphic functional calculus self-consistently maps the eigenvalue as 0⁻¹ = 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.

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 develop 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 show how the general results specialize in several common cases of interest. After these, we turn to examples and applications.

To obtain the projection operators associated with each index-one eigenvalue λ ∈ {ζ ∈ Λ_A:ν_ζ = 1}, we apply the functional 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.

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.

Unitary, Hermitian, skew-Hermitian, orthogonal, symmetric, and skew-symmetric matrices are all special cases of normal matrices. As noted, normal matrices are those that commute with their Hermitian adjoint (complex-conjugate transpose): AA^† = A^†A. Moreover, a matrix is normal if and only if it can be diagonalized by a unitary transformation: A = UΛU^†, where the columns of the unitary matrix U are the orthonormal right eigenvectors of A corresponding to the eigenvalues ordered along the diagonal matrix Λ. For an M-by-M matrix A, the eigenvalues in Λ_A are ordered and enumerated according to the possibly degenerate M-tuple (Λ_A) = (λ₁, …, λ_M). 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 right-eigenvectors 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 only applies 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⁻¹ = U^† are A’s left eigenvectors. In this case, they are simply the conjugate transpose of the right eigenvectors. Note that conjugate transposition is the familiar transformation rule between ket and bra spaces in quantum mechanics (see, e.g., Ref. 45)—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.¹⁵ Implicitly in quantum mechanics, the basis is taken to be the basis of eigenstates of any Hermitian operator, supposedly 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 dot-product 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 inner-product-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.

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⁻¹, where the columns of the square matrix P are the not-necessarily-orthogonal right eigenvectors of A corresponding to the eigenvalues ordered along the diagonal matrix Λ and where the rows of P⁻¹ 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.

Given the M-tuple of possibly-degenerate eigenvalues (Λ_A) = (λ₁, λ₂, …, λ_M), there is a corresponding M-tuple of linearly-independent right-eigenvectors (|λ₁⟩, |λ₂⟩, …, |λ_M⟩) and a corresponding M-tuple of linearly-independent left-eigenvectors (⟨λ₁|, ⟨λ₂|, …, ⟨λ_M|) such that:

and:

with the orthonormality condition that:

To avoid misinterpretation, we stress that the bras and kets that appear above are the left and right eigenvectors, respectively, and typically do not correspond to complex-conjugate transposition.

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.

Not all matrices can be diagonalized, but all square matrices can be put into Jordan canonical form via the transformation: A = Y JY⁻¹.¹⁷ 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⁻¹ 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₁, J₂, …, J_n), with 1 ≤ n ≤ M. The k^th Jordan block J_k has dimension m_k-by-m_k:

such that: $∑k=1nmk=M$⁠.

Note that eigenvalue λ ∈ Λ_A corresponds to g_λ different Jordan blocks, where g_λ is the geometric multiplicity of the eigenvalue λ. Indeed: $n=∑λ∈ΛAgλ$⁠. 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 M-dimensional vector space. However, the set of M generalized eigenvectors does form a basis for the vector space.⁴⁶ 

Given the n-tuple of possibly-degenerate eigenvalues (Λ_A) = (λ₁, λ₂, …, λ_n) where each eigenvalue λ ∈ Λ_A is listed g_λ times, there is a corresponding n-tuple of m_k-tuples of linearly-independent generalized right-eigenvectors:

where, for each λ_k ∈ (Λ_A):

and a corresponding n-tuple of m_k-tuples of linearly-independent generalized left-eigenvectors:

where:

such that:

and:

for 0 ≤ m ≤ m_k − 1, where $|λj(0)〉=0→$ and $〈λj(0)|=0→$⁠. Specifically, $|λk(1)〉$ and $〈λk(1)|$ are conventional right and left eigenvectors, respectively.