Подтвердить что ты не робот

Алгоритм для индексных чисел треугольных матричных коэффициентов

Я думаю, что это должно быть просто, но я не могу понять это правильно.

У меня есть треугольная матрица MxM, коэффициенты которой хранятся в векторе, строка за строкой. Например:

M =   [ m00 m01 m02 m03 ] 
      [     m11 m12 m12 ]
      [         m22 m23 ]
      [             m33 ]

хранится как

coef[ m00 m01 m02 m03 m11 m12 m13 m22 m23 m33 ]

Теперь я ищу нерекурсивный алгоритм, который дает мне размер матрицы M и индекс массива коэффициентов i

unsigned int row_index(i,M)

и

unsigned int column_index(i,M)

элемента матрицы, на который он ссылается. Так, row_index(9,4) == 3, column_index(7,4) == 2 и т.д., если подсчет индекса основан на нуле.

EDIT: несколько ответов с использованием итерации. Кто-нибудь знает об алгебраических выражениях?

4b9b3361

Ответ 1

Здесь алгебраическое (в основном) решение:

unsigned int row_index( unsigned int i, unsigned int M ){
    double m = M;
    double row = (-2*m - 1 + sqrt( (4*m*(m+1) - 8*(double)i - 7) )) / -2;
    if( row == (double)(int) row ) row -= 1;
    return (unsigned int) row;
}


unsigned int column_index( unsigned int i, unsigned int M ){
    unsigned int row = row_index( i, M);
    return  i - M * row + row*(row+1) / 2;
}

EDIT: fixed row_index()

Ответ 2

Одиночные вкладыши в конце этого ответа, объяснение следует: -)

В массиве коэффициентов есть: M элементов для первой строки (строка 0, в вашей индексации), (M-1) для второго (строка 1) и т.д. для всего M + (M-1) + & hellip; + 1 = M (M + 1)/2 элементов.

Немного легче работать с конца, потому что массив коэффициентов всегда имеет 1 элемент для последней строки (строка M-1), 2 элемента для второй строки (строка M-2), 3 элемента для строки М-3 и т.д. Последние K-строки занимают последние K (K + 1)/2 позиции массива коэффициентов.

Теперь предположим, что вам предоставляется индекс я в массиве коэффициентов. Есть элементы M (M + 1)/2-1-i в положениях, больших i. Назовите это число я '; вы хотите найти наибольшее k такое, что k (k + 1)/2 & le; я '. Форма этой проблемы является источником ваших страданий - насколько я вижу, вы не можете избежать квадратных корней: -)

Хорошо, пусть это все равно: k (k + 1) & le; 2i 'означает (k + 1/2) 2 - 1/4 & le; 2i 'или эквивалентно k & le; (SQRT (8i '+ 1) -1)/2. Позвольте мне назвать наибольший такой k, как K, то есть K строк, которые позже (из общего числа M строк), поэтому row_index (i, M) является M-1-K. Что касается индекса столбца, то K (K + 1)/2 элементов из я 'находятся в более поздних строках, поэтому в этой строке есть j' = i'-K (K + 1)/2, более поздние элементы K + 1 элементов), поэтому индекс столбца K-j '. [Или, что то же самое, эта строка начинается с (K + 1) (K + 2)/2 с конца, поэтому нам просто нужно вычесть начальную позицию этой строки из i: i- [M (M + 1)/2 - (К + 1) (К + 2)/2]. Отрадно, что оба выражения дают один и тот же ответ!]

(Pseudo-) Code [используя ii вместо я ', поскольку некоторые языки не допускают переменные с именами, такими как i'; -)]:

row_index(i, M):
    ii = M(M+1)/2-1-i
    K = floor((sqrt(8ii+1)-1)/2)
    return M-1-K

column_index(i, M):
    ii = M(M+1)/2-1-i
    K = floor((sqrt(8ii+1)-1)/2)
    return i - M(M+1)/2 + (K+1)(K+2)/2

Конечно, вы можете превратить их в однострочные, заменив выражения для ii и K. Возможно, вам придется быть осторожными с ошибками точности, но есть способы найти целочисленный квадратный корень, который не требует вычисления с плавающей запятой, Кроме того, процитировать Кнута: "Остерегайтесь ошибок в приведенном выше коде, я только доказал это правильно, не пробовал".

Если я могу рискнуть сделать еще одно замечание: просто сохраняя все значения в M & times; массив M будет занимать не более двух раз больше памяти, и в зависимости от вашей проблемы коэффициент 2 может быть незначителен по сравнению с алгоритмическими улучшениями, и, возможно, стоит торговать, возможно, дорогостоящим вычислением квадратного корня для более простых выражений, которые у вас будут.

