Скорость работы с __constant__ памятью

Доброго времени суток.
Решил в порядке изучения предмета выполнить свертку двух функций (u и g), заданных одномерными массивами, "в лоб". Написал три ядра, различающихся только тем, в какой памяти лежат исходные данные (результат во всех случаях заносится в глобальную память).

1) оба массива в глобальной
__global__ void DirConv1(short* u, short* g, int len_g, float* conv)
//u - массив отсчетов сигнала, len_u - число отсчетов u, g - массив опорной ф-ции,len_g - число отсчетов g
//conv - массив результата свертки . Все массивы выделяются вызывающей функцией
{
int idx = threadIdx.x + blockIdx.x * blockDim.x; //глобальный индекс нити

float sum=0;
int j;
for (j=0; j < len_g; j++)
{
sum = sum + (long int)g[j] * u[idx+j];
}
conv[idx] = sum /len_g;
}
/////////////////////////////////////////////////////////////////////////////////////////////////////////////
2) Один массив (меньшего размера) - в глобальной, другой - в константной
__global__ void DirConv1C(short* u, int len_g, float* conv)
{
int idx = threadIdx.x + blockIdx.x * blockDim.x; //глобальный индекс нити

float sum=0;
int j;
for (j=0; j < len_g; j++)
{
sum = sum + (long int)const_g[j] * u[idx+j];
}
conv[idx] = sum /len_g;
}
/////////////////////////////////////////////////////////////////////////////////////////////////////////////
3) Оба - в константной
__global__ void DirConv1CC(int len_g, float* conv)
{
int idx = threadIdx.x + blockIdx.x * blockDim.x; //глобальный индекс нити

float sum=0;
int j;
for (j=0; j < len_g; j++)
{
sum = sum + (long int)const_g[j] * const_u[idx+j];
}
conv[idx] = sum /len_g;
}
////////////////////////////////////////////////////////////////////
///Объявления массивов:
__constant__ short const_g[2048];
__constant__ short const_u[16384 + 127*4];
//Выделение памяти для GPU
cutilSafeCall(cudaMalloc((void**)&gpuConv, N * sizeof(float) ) );
cutilSafeCall(cudaMalloc((void**)&gpu_u, len_total * sizeof(short) ) );
cutilSafeCall(cudaMalloc((void**)&gpu_g, len_g * sizeof(short) ) );

//копируем данные в GPU
cutilSafeCall(cudaMemcpy(gpu_u, u,len_total * sizeof(short), cudaMemcpyHostToDevice ) );
cutilSafeCall(cudaMemcpy(gpu_g, g,len_g * sizeof(short), cudaMemcpyHostToDevice ) );

cutilSafeCall(cudaMemcpyToSymbol(const_g, g, len_g * sizeof(short), 0, cudaMemcpyHostToDevice ) );
cutilSafeCall(cudaMemcpyToSymbol(const_u, u, len_total * sizeof(short), 0, cudaMemcpyHostToDevice ) );

Результаты времени выполнения ядер (GTX580, CUDAToolkit 3.2, XP x64 SP2):
1.Processing time: 0.186207 (ms)
2.Processing time: 0.149459 (ms)
3.Processing time: 0.724911 (ms)

Результаты вычислений, как ни странно ;), одинаковые
Вопрос: почему в последнем случае, когда ОБА массива размещены в константной памяти, происходит падение производительности? Не хватает объема кеша? Или в подобных случаях вообще лучше "разносить" данные по разной памяти?

Forums: 

Падение производительности

Скорее всего в последнем случае такая штука из-за того, что каждая нить обращается к константной памяти и читает из нее 2 числа, расположенных не подряд. В результате выборка из памяти в пределах варпа не может объединиться и разбивается на множество запросов доступа к ней.

