Ticket #10848: trac_10848-hermitian-matrices-v4.patch

File trac_10848-hermitian-matrices-v4.patch, 7.0 KB (added by rbeezer, 9 years ago)
  • sage/matrix/matrix0.pyx

    # HG changeset patch
    # User Rob Beezer <beezer@ups.edu>
    # Date 1298574394 28800
    # Node ID 07a2e975404a08cddf524ef0ad6ac59217519672
    # Parent  aff88d266c1d2bca6b7be97d579e198ebf466863
    10848: add is_hermitian for matrices
    
    diff -r aff88d266c1d -r 07a2e975404a sage/matrix/matrix0.pyx
    a b  
    28962896                    return False
    28972897        return True
    28982898
     2899    def is_hermitian(self):
     2900        r"""
     2901        Returns ``True`` if the matrix is equal to its conjugate-transpose.
     2902
     2903        OUTPUT:
     2904
     2905        ``True`` if the matrix is square and equal to the transpose
     2906        with every entry conjugated, and ``False`` otherwise.
     2907
     2908        Note that if conjugation has no effect on elements of the base
     2909        ring (such as for integers), then the :meth:`is_symmetric`
     2910        method is equivalent and faster.
     2911
     2912        This routine is for matrices over exact rings and so may not
     2913        work properly for matrices over ``RR`` or ``CC``.  For matrices with
     2914        approximate entries, the rings of double-precision floating-point
     2915        numbers, ``RDF`` and ``CDF``, are a better choice since the
     2916        :meth:`sage.matrix.matrix_double_dense.Matrix_double_dense.is_hermitian`
     2917        method has a tolerance parameter.  This provides control over
     2918        allowing for minor discrepancies between entries when checking
     2919        equality.
     2920
     2921        EXAMPLES::
     2922
     2923            sage: A = matrix(QQbar, [[ 1 + I,  1 - 6*I, -1 - I],
     2924            ...                      [-3 - I,     -4*I,     -2],
     2925            ...                      [-1 + I, -2 - 8*I,  2 + I]])
     2926            sage: A.is_hermitian()
     2927            False
     2928            sage: B = A*A.conjugate_transpose()
     2929            sage: B.is_hermitian()
     2930            True
     2931
     2932        Sage has several fields besides the entire complex numbers
     2933        where conjugation is non-trivial. ::
     2934
     2935            sage: F.<b> = QuadraticField(-7)
     2936            sage: C = matrix(F, [[-2*b - 3,  7*b - 6, -b + 3],
     2937            ...                  [-2*b - 3, -3*b + 2,   -2*b],
     2938            ...                  [   b + 1,        0,     -2]])
     2939            sage: C.is_hermitian()
     2940            False
     2941            sage: C = C*C.conjugate_transpose()
     2942            sage: C.is_hermitian()
     2943            True
     2944
     2945        A matrix that is nearly Hermitian, but for a non-real
     2946        diagonal entry. ::
     2947
     2948            sage: A = matrix(QQbar, [[    2,   2-I, 1+4*I],
     2949            ...                      [  2+I,   3+I, 2-6*I],
     2950            ...                      [1-4*I, 2+6*I,     5]])
     2951            sage: A.is_hermitian()
     2952            False
     2953            sage: A[1,1] = 132
     2954            sage: A.is_hermitian()
     2955            True
     2956
     2957        Rectangular matrices are never Hermitian.  ::
     2958
     2959            sage: A = matrix(QQbar, 3, 4)
     2960            sage: A.is_hermitian()
     2961            False
     2962        """
     2963        if not self.is_square():
     2964            return False
     2965
     2966        cdef Py_ssize_t i,j
     2967
     2968        for i from 0 <= i < self._nrows:
     2969            for j from 0 <= j <= i:
     2970                if self.get_unsafe(i,j) != self.get_unsafe(j,i).conjugate():
     2971                    return False
     2972        return True
     2973
    28992974    def is_skew_symmetric(self):
    29002975        """
    29012976        Returns True if this is a skew symmetric matrix.
  • sage/matrix/matrix_double_dense.pyx

    diff -r aff88d266c1d -r 07a2e975404a sage/matrix/matrix_double_dense.pyx
    a b  
    13341334        self.cache(key, b)
    13351335        return b
    13361336
     1337    def is_hermitian(self, tol = 1e-12):
     1338        r"""
     1339        Return ``True`` if the matrix is equal to its conjugate-transpose.
     1340
     1341        INPUT:
     1342
     1343        - ``tol`` - default: ``1e-12`` - the largest value of the
     1344          absolute value of the difference between two matrix entries
     1345          for which they will still be considered equal.
     1346
     1347        OUTPUT:
     1348
     1349        ``True`` if the matrix is square and equal to the transpose
     1350        with every entry conjugated, and ``False`` otherwise.
     1351
     1352        Note that if conjugation has no effect on elements of the base
     1353        ring (such as for integers), then the :meth:`is_symmetric`
     1354        method is equivalent and faster.
     1355
     1356        The tolerance parameter is used to allow for numerical values
     1357        to be equal if there is a slight difference due to round-off
     1358        and other imprecisions.
     1359
     1360        The result is cached, on a per-tolerance basis.
     1361
     1362        EXAMPLES::
     1363
     1364            sage: A = matrix(CDF, [[ 1 + I,  1 - 6*I, -1 - I],
     1365            ...                    [-3 - I,     -4*I,     -2],
     1366            ...                    [-1 + I, -2 - 8*I,  2 + I]])
     1367            sage: A.is_hermitian()
     1368            False
     1369            sage: B = A*A.conjugate_transpose()
     1370            sage: B.is_hermitian()
     1371            True
     1372
     1373        We get a unitary matrix from the SVD routine and use this
     1374        numerical matrix to create a matrix that should be Hermitian
     1375        (indeed it should be the identity matrix), but with some
     1376        imprecision.  We use this to illustrate that if the tolerance
     1377        is set too small, then we can be too strict about the equality
     1378        of entries and achieve the wrong result.  ::
     1379
     1380            sage: A = matrix(CDF, [[ 1 + I,  1 - 6*I, -1 - I],
     1381            ...                    [-3 - I,     -4*I,     -2],
     1382            ...                    [-1 + I, -2 - 8*I,  2 + I]])
     1383            sage: U, _, _ = A.SVD()
     1384            sage: B=U*U.conjugate_transpose()
     1385            sage: B.is_hermitian()
     1386            True
     1387            sage: B.is_hermitian(tol=5.0e-17)
     1388            True
     1389            sage: B.is_hermitian(tol=3.0e-17)
     1390            False
     1391
     1392        A matrix that is nearly Hermitian, but for a non-real
     1393        diagonal entry. ::
     1394
     1395            sage: A = matrix(CDF, [[    2,   2-I, 1+4*I],
     1396            ...                    [  2+I,   3+I, 2-6*I],
     1397            ...                    [1-4*I, 2+6*I,     5]])
     1398            sage: A.is_hermitian()
     1399            False
     1400            sage: A[1,1] = 132
     1401            sage: A.is_hermitian()
     1402            True
     1403
     1404        Rectangular matrices are never Hermitian.  ::
     1405
     1406            sage: A = matrix(CDF, 3, 4)
     1407            sage: A.is_hermitian()
     1408            False
     1409        """
     1410        global numpy
     1411        tol = float(tol)
     1412        key = 'hermitian_%s'%tol
     1413        b = self.fetch(key)
     1414        if not b is None:
     1415            return b
     1416        if not self.is_square():
     1417            self.cache(key, False)
     1418            return False
     1419        if numpy is None:
     1420            import numpy
     1421        cdef Py_ssize_t i, j
     1422        hermitian = True
     1423        for i from 0 < i < self._nrows:
     1424            for j from 0 <= j <= i:
     1425                if numpy.absolute(self.get_unsafe(i,j) - self.get_unsafe(j,i).conjugate()) > tol:
     1426                    hermitian = False
     1427                    break
     1428        self.cache(key, hermitian)
     1429        return hermitian
     1430
    13371431    def schur(self, base_ring=None):
    13381432        r"""
    13391433        Returns the Schur decomposition of the matrix.