Доступ к подматрицам с использованием LAPACK

Есть ли функция в LAPACK, которая даст мне элементы конкретной подматрицы? Если да, то каков синтаксис в С++?

Или мне нужно его закодировать?

Ответы

Ответ 1

Нет функции для доступа к подматрице. Однако из-за того, что данные матрицы хранятся в процедурах LAPACK, вам не нужен. Это экономит много копий, и расположение данных было (частично) выбрано по этой причине:

Напомним, что плотная (т.е. неполосная, треугольная, эрмитова и т.д.) матрица в LAPACK определяется четырьмя значениями:

  • указатель на верхний левый угол матрицы
  • количество строк в матрице
  • количество столбцов в матрице
  • "ведущий размер" матрицы; обычно это расстояние в памяти между соседними элементами строки.

В большинстве случаев большинство людей используют только ведущее измерение, равное количеству строк; матрица 3x3 обычно хранится так:

a[0] a[3] a[6] 
a[1] a[4] a[7]
a[2] a[5] a[8]

Предположим, что нам нужна подматрица 3x3 огромной матрицы с ведущим размером lda. Предположим, мы специально хотим, чтобы подматрица 3x3, верхний левый угол которой расположен в точке (15,42):

         .             .            .
         .             .            .
... a[15+42*lda] a[15+43*lda] a[15+44*lda] ...
... a[16+42*lda] a[16+43*lda] a[16+44*lda] ...
... a[17+42*lda] a[17+43*lda] a[17+44*lda] ...
         .             .            .
         .             .            .

Мы могли бы скопировать эту матрицу 3x3 в непрерывное хранилище, но если мы хотим передать ее как входную (или выходную) матрицу в процедуру LAPACK, нам не нужно; нам нужно только определить параметры соответствующим образом. Позвольте называть эту подматрицу b; затем определим:

// pointer to the top-left corner of b:
float *b = &a[15 + 42*lda];
// number of rows in b:
const int nb = 3;
// number of columns in b:
const int mb = 3;
// leading dimension of b:
const int ldb = lda;

Единственное, что может быть удивительно, это значение ldb; используя значение lda "большой матрицы", мы можем обращаться к подматрице без копирования и работать на ней на месте.

Однако Я солгал (вроде). Иногда вы действительно не можете работать на подматрице на месте и действительно должны ее скопировать. Я не хотел об этом говорить, потому что это редко, и вы должны использовать операции на месте, когда это возможно, но я бы плохо себя чувствовал, не говоря вам, что это возможно. Подпрограмма:

SLACPY(UPLO,M,N,A,LDA,B,LDB)

копирует матрицу M x N, верхний левый угол которой A и сохраняется с ведущим размером lda в матрице M x N, верхний левый угол которой равен b и имеет ведущий размер ldb. Параметр UPLO указывает, следует ли копировать верхний треугольник, нижний треугольник или всю матрицу.

В примере, который я привел выше, вы бы использовали его так (при условии привязки клапека):

...
const int m = 3;
const int n = 3;
float b[9];
const int ldb = 3;
slacpy("A",  // anything except "U" or "L" means "copy everything"
       &m,   // number of rows to copy
       &n,   // number of columns to copy
       &a[15 + 42*lda], // pointer to top-left element to copy
       lda,  // leading dimension of a (something huge)
       b,    // pointer to top-left element of destination
       ldb); // leading dimension of b (== m, so storage is dense)
...