Умножение большого вектора комплексного числа на Scalar эффективно C ++ - PullRequest
Новая
       15

Умножение большого вектора комплексного числа на Scalar эффективно C ++

4 голосов
/ 28 июля 2011

В настоящее время я пытаюсь наиболее эффективно выполнить на месте умножение массива комплексных чисел (память выровнена так же, как в std :: complex, но в настоящее время мы используем нашу собственную ADT) на массив скалярных значений, которыеимеет тот же размер, что и массив комплексных чисел.

Алгоритм уже распараллелен, то есть вызывающий объект разбивает работу на потоки.Этот расчет выполняется для массивов в сотнях миллионов, поэтому для его завершения может потребоваться некоторое время.CUDA не является решением для этого продукта, хотя я бы хотел, чтобы это было.У меня есть доступ к надстройке, и поэтому у меня есть некоторый потенциал для использования BLAS / uBLAS.

Я думаю, однако, что SIMD может дать гораздо лучшие результаты, но я недостаточно знаком с тем, как это сделать.с комплексными числами.У меня есть следующий код (помните, что он разбит на потоки, соответствующие количеству ядер на целевой машине).Целевая машина также неизвестна.Таким образом, общий подход, вероятно, лучше всего. 

void cmult_scalar_inplace(fcomplex *values, const int start, const int end, const float *scalar)
{
    for (register int idx = start; idx < end; ++idx)
    {
        values[idx].real *= scalar[idx];
        values[idx].imag *= scalar[idx];
    }
}

fcomplex определяется следующим образом: 

struct fcomplex
{
    float real;
    float imag;
};

Я пробовал вручную развернуть цикл, так как мой счетчик цикла finally всегда будетстепень 2, но компилятор уже делает это для меня (я развернул до 32).Я попробовал ссылаться на скаляр const float - думая, что я сохраню один доступ - и это оказалось равным тому, что компилятор уже делал.Я пробовал STL и Transform, результаты которых близки к игре, но все еще хуже.Я также попытался привести к std :: complex и разрешить ему использовать перегруженный оператор для скалярного * complex для умножения, но в конечном итоге это дало те же результаты.

Итак, у кого-нибудь есть идеи?Огромная благодарность за то, что вы уделили время этому размышлению!Целевой платформой является Windows.Я использую Visual Studio 2008. Продукт также не может содержать код GPL!Большое спасибо.

Ответы [ 4 ]

1 голос
/ 28 июля 2011

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

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

То, что я хотел бы рассмотреть, - это использовать разные размеры развертки, а также играть с выравниванием scalar и values (шаблон доступа к памяти может оказывать большое влияние на эффекты кэширования).

Для решения проблемы нежелательной сериализации можно посмотреть, что является сгенерированным кодом для чего-то вроде 

float r0 = values[i].real, i0 = values[i].imag, s0 = scalar[i];
float r1 = values[i+1].real, i1 = values[i+1].imag, s1 = scalar[i+1];
float r2 = values[i+2].real, i2 = values[i+2].imag, s2 = scalar[i+2];
values[i].real = r0*s0; values[i].imag = i0*s0;
values[i+1].real = r1*s1; values[i+1].imag = i1*s1;
values[i+2].real = r2*s2; values[i+2].imag = i2*s2;

потому что здесь у оптимизатора в теории немного больше свободы.

1 голос
/ 28 июля 2011

Лучше всего будет использовать оптимизированный BLAS, который будет использовать все, что доступно на вашей целевой платформе.

1 голос
/ 28 июля 2011

Вы можете сделать это довольно легко с SSE, например, 

void cmult_scalar_inplace(fcomplex *values, const int start, const int end, const float *scalar)
{
    for (int idx = start; idx < end; idx += 2)
    {
        __m128 vc = _mm_load_ps((float *)&values[idx]);
        __m128 vk = _mm_set_ps(scalar[idx + 1], scalar[idx + 1], scalar[idx], scalar[idx]);
        vc = _mm_mul_ps(vc, vk);
        _mm_store_ps((float *)&values[idx], vc);
    }
}

Обратите внимание, что values и scalar должны быть выровнены на 16 байт.

Или вы можете просто использоватьКомпилятор Intel ICC и пусть он сделает за вас тяжелую работу.

ОБНОВЛЕНИЕ

Вот улучшенная версия, которая развертывает цикл с коэффициентом 2 ииспользует одну инструкцию загрузки, чтобы получить 4 скалярных значения, которые затем распаковываются в два вектора: 

void cmult_scalar_inplace(fcomplex *values, const int start, const int end, const float *scalar)
{
    for (int idx = start; idx < end; idx += 4)
    {
        __m128 vc0 = _mm_load_ps((float *)&values[idx]);
        __m128 vc1 = _mm_load_ps((float *)&values[idx + 2]);
        __m128 vk = _mm_load_ps(&scalar[idx]);
        __m128 vk0 = _mm_shuffle_ps(vk, vk, 0x50);
        __m128 vk1 = _mm_shuffle_ps(vk, vk, 0xfa);
        vc0 = _mm_mul_ps(vc0, vk0);
        vc1 = _mm_mul_ps(vc1, vk1);
        _mm_store_ps((float *)&values[idx], vc0);
        _mm_store_ps((float *)&values[idx + 2], vc1);
    }
}
0 голосов
/ 28 июля 2011

Есть ли у вас доступ к интегрированным примитивам производительности Intel? Интегрированные примитивы производительности У них есть ряд функций, которые обрабатывают подобные случаи с довольно приличной производительностью. Возможно, вам удастся решить вашу конкретную проблему, но я не удивлюсь, если ваш компилятор уже выполнит достойную работу по оптимизации кода.

...