@verbose
def _tf_mixed_norm_solver_bcd_(M, G, Z, active_set, alpha_space, alpha_time,
                               lipschitz_constant, phi, phiT,
                               wsize=64, tstep=4, n_orient=1,
                               maxit=200, tol=1e-8, log_objective=True,
                               perc=None, verbose=None):
    # First make G fortran for faster access to blocks of columns
    G = np.asfortranarray(G)

    n_sensors, n_times = M.shape
    n_sources = G.shape[1]
    n_positions = n_sources // n_orient

    n_step = int(ceil(n_times / float(tstep)))
    n_freq = wsize // 2 + 1
    shape = (-1, n_freq, n_step)

    G = dict(zip(np.arange(n_positions), np.hsplit(G, n_positions)))
    R = M.copy()  # residual
    active = np.where(active_set)[0][::n_orient] // n_orient
    for idx in active:
        R -= np.dot(G[idx], phiT(Z[idx]))

    E = []  # track cost function

    alpha_time_lc = alpha_time / lipschitz_constant
    alpha_space_lc = alpha_space / lipschitz_constant

    converged = False

    for i in range(maxit):
        val_norm_l21_tf = 0.0
        val_norm_l1_tf = 0.0
        max_diff = 0.0
        active_set_0 = active_set.copy()
        for j in range(n_positions):
            ids = j * n_orient
            ide = ids + n_orient

            G_j = G[j]
            Z_j = Z[j]
            active_set_j = active_set[ids:ide]

            Z0 = deepcopy(Z_j)

            was_active = np.any(active_set_j)

            # gradient step
            GTR = np.dot(G_j.T, R) / lipschitz_constant[j]
            X_j_new = GTR.copy()

            if was_active:
                X_j = phiT(Z_j)
                R += np.dot(G_j, X_j)
                X_j_new += X_j

            rows_norm = linalg.norm(X_j_new, 'fro')
            if rows_norm <= alpha_space_lc[j]:
                if was_active:
                    Z[j] = 0.0
                    active_set_j[:] = False
            else:
                if was_active:
                    Z_j_new = Z_j + phi(GTR)
                else:
                    Z_j_new = phi(GTR)

                col_norm = np.sqrt(np.sum(np.abs(Z_j_new) ** 2, axis=0))

                if np.all(col_norm <= alpha_time_lc[j]):
                    Z[j] = 0.0
                    active_set_j[:] = False
                else:
                    # l1
                    shrink = np.maximum(1.0 - alpha_time_lc[j] / np.maximum(
                                        col_norm, alpha_time_lc[j]), 0.0)
                    Z_j_new *= shrink[np.newaxis, :]

                    # l21
                    shape_init = Z_j_new.shape
                    Z_j_new = Z_j_new.reshape(*shape)
                    row_norm = np.sqrt(stft_norm2(Z_j_new).sum())
                    if row_norm <= alpha_space_lc[j]:
                        Z[j] = 0.0
                        active_set_j[:] = False
                    else:
                        shrink = np.maximum(1.0 - alpha_space_lc[j] /
                                            np.maximum(row_norm,
                                            alpha_space_lc[j]), 0.0)
                        Z_j_new *= shrink
                        Z[j] = Z_j_new.reshape(-1, *shape_init[1:]).copy()
                        active_set_j[:] = True
                        R -= np.dot(G_j, phiT(Z[j]))

                        if log_objective:
                            val_norm_l21_tf += norm_l21_tf(
                                Z[j], shape, n_orient)
                            val_norm_l1_tf += norm_l1_tf(
                                Z[j], shape, n_orient)

            max_diff = np.maximum(max_diff, np.max(np.abs(Z[j] - Z0)))

        if log_objective:  # log cost function value
            pobj = (0.5 * (R ** 2.).sum() + alpha_space * val_norm_l21_tf +
                    alpha_time * val_norm_l1_tf)
            E.append(pobj)
            logger.info("Iteration %d :: pobj %f :: n_active %d" % (i + 1,
                        pobj, np.sum(active_set) / n_orient))
        else:
            logger.info("Iteration %d" % (i + 1))

        if perc is not None:
            if perc * n_positions <= np.sum(active_set) / float(n_orient):
                break

        if np.array_equal(active_set, active_set_0):
            if max_diff < tol:
                logger.info("Convergence reached !")
                converged = True
                break

    return Z, active_set, E, converged
dataset/ETHPy150Open mne-tools/mne-python/mne/inverse_sparse/mxne_optim.py/_tf_mixed_norm_solver_bcd_/swapped operands of `<=`