ps а вообще nvidia каки.... маловато константной памяти. Хотел применять к цифровой обработке сигналов - на борту 4 гига оперы (quadro 5800FX) а константной уже очень давно 64k.... да и копирование туда-сюда долговато.... в итоге два xeon 6ти ядерных используя tbb и ipp уделывают квадру почти по всему

По формальным попугаям 5800FX

По формальным попугаям 5800FX если и быстрее 12-ти ядер на 3 гигагерца (частота условная, для примера), то не в большие разы. Ну может раза в полтора-два если без MAD. Это для float, для двойной точности CPU будет заметно быстрее.

Поэтому если у вас пересылки занимают заметное время в сравнении с вычислениями, нужно время на пересылки прятать: пока что-то считается, другая порция пересылается. Т.к. пересылка по DMA, процессор свободен, можно в Angry Birds играть.

Второй момент, на котором видеокарты блистают - это тригонометрия, экспоненты и логарифмы. В SSE таких инструкций не завезли, или на FP87 (безумно медленно) или на SSE "приближенные аналоги", все едино довольно долгие. Если точности "приближенных" sin/cos/exp/ln на видеокарте достаточно, то вот в этом месте выигрыш очень большой, порядок-полтора запросто.

Про попугаев

Насчет попугаев - да. 5800FX - 90 терафлопс а один xeon на 3.3 гигаерца - примерно 1 терафлоп. Но тут опятьже маркетинг - 90 у квадры - это говорит нвидиа, на определенной смеси команд. На всех тестах из sdk 5800 дает максимум 65 терафлоп (умножение с матрицами).

Это то понятно что пересылку можно делать и не часто - вход и выход чисто + дробление. Но когда много обработки то все валится. У меня работа с векторами размерностью под 200 миллионов флоатов - кернелы все быстро делают но потом то все равно результат брать придется из видяхи....

Посмотрим что будет в куде 4.0 =)

Какие-то странные у вас

Какие-то странные у вас попугаи, неизвестной мне системы.

На практике я вижу так
1) 5800FX - это практически GTX280, которая вот у меня есть. С точностью до количества и скорости памяти и клоков.
2) У GTX280 писали про терафлопс (один), но откуда получена эта цифра я не понимаю. На практике я вижу ~400 GFLOP/s на задаче вроде nbody и 370 GFLOP/s на умножении матриц сам мерял, см. на этом сайте.
Пример из SDK про умножение матриц - негодный.
3) С double - порядка 75 GFLOP/s на dgemm, тоже сам мерял.

Это все для GTX280, для квадры надо поправлять на скорости памяти/клоки.

Теперь про CPU:
4) Одно ядро фигачит 8 операций на такт на SSE в теории (float) и 4 в double. Это простой арифметики, эти операции реально быстрые.
5) 12 ядер x 3.3 гигагерца x 8 операций = 316.8 GFLOP/s в float и вдвое меньше для double
6) На практике процентов 80 от этого получаются достаточно легко (не всегда, но часто), процентов 90 - иногда тоже выходит.
На 6 ядрах и несложных задачах (вроде сложения пары векторов) реальный упор будет, скорее всего, в bandwidth памяти.

Что касается скорости пересылки: если у вас есть какие-то еще PCIe-устройства то нужно очень внимательно смотреть в какие слоты что воткнуто. Одно неловкое движение - и видеокарта работает на x8 вместо x16.
Не знаю про квадру, а обычные "бытовые" видеокарты в x8 работают на пересылках со скоростью гигабайта 2.5 в секунду (хотя, казалось бы, на PCIe 2.0 должно быть 4), а на x16 - у меня получается 4.5, иногда, очень редко (и от драйверов зависит) доходит до 5.

Что касается CUDA4, то железо то не меняется, откуда взяться особым чудесам?

Попугаи странные ибо это я из

