def symmetric_rescaling(A, copy=True):
    """Scale the matrix symmetrically.

        A = D^{-1/2} A D^{-1/2}

    where D=diag(A).

    The left multiplication is accomplished through scale_rows and the right
    multiplication is done through scale columns.

    Parameters
    ----------
    A : sparse matrix
        Sparse matrix with N rows
    copy : {True,False}
        - If copy=True, then the matrix is copied to a new and different return
          matrix (e.g. B=symmetric_rescaling(A))
        - If copy=False, then the matrix is overwritten deeply (e.g.
          symmetric_rescaling(A,copy=False) overwrites A)

    Returns
    -------
    D_sqrt : array
        Array of sqrt(diag(A))
    D_sqrt_inv : array
        Array of 1/sqrt(diag(A))
    DAD    : csr_matrix
        Symmetrically scaled A

    Notes
    -----
    - if A is not csr, it is converted to csr and sent to scale_rows

    Examples
    --------
    >>> import numpy as np
    >>> from scipy.sparse import spdiags
    >>> from pyamg.util.utils import symmetric_rescaling
    >>> n=5
    >>> e = np.ones((n,1)).ravel()
    >>> data = [ -1*e, 2*e, -1*e ]
    >>> A = spdiags(data,[-1,0,1],n,n).tocsr()
    >>> Ds, Dsi, DAD = symmetric_rescaling(A)
    >>> print DAD.todense()
    [[ 1.  -0.5  0.   0.   0. ]
     [-0.5  1.  -0.5  0.   0. ]
     [ 0.  -0.5  1.  -0.5  0. ]
     [ 0.   0.  -0.5  1.  -0.5]
     [ 0.   0.   0.  -0.5  1. ]]

    """
    if isspmatrix_csr(A) or isspmatrix_csc(A) or isspmatrix_bsr(A):
        if A.shape[0] != A.shape[1]:
            raise ValueError('expected square matrix')

        D = diag_sparse(A)
        mask = (D != 0)

        if A.dtype != complex:
            D_sqrt = np.sqrt(abs(D))
        else:
            # We can take square roots of negative numbers
            D_sqrt = np.sqrt(D)

        D_sqrt_inv = np.zeros_like(D_sqrt)
        D_sqrt_inv[mask] = 1.0/D_sqrt[mask]

        DAD = scale_rows(A, D_sqrt_inv, copy=copy)
        DAD = scale_columns(DAD, D_sqrt_inv, copy=False)

        return D_sqrt, D_sqrt_inv, DAD

    else:
        return symmetric_rescaling(csr_matrix(A))
Scale the matrix symmetrically.

        A = D^{-1/2} A D^{-1/2}

    where D=diag(A).

    The left multiplication is accomplished through scale_rows and the right
    multiplication is done through scale columns.

    Parameters
    ----------
    A : sparse matrix
        Sparse matrix with N rows
    copy : {True,False}
        - If copy=True, then the matrix is copied to a new and different return
          matrix (e.g. B=symmetric_rescaling(A))
        - If copy=False, then the matrix is overwritten deeply (e.g.
          symmetric_rescaling(A,copy=False) overwrites A)

    Returns
    -------
    D_sqrt : array
        Array of sqrt(diag(A))
    D_sqrt_inv : array
        Array of 1/sqrt(diag(A))
    DAD    : csr_matrix
        Symmetrically scaled A

    Notes
    -----
    - if A is not csr, it is converted to csr and sent to scale_rows

    Examples
    --------
    >>> import numpy as np
    >>> from scipy.sparse import spdiags
    >>> from pyamg.util.utils import symmetric_rescaling
    >>> n=5
    >>> e = np.ones((n,1)).ravel()
    >>> data = [ -1*e, 2*e, -1*e ]
    >>> A = spdiags(data,[-1,0,1],n,n).tocsr()
    >>> Ds, Dsi, DAD = symmetric_rescaling(A)
    >>> print DAD.todense()
    [[ 1.  -0.5  0.   0.   0. ]
     [-0.5  1.  -0.5  0.   0. ]
     [ 0.  -0.5  1.  -0.5  0. ]
     [ 0.   0.  -0.5  1.  -0.5]
     [ 0.   0.   0.  -0.5  1. ]]