Scipy - найти основания пространства столбцов матрицы

Я пытаюсь написать простой алгоритм Simplex, первым шагом которого является поиск базового возможного решения:

1. Выберите набор B линейно независимых столбцов A
2. Установите все компоненты x, соответствующие столбцам, не входящим в B, равными нулю.
3. Решите получившиеся m уравнений, чтобы определить компоненты x. Это основные переменные.

Я знаю, что решение будет включать использование scipy.linalg.svd (или scipy.linalg.lu) и некоторой магии numpy.argwhere / numpy.where, но я не уверен, как именно.

Есть ли у кого-нибудь чистая реализация Numpy/Scipy для поиска основы (шаг 1) или, что еще лучше, всего вышеперечисленного?

Пример:

>>> A
array([[1, 1, 1, 1, 0, 0, 0],
       [1, 0, 0, 0, 1, 0, 0],
       [0, 0, 1, 0, 0, 1, 0],
       [0, 3, 1, 0, 0, 0, 1]])

>>> u, s, v = scipy.linalg.svd(A)
>>> non_zero, = numpy.where(s > 1e-7)
>>> rank = len(non_zero)
>>> rank
4
>>> for basis in some_unknown_function(A):
...     print(basis)
{3, 4, 5, 6}
{1, 4, 5, 6}

и так далее.