0
$\begingroup$

In the mathematical field of linear algebra and convex analysis, the numerical range or field of values of a complex

$n\times n$ matrix $A$ is the set

$$W(A)=\left\{{\frac {{\mathbf {x}}^{*}A{\mathbf {x}}}{{\mathbf {x}}^{*}{\mathbf {x}}}}\mid {\mathbf {x}}\in {\mathbb {C}}^{n},\ x\not =0\right\}$$

where $\mathbf{x}^*$ denotes the conjugate transpose of the vector $\mathbf {x}$ .

In engineering, numerical ranges are used as a rough estimate of eigenvalues of A. Recently, generalizations of the numerical range are used to study quantum computing.

A related concept is the numerical radius, which is the largest absolute value of the numbers in the numerical range, i.e.

$$r(A)=\sup\{|\lambda |:\lambda \in W(A)\}=\sup _{{\|x\|=1}}|\langle Ax,x\rangle |.$$

I want to implement the numerical range of a matrix in Python. I found the following code but contains a function command mlab.find which no longer exists. I also found how to replace it but the output is not the desired one. Any help?

import numpy as np
import matplotlib.pyplot as plt


def find(condition):
    res,  =  np.nonzero(np.ravel(condition))
    return res


def numerical_range(A, resolution=0.01):
    A = np.asmatrix(A)
    th = np.arange(0, 2*np.pi + resolution, resolution)
    k = 0
    w = []
    for j in th:
        Ath = np.exp(-1j*j)*A
        Hth = (Ath + Ath.H)/2
        e,r = np.linalg.eigh(Hth)
        r = np.matrix(r)
        e = np.real(e)
        m = e.max()
        s = find(e == m)
        if np.size(s) == 1:
            w.append(np.matrix.item(r[:,s].H*A*r[:,s]))
        else:
            Kth = 1j*(Hth - Ath)
            pKp = r[:,s].H*Kth*r[:,s]
            ee,rr = np.linalg.eigh(pKp)
            rr = np.matrix(rr)
            ee = np.real(ee)
            mm = ee.min()
            sm = find(ee == mm)
            temp = rr[:,sm[0]].H*r[:,s].H*A*r[:,s]*rr[:,sm[0]]
            w.append(temp[0,0])
            k += 1
            mM = ee.max()
            sM = find(ee == mM)
            temp = rr[:,sM[0]].H*r[:,s].H*A*r[:,s]*rr[:,sM[0]]
            w.append(temp[0,0])
        k += 1
    plt.plot(w)
    return None


A = np.matrix([[.6,.5], [.6,.5]])
numerical_range(A, resolution=0.01)
plt.show()
```
$\endgroup$
9
  • $\begingroup$ Welcome to SciComp.SE. As it is right now your question is closer to a programming question than a computational Science one. Does "range" the to the "rank" of the Matrix? If that's the case, why but using numpy.linalg.matrix_rank? Also, what does the name mlab refers to? $\endgroup$ – nicoguaro Apr 1 at 5:21
  • $\begingroup$ Mlab is from Matlab.it refers to the spectrum of a matrix which is an ellipse with foci the eigenvalues of the matrix $\endgroup$ – Tsiantakis Apr 1 at 5:41
  • 1
    $\begingroup$ Separately, questions of the form "the output is not the desired one. Any help?" after showing 50 lines of undocumented code are quite unpopular here. We have no idea what your code is supposed to do, and we have no idea what you think the output should be. You are unlikely going to get answers without providing us a bit more information. $\endgroup$ – Wolfgang Bangerth Apr 1 at 13:40
  • 1
    $\begingroup$ Oh, I see mlab is the deprecated module in matplotlib, I think that you could use np.argwhere for that. Regarding your question, I still don't get what the range is, but maybe somebody does. It seems to be a set of scalars, although your definition does not describe how to choose $\mathbf{x}$. $\endgroup$ – nicoguaro Apr 1 at 15:12
  • 2
    $\begingroup$ @nicoguaro The numerical range of a matrix $A$ is meant to be a better way to assess the stability of the dynamical system $\dot x = Ax$ that works just as well for aggressively non-normal or non-diagonalizable matrices as it does for normal matrices. The motivation for it is that the usual eigenvalue stability analysis that everyone gets taught in undergrad can predict that a system is stable at $t = \infty$ but fails to show that the system might exhibit huge departures in finite time. See also the $\epsilon$-pseudospectrum, field of values, etc. $\endgroup$ – Daniel Shapero Apr 1 at 16:51

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.