Попугаи странные ибо это я из пояснительной записки взял где народ подбирал оборудование под наши задачи и приводил официальные цифры - а такие попугаи обычно на спец смеси команд.
Тест матричный и правда не очень....
А вот насчет скорости пересылки - и правда интересно. У меня больше 2,5 гигов в секунду не выходило... хотя показывает что 16x pcie, но ничего больше не стоит.... мб опять мамку глючную дали (заказывали на работе граф станции у РАМЕК - такое фуфло собирают за бешенные бабки, да еще уже 2 станции из 4 отвозили на гарантийный ремонт как раз из за pci-e)

У меня в боевой рабочей

У меня в боевой рабочей станции стоят ATI 5870 и GTX480.

На старой мамке (на X58 + nForce 200) было ~4Gb/sec у ATI (измерение OpenCL, через oclBandwidth имени NVidia) и около 5Gb/sec у NVidia. Это оба больше чем x8

С тех пор прошло время, сменилась материнка (стала на P67+NForce 200) и сменились драйвера. На новой материнке при неудачном подборе слотов получалось ~2.9Gb у ATI и 2.5 у NVidia, а если немного подумать, то 2.9 у ATI и 4.2 у NVidia.
Но материнка устроена так, что 2 по x16 можно получить только если никаких еще PCIe устройств нет, а у меня еще RAID. Т.е. 4+8+16 - это максимум, который на ней получается, увы.

>В результате выборка из

>В результате выборка из памяти в пределах варпа не может объединиться и разбивается на множество запросов доступа к ней.
Похоже, что-то такое и происходит. Но о-о-о-о-чччень ме-е-е-ее-дленно ;)

Как мне кажется, вы выяснили,

Как мне кажется, вы выяснили, что доступ к кэш-памяти быстрее, чем доступ к константам.

Кроме того, размер блока/грида мне непонятен, оттого сложно представить, что на самом деле происходит.

>размер блока/грида мне

>размер блока/грида мне непонятен

int N = 16384;
blockSize = 1024; numBlocks = N/blockSize;
DirConv1 <<< numBlocks, blockSize >>>(gpu_u, gpu_g, len_g, gpuConv);

Как-то не понимаю (не чувствую) я глубокий смысл в выборе размерности блока : одномерный, двумерный или трехмерный. По условиям задачи обрабатывется одномерный массив, так значит и размерность блока можно взять одномерную (если это обеспечивает нужное количество нитей - 16384) ?
Опять же, что лучше - увеличивать размер блока или количество блоков? Увеличение размера блока, как я понимаю, (помимо аппаратных пределов) ограничивается уменьшением доступной для каждого ядра разделяемой и регистровой памятью. Количество блоков имеет чисто аппаратное ограничение. В чем я неправ?

Про размер блока и их

Про размер блока и их количество вполне внятно написано в документации. Там, правда, речь идет о повышении occupancy т.е. о том, что надо запихать как можно больше потоков (warp-ов) на процессор.

Смысл в том, что если варпов много, то вся latency по доступу к памяти разных видов - эффективно прячется (пока варп ждет своей очереди - ему из памяти все прочитают и принесут)

Есть альтернативный подход Василия Волкова: самая быстрая память это регистры (быстрее shared), поэтому будем делать блоки с большим потреблением регистров, а тредов в блоке немного, 192 или 256. Occupancy будет маленькая, а общее быстродействие вырастет.

Если целиться в оптимальную исполняемость на любом чипе от G80 до ферми, то оптимальный размер блока получается 256 (или кратно меньше) - на G80 влезет три таких блока, на 1.3 - 4, на ферми - 6.

И для Fermi в-общем намекают,

И для Fermi в-общем намекают, что нужно порядка 22 варпов (700 тредов на SM) чтобы спрятать латентность инструкций и порядка 40 варпов чтобы спрятать латентность глобальной памяти (Programming Guide от 4.0RC2, раздел 5.2.3)

Т.е. блоки по 256 (и 6 штук на ферми) - хороший выбор в общем случае. Если регистров/shared хватит

Ну не укладываются в моей

