def _svdgram(A, tol=None, tol2=1e-7):
    ''' Highly experimental acceleration of SVD/QR using Gram matrix.
        Use with caution for m>>n only!
        function [u,s,r]=_svdgram(A,[tol])
        u is the left singular factor of A,
        s is the singular values (vector!),
        r has the meaning of diag(s)*v'.
        if tol is given, performs the truncation with Fro-threshold.
    '''

    R2 = _np.dot(_tconj(A), A)
    [u, s, vt] = _np.linalg.svd(R2, full_matrices=False)
    u = _np.dot(A, _tconj(vt))
    s = (u**2).sum(axis=0)
    s = s ** 0.5
    if tol is not None:
        p = _my_chop2(s, _np.linalg.norm(s) * tol)
        u = u[:, :p]
        s = s[:p]
        vt = vt[:p, :]

    tmp = _spdiags(1. / s, 0, len(s), len(s)).todense()
    tmp = _np.array(tmp)
    u = _np.dot(u, tmp)
    r = _np.dot(_np.diag(s), vt)

    # Run chol for reortogonalization.
    # It will stop if the matrix will be singular.
    # Fortunately, it means rank truncation with eps in our business.

    if (s[0] / s[-1] > 1. / tol2):
        p = 1
        while (p > 0):
            R2 = _np.dot(_tconj(u), u)
            #[u_r2, s_r2, vt_r2] = _np.linalg.svd(R2) # in matlab [R, p] = chol(a) - here is a *dirty* patch
            #p = s_r2[s_r2 > eps].size
            #R2 = R2[:p, :p]
            R = _cholesky(R2, lower=True)
            if (p > 0):
                u = u[:, :p]
                s = s[:p]
                r = r[:p, :]
            iR = _np.linalg.inv(R)
            u = _np.dot(u, iR)
            r = _np.dot(R, r)
    return u, s, r
Highly experimental acceleration of SVD/QR using Gram matrix.
        Use with caution for m>>n only!
        function [u,s,r]=_svdgram(A,[tol])
        u is the left singular factor of A,
        s is the singular values (vector!),
        r has the meaning of diag(s)*v'.
        if tol is given, performs the truncation with Fro-threshold.