[Edit: BTW, вы можете доказать, что floor ((sqrt (8ii + 1) -1)/2) = (isqrt (8ii + 1) -1)/2 где isqrt (x) = floor (sqrt ( x)) представляет собой целочисленный квадратный корень, а деление - целочисленное деление (усечение, значение по умолчанию в C/С++/Java и т.д.) - поэтому, если вы беспокоитесь о проблемах с точностью, вам нужно только беспокоиться о реализации правильного целочисленного квадрата корень.]

Ответ 3

Должно быть, что

i == col + row*(M-1)-row*(row-1)/2

Итак, один способ найти col и row - перебирать возможные значения строки:

for(row = 0; row < M; row++){
  col = i - row*(M-1)-row*(row-1)/2
  if (row <= col < M) return (row,column);
}

Это, по крайней мере, нерекурсивный, я не знаю, можете ли вы сделать это без итерации.

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

Ответ 4

Объяснение ShreevatsaR превосходно и помогло мне решить мою проблему. Однако объяснение и код, предоставленные для индекса столбца, дают немного другой ответ, чем тот, который задает вопрос.

Чтобы повторить, в строке после я есть j '= i' - K (K + 1)/2 элементов. Но строка, как и любая другая строка, имеет M элементов. Следовательно, индекс столбца (на основе нуля) равен y = M - 1 - j '.

Соответствующий псевдокод:

column_index(i, M):
  ii = M(M+1)/2-1-i
  K = floor((sqrt(8ii+1)-1)/2)
  jj = ii - K(K+1)/2
  return M-1-jj

Ответ ShreevatsaR, K - j ', является индексом столбца при начале отсчета (с нулем) по диагонали матрицы. Следовательно, его вычисление дает column_index (7,4) = 0 и не, как указано в вопросе, column_index (7,4) = 2.

Ответ 5

Для них может быть умный один вкладыш, но (минус любая проверка ошибок):

unsigned int row_index( unsigned int i, unsigned int M ){
    unsigned int row = 0;
    unsigned int delta = M - 1;
    for( unsigned int x = delta; x < i; x += delta-- ){
        row++;
    }
    return row;
}

unsigned int column_index( unsigned int i, unsigned int M ){
    unsigned int row = 0;
    unsigned int delta = M - 1;
    unsigned int x;
    for( x = delta; x < i; x += delta-- ){
        row++;
    }
    return M + i - x - 1;
}

Ответ 6

Я подумал немного, и я получил следующий результат. Обратите внимание, что вы получаете как строки, так и столбцы на один снимок.

Предположения: строки начинаются с 0. Начало столбцов начинается с нуля. Начало индекса при 0

нотация

N = размер матрицы (был M в исходной задаче)

m = Индекс элемента

Код psuedo

function ind2subTriu(m,N)
{
  d = 0;
  i = -1;
  while d < m
  {
    i = i + 1
    d = i*(N-1) - i*(i-1)/2
  }
  i0 = i-1;
  j0 = m - i0*(N-1) + i0*(i0-1)/2 + i0 + 1;
  return i0,j0
}

И некоторый код октавы /matlab

function [i0 j0]= ind2subTriu(m,N)
 I = 0:N-2;
 d = I*(N-1)-I.*(I-1)/2;
 i0 = I(find (d < m,1,'last'));
 j0 = m - d(i0+1) + i0 + 1;

Как вы думаете?

По состоянию на декабрь 2011 года действительно хороший код для этого был добавлен в GNU/Octave. Потенциально они будут расширять sub2ind и ind2sub. Код может быть найден на данный момент как частные функции ind2sub_tril и sub2ind_tril

Ответ 7

Пришло время понять, что вам нужно!:)

unsigned int row_index(int i, int m)
{
    int iCurrentRow = 0;
    int iTotalItems = 0;
    for(int j = m; j > 0; j--)
    {
        iTotalItems += j;

        if( (i+1) <= iTotalItems)
            return iCurrentRow;

        iCurrentRow ++;
    }

    return -1; // Not checking if "i" can be in a MxM matrix.
}

Извините, забыл другую функцию.....

unsigned int column_index(int i, int m)
{
    int iTotalItems = 0;
    for(int j = m; j > 0; j--)
    {
        iTotalItems += j;

        if( (i+1) <= iTotalItems)
            return m - (iTotalItems - i);
    }

    return -1; // Not checking if "i" can be in a MxM matrix.
}