Ну не укладываются в моей голове новые концепции параллельного программирования :(
Из чтения документации я вынес следующее.
Кол-во варпов зависит от кол-ва тридов, которое, в свою очередь задается алгоритмом задачи.
Но вот как все эти абстракции grid - blocks - threads связаны с реальными мультипроцессорами - хоть убей понять не могу.
В моей задаче imho как блоки ни комбинируй - все равно задействуются все 512 SM и не один раз

Ну да. Число варпов равно

Ну да. Число варпов равно числу потоков / 32 (округляя вверх).

А дальше принцип очень простой
- SM - это группа потоковых процессоров. В случае Fermi их 32 или 48 в одном SM. В вашей 580й - 16 SM по 32 'cuda core' (+SFU). В рамках SM есть регистры и разделяемая память, общие для этих самых core/SFU
- блок - это нечто, что исполняется на одном SM одновременно.
- в блоке есть треды, которые пилятся по 32 штуки (warp) и каждый warp исполняется SIMD (одна инструкция, разные данные). Номера тредов в варпе идут подряд (для двумерной нумерации тредов в блоке - не задумывался т.к. никогда не использовал)
- на одном SM одновременно могут исполняться несколько блоков, количество ограничено общим числом тредов, количеством регистров, количеством shared memory
- для блоков не гарантируется порядок исполнения, расположение на SM и так далее.

Как-то так. Если непонятно - спрашивайте, это - базовая вещь для классического подхода к CUDA ("обеспечиваем максимальное occupancy")

Есть альтернативный подход

Есть альтернативный подход Василия Волкова: самая быстрая память это регистры (быстрее shared), поэтому будем делать блоки с большим потреблением регистров, а тредов в блоке немного, 192 или 256. Occupancy будет маленькая, а общее быстродействие вырастет

Если речь идёт об GEMM-подобных задачах (где считанные данные могут использоваться много раз), то этот подход называется "Register blocking", в русской литературе я встречал название "блокирование регистров". Есть ещё cache blocking, можно делать многоуровневое блокирование и т.п.
Ещё в 1991г. ( http://www.cse.iitd.ernet.in/~nagaraju/TA/2004-2005/SDS-04/project/loop_... ) писали: "Blocking is a well-known optimization technique for improving the effectiveness of memory hierarchies"

Не обязательно GEMM,

Не обязательно GEMM, собственно вот статья: http://www.cs.berkeley.edu/~volkov/volkov10-PMAA.pdf

Для GEMM как раз все понятно, там увеличивая размер (квадратного) кэша в N раз получаешь выигрыш по отношению ops/memory access в те же N раз (т.к. фетчей N^2, а операций с этим блоком N^3).

Есть альтернативный подход

Есть альтернативный подход Василия Волкова: самая быстрая память это регистры (быстрее shared), поэтому будем делать блоки с большим потреблением регистров, а тредов в блоке немного, 192 или 256. Occupancy будет маленькая, а общее быстродействие вырастет

Просмотрел статью: имхо там суть отнюдь не в регистрах, а в том, что можно съесть все lantency двумя варпами(двумя, как я понял потому, что один SM реально одновременно выполняет только два варпа) за счёт Instruction Level Parallelism(то есть за счёт того, что варп активен, когда на очереди инструкция данные для которой уже готовы). Большее число доступных регистров на поток это уже побочное следствие из возможности съедать lantency двумя варпами на SM.
Если же съесть все latency небольшим количеством варпов на SM не получается, соответственно при большом использовании регистров, то может быть ухудшение производительности..

Тут где у этой лошади голова,

Тут где у этой лошади голова, а где - хвост - зависит от интерпретации.

Латентность можно съесть двумя варпами. Отлично.
А зачем это захотелось бы делать? А затем, что в каждом варпе регистров тогда будет доступно просто дохрена.

Если бы плюсом было только

Если бы плюсом было только большее количество регистров, то да, зависит только от интерпретации.
Но, имхо, могут быть ситуации где выигрыш будет не только из-за большего количества регистров. Например если есть задача с небольшим количеством потоков(прям 15*32*2 на gtx480), которые разумно разделить не получается: едят эти потоки много данных(достаточное для скрытия латентности), каждый элемент входных данных используется только для вычислений промежуточных результатов(за счёт этого, требования к количеству регистров не большое), потоки вычислений сильно переплетаются (за счёт этого переплетения разбивать потоки не разумно). Применить статью Волкова с позиции register blocking не получится - поэтому лучше помнить с позиции сокрытия латентности. Если не применить эту статью, и следовать идеи повышения occupancy, то часть потоков будет считывать данные в shared, для основных вычислительных потоков. В то время, как применяя идеи статьи, можно смело использовать небольшое количество потоков, без shared памяти.
Но, конечно, пример этот сильно надуманный...

Представим, что все как вы

Представим, что все как вы написали
- 15 блоков по 64 потока
- регистров надо мало
- shared надо мало
- латентность прячется в данных
Т.е. задача - CPU bound в наичистейшем виде.

Получится забавно: нам надо, чтобы планировщик эти блоки распихал по одному на SM. Но на ферми он делать этого не будет, наоборот запихает все на один, потому что вдруг придет еще один kernel, тут его параллельно и выполнят.

Т.е. чтобы было счастье - нужно фейково затребовать больше shared (или регистров) чем на самом деле надо, чтобы планировщик не чудил и пихал по одному блоку на SM.

Т.е. чтобы было счастье -

Т.е. чтобы было счастье - нужно фейково затребовать больше shared (или регистров) чем на самом деле надо, чтобы планировщик не чудил и пихал по одному блоку на SM

Да, у Волкова это называется "Allocate shared memory dynamiclly to control occupancy".
Получается неявный контроль планировщика, что можно назвать workaround'ом.. также желательно как-то проверять количество shared/регистров на SM, а то вдруг на новых Compute Capability больше будет..

Т.е. задача - CPU bound в

Т.е. задача - CPU bound в наичистейшем виде
А почему не может быть Memory bound при этих условиях?

Что вы понимаете под latency?

Что вы понимаете под latency? как оно соотносится с throughput/bandwidth?
Что Волков скрывал(всмысле memory latency,) во втором примере (со страницы 11)? Что он получил? Разве этот пример перестал быть memory bound?

Ну я скорее о том, что memcpy

Ну я скорее о том, что memcpy - без разницы как делать. Можно в 64 потока и 8 копирований в потоке, можно в 320 потоков и одно копирование, все одно упрешься в память (bandwidth)

Если другие ресурсы несущественны (регистры, shared), то разницы никакой.

Если другие ресурсы

Если другие ресурсы несущественны (регистры, shared), то разницы никакой

Я как раз пытаюсь привести некий сферический пример для которого есть разница (и для которого и регистры, и shared несущественны):

Пусть каждый поток считывает достаточное количество данных "n" для скрытия латентности всего 64мя потоками на SM.

Алгоритм по которому поток работает с данными не позволяет эффективного распараллеливания, например что-то типа
a(... b(c(d(e(f(X1,X2),X3),X4),X5),X6)... ,Xn)

Вычисления в любом случае будут идти в одном потоке.
Вопрос в том, будут ли эти "n" данных считываться только одним потоком (который и будет делать вычисления) или же их будут считывать несколько потоков и кормить вычислительный поток через shared память.

Ну то есть такой

Ну то есть такой искусственный случай, когда распараллелиться на 64*числоSM можно, а дальше - нет?

Ну, наверное, да.

Ну то есть такой

Ну то есть такой искусственный случай, когда распараллелиться на 64*числоSM можно, а дальше - нет?

Да, именно (но только в этом искусственном случае ещё важное свойство: данных одному потоку нужно много). В таких случаях полезно помнить об ILP, описанном в статье. Кстати, хорошая статья, спасибо за ссылку!