def pinv_array(a, cond=None):
    """Calculate the Moore-Penrose pseudo inverse of each block of the three dimensional array a.

    Parameters
    ----------
    a   : {dense array}
        Is of size (n, m, m)
    cond : {float}
        Used by gelss to filter numerically zeros singular values.
        If None, a suitable value is chosen for you.

    Returns
    -------
    Nothing, a is modified in place so that a[k] holds the pseudoinverse
    of that block.

    Notes
    -----
    By using lapack wrappers, this can be much faster for large n, than
    directly calling pinv2

    Examples
    --------
    >>> import numpy as np
    >>> from pyamg.util.linalg import pinv_array
    >>> a = np.array([[[1.,2.],[1.,1.]], [[1.,1.],[3.,3.]]])
    >>> ac = a.copy()
    >>> # each block of a is inverted in-place
    >>> pinv_array(a)

    """
    n = a.shape[0]
    m = a.shape[1]

    if m == 1:
        # Pseudo-inverse of 1 x 1 matrices is trivial
        zero_entries = (a == 0.0).nonzero()[0]
        a[zero_entries] = 1.0
        a[:] = 1.0/a
        a[zero_entries] = 0.0
        del zero_entries

    else:
        # The block size is greater than 1

        # Create necessary arrays and function pointers for calculating pinv
        gelss, gelss_lwork = get_lapack_funcs(('gelss', 'gelss_lwork'),
                                              (np.ones((1,), dtype=a.dtype)))
        RHS = np.eye(m, dtype=a.dtype)
        lwork = _compute_lwork(gelss_lwork, m, m, m)

        # Choose tolerance for which singular values are zero in *gelss below
        if cond is None:
            t = a.dtype.char
            eps = np.finfo(np.float).eps
            feps = np.finfo(np.single).eps
            geps = np.finfo(np.longfloat).eps
            _array_precision = {'f': 0, 'd': 1, 'g': 2, 'F': 0, 'D': 1, 'G': 2}
            cond = {0: feps*1e3, 1: eps*1e6, 2: geps*1e6}[_array_precision[t]]

        # Invert each block of a
        for kk in range(n):
            gelssoutput = gelss(a[kk], RHS, cond=cond, lwork=lwork,
                                overwrite_a=True, overwrite_b=False)
            a[kk] = gelssoutput[1]
def pinv_array ( a , cond = None ) : n = a . shape [ 0 ] m = a . shape [ 1 ] if m == 1 : zero_entries = ( a == 0.0 ) . nonzero ( ) [ 0 ] a [ zero_entries ] = 1.0 a [ : ] = 1.0 / a a [ zero_entries ] = 0.0 del zero_entries else : gelss , gelss_lwork = get_lapack_funcs ( ( 'gelss' , 'gelss_lwork' ) , ( np . ones ( ( 1 , ) , dtype = a . dtype ) ) ) RHS = np . eye ( m , dtype = a . dtype ) lwork = _compute_lwork ( gelss_lwork , m , m , m ) if cond is None : t = a . dtype . char eps = np . finfo ( np . float ) . eps feps = np . finfo ( np . single ) . eps geps = np . finfo ( np . longfloat ) . eps _array_precision = { 'f' : 0 , 'd' : 1 , 'g' : 2 , 'F' : 0 , 'D' : 1 , 'G' : 2 } cond = { 0 : feps * 1e3 , 1 : eps * 1e6 , 2 : geps * 1e6 } [ _array_precision [ t ] ] for kk in range ( n ) : gelssoutput = gelss ( a [ kk ] , RHS , cond = cond , lwork = lwork , overwrite_a = True , overwrite_b = False ) a [ kk ] = gelssoutput [ 1 ]
Calculate the Moore-Penrose pseudo inverse of each block of the three dimensional array a.

    Parameters
    ----------
    a   : {dense array}
        Is of size (n, m, m)
    cond : {float}
        Used by gelss to filter numerically zeros singular values.
        If None, a suitable value is chosen for you.

    Returns
    -------
    Nothing, a is modified in place so that a[k] holds the pseudoinverse
    of that block.

    Notes
    -----
    By using lapack wrappers, this can be much faster for large n, than
    directly calling pinv2

    Examples
    --------
    >>> import numpy as np
    >>> from pyamg.util.linalg import pinv_array
    >>> a = np.array([[[1.,2.],[1.,1.]], [[1.,1.],[3.,3.]]])
    >>> ac = a.copy()
    >>> # each block of a is inverted in-place
    >>> pinv_array(a)