Оптимизация циклов C для получения диагонали массива - PullRequest
7 голосов
/ 26 февраля 2011

Великий Бог, Google не предоставил мне объяснения некоторых проблем оптимизации цикла. Итак, из-за печали, что у меня недостаточно Google-фу, я обращаюсь к вам StackOverflow.

Я оптимизирую программу на C для решения конкретной системы дифференциальных уравнений. В процессе поиска численного решения я вызываю функцию, которая устанавливает линейную систему уравнений, а затем функцию, которая ее решает.

Функция решения изначально имела узкое место при доступе к элементам по диагонали массива, который определяет линейную систему. Поэтому я включил одномерный массив, который устанавливается во время инициализации системы и содержит значения по диагонали массива.

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

Примечание: Я собрал все версии, которые я пробовал, в одну функцию и профилировал эту функцию, чтобы увидеть, где было потрачено время. Я сообщу время выполнения для версии в процентах от общего времени в функции. Функция была оценена несколько миллионов раз. Чем меньше номер, тем лучше.

Соответствующие объявления данных, используемых в коде:

/* quick definitions of the relevant variables, data is a struct */

static const int sp_diag_ind[98] = {2,12,23,76,120,129,137,142,.../* long list */};

double *spJ = &(data->spJ[0]);
/* data has double spJ[908] that represents a sparse matrix stored in triplet
*  form, I grab the pointer because I've found it to be more 
*  efficient than referencing data->spJ[x] each time I need it
*/

int iter,jter;
double *diag_data = NV_DATA_S(data->J_diag);
/* data->J_diag has a content field that has an array double diag_data[150]
*  NV_DATA_S is a macro to return the pointer to the relevant data
*/

Мой оригинальный цикл для инициализации diag_data. Время составило 16,1% от оценки (см. Примечание).

/* try 1 */
for (iter = 0; iter<3; iter++) {
    diag_data[iter] = 0; 
}
jter = 0;
for (iter = 3; iter<101; iter++) { // unaligned loop start
    diag_data[iter] = spJ[sp_diag_ind[jter]];
    jter++; // heavy line for loop
}

for (iter = 101; iter<150; iter++) {
    diag_data[iter] = 0; 
}

Подводя итог, мы берем указатель на диагональ, устанавливаем некоторые компоненты на ноль (это не является обязательным в зависимости от используемого мной алгоритма), затем получаем значения, находящиеся на диагонали представленного «массива». в редкой форме от spJ. Поскольку spJ является одномерным массивом из 908 ненулей (в основном нулевого) массива 150x150, мы должны использовать поиск, чтобы найти позиции диагональных элементов в spJ. Этот поиск определяется массивом 98 элементов sp_diag_ind.

Я попытался отказаться от использования jter, потому что он оказался несвободным для увеличения. Средняя петля моя вторая попытка:

for (iter = 0; iter<98; iter++) { // unaligned loop start
    diag_data[iter+3] = spJ[sp_diag_ind[iter]];
}

Это немного улучшило ситуацию. Время было 15,6% для этой версии. Но когда я смотрю на анализ Shark этого кода (инструмент, который поставляется с XCode на Mac), он предупреждает меня, что это цикл без выравнивания.

Третья попытка для улучшения заключалась в удалении циклов "обнуления" и использовании memset для обнуления diag_data:

memset(diag_data, '\0', sizeof(diag_data));

for (iter = 0; iter<98; iter++) { // unaligned loop start
    diag_data[iter+3] = spJ[sp_diag_ind[iter]];
}

Это было приурочено к 14,9%. Не зная, что такое не выровненный цикл, я продолжал возиться. Я нашел улучшенную четвертую реализацию , в которой выполняется смещение выравнивания между diag_data и spJ [сумасшедший индекс] с указателем:

realtype * diag_mask = &diag_data[3];
for (iter = 0; iter<98; iter++) { // unaligned loop start
    diag_mask[iter] = spJ[sp_diag_ind[iter]];
}

Использование diag_mask позволило немного улучшить скорость. Пришло 13,1%.

Редактировать: Оказалось, что этот раздел был глупее, чем я первоначально думал. Использование iter не определено. Подсказки к @caf и @rlibby для его ловли .

Наконец, я попробовал то, что мне показалось глупым:

memset(diag_data, '\0', sizeof(diag_data));

for (iter = 0; iter<98;) {
    diag_mask[iter] = spJ[sp_diag_ind[iter++]];
}

Это было рассчитано на 10,9%. Кроме того, Shark не выдает предупреждение о невыровненном цикле, когда я смотрю на аннотированный исходный код. Конец глупой секции

Итак, мои вопросы:

  1. Что такое невыровненный цикл?
  2. Почему пятая реализация выровнена, а четвертая нет?
  3. Имеет ли выровненный цикл ответственность за улучшение скорости выполнения между моей четвертой и пятой реализациями или встраивает шаг приращения в поиск значения sp_diag_ind ответственного?
  4. Видите ли вы какие-либо улучшения, которые я могу сделать?

Спасибо за помощь.

- Andrew

Ответы [ 3 ]

2 голосов
/ 26 февраля 2011

Цикл без выравнивания - это цикл, в котором первая инструкция не начинается на определенной границе (кратной 16 или 32).Должен быть флаг компилятора для выравнивания циклов;это может или не может помочь производительности.Независимо от того, выровнен ли цикл или нет без флага, все зависит от того, какие инструкции стоят перед ним, поэтому он не предсказуем.Еще одна оптимизация, которую вы можете попробовать - пометить diag_mask, spJ и sp_diag_ind как restrict (функция C99).Это означает, что они не являются псевдонимами и могут помочь компилятору лучше оптимизировать цикл.Счетчик 98, вероятно, будет слишком мал, чтобы увидеть какой-либо эффект.

1 голос
/ 26 февраля 2011

Видите ли вы какие-либо другие улучшения, которые я могу сделать?

Вы настраиваете дневное освещение на что-то, что использует около 11% времени.Разве в остальных 89% нет ничего, что можно было бы оптимизировать?

1 голос
/ 26 февраля 2011

Ваша пятая версия неверна - она ​​имеет неопределенное поведение , потому что она и изменяет iter, и ссылается на ее значение, для цели, отличной от вычисления нового значения, без промежуточной точки последовательности. *

Вы пытались сохранить фактические значения диагоналей, а не их индексы в пределах spJ, в точке, которую вы вычисляете sp_diag_ind[]? Тогда вы можете просто скопировать их прямо в diag_data (или, что еще лучше, напрямую использовать вектор диагоналей).


Соответствующая часть стандарта C - §6.5 Выражения:

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

Это относится к объекту iter в вашем выражении. Нарушение «обязательного» ограничения - неопределенное поведение.

gcc (протестировано с версией 4.4.5) даже предупреждает о вашем выражении:

x.c:16: warning: operation on ‘iter’ may be undefined
...