"38 попугаев" за 6000 евро - очередной супер-комп

http://www.regard.ru/catalog/tovar1718.htm
Оказывается это можно купить. И ответ на мой вопрос - действительно 7 PCIe 16x.

Что касается статей, я два года назад уже запускал перемножение матриц. Но тогда СUDA не было - я использовал продукт от rapidmind. Результат был не впечатляющий, но по-моему там всего 8 конвееров использовалось, хоть и "коммерческая" версия была. Потом были другие дела и забросил.

 

Нашел я подробности : http://fastra2.ua.ac.be/?page_id=214

А вот тут интересно :
http://www.dvhardware.net/articles25_fastra_2_desktop_supercomputer.html

 

W4FhLF

Что-то меня заинтересовало перемножение матриц, решил поумножать матрицы типа float на процессоре C2D T9400 2.53ГГц(L2 - 6Мб).
Исходный алгорит, откомпилирован в 2008 студии с оптимизацией по скорости:

Код (Text):

1. void line_mul3(float* d,float* s,float k)
2. {
3.     for(int i=0;i<N;i++)
4.         *d++ +=*s++*k;
5. }
6.  
7. void matrix_mul3(float* m3,float* m1,float* m2){
8.     memset(m3,0,N*N*sizeof(*m3));
9.     for(int i=0;i<N;i++){
10.         float* l1=m1;
11.         float* l2=m2;
12.         for(int k=0;k<N;k++){
13.             line_mul3(m3,l2,*l1);
14.             l1++;
15.             l2+=N;
16.         }
17.         m1+=N;
18.         m3+=N;
19.     }
20. }

Результаты для матриц 512x512 - 120 мс, 1024x1024 - 950 мс, 2048x2048 - 9 с, 4096x4096 - 71 c.

А вот алгоритм с привлечением SSE и разбиением внешних циклов(без разбиения начинает раза в два тормозить когда матрица не лезет в кеш):

Код (Text):

1. #define CACHE 0x80000
2. #define nk (CACHE/N)
3. #define ni (CACHE/N)
4.  
5. void matrix_mul5(float* m3,float* m1,float* m2){
6.     memset(m3,0,N*N*sizeof(*m3));
7.     __asm{
8.         mov edi,[m3];
9.         mov esi,[m1];
10.         mov ebx,[m2];
11.         mov ecx,N/ni;
12.     };
13. l1:
14.     __asm{
15.         push ebx;
16.         push ecx;
17.         push esi;
18.         mov ecx,N/nk;
19.     };
20. l2:
21.     __asm{
22.         push ecx;
23.         push esi;
24.         push edi;
25.         mov ecx,ni;
26.     };
27. l3:
28.     __asm{
29.         push ebx
30.         push ecx
31.         push esi
32.         mov ecx,nk/4;
33.     };
34. l5:
35.     __asm{
36.         movdqa xmm0,[esi]
37.         pshufd xmm3,xmm0,255
38.         pshufd xmm2,xmm0,170
39.         pshufd xmm1,xmm0,85
40.         pshufd xmm0,xmm0,0
41.         push ebx
42.         push ecx
43.         push edi
44.         mov ecx,-N*4
45.         sub ebx,ecx
46.         sub edi,ecx
47.     };
48. l4:
49.     __asm{
50.         movdqa xmm4,[ebx+ecx]
51.         movdqa xmm5,[ebx+ecx+N*4]
52.         mulps xmm4,xmm0
53.         movdqa xmm6,[ebx+ecx+2*N*4]
54.         mulps xmm5,xmm1
55.         movdqa xmm7,[ebx+ecx+3*N*4]
56.         mulps xmm6,xmm2
57.         mulps xmm7,xmm3
58.         addps xmm5,xmm4
59.         movdqa xmm4,[edi+ecx]
60.         addps xmm7,xmm6
61.         addps xmm5,xmm4
62.         addps xmm7,xmm5
63.         movdqa [edi+ecx],xmm7
64.         add ecx,4*4
65.         js l4
66.         pop edi
67.         pop ecx
68.         pop ebx
69.         add ebx,4*N*4
70.         add esi,4*4
71.         loop l5
72.         pop esi
73.         pop ecx
74.         pop ebx
75.         add esi,N*4;
76.         add edi,N*4;
77.         dec ecx;
78.         ja l3;
79.         pop edi;
80.         pop esi;
81.         pop ecx;
82.         add ebx,nk*N*4
83.         add esi,nk*4
84.         dec ecx
85.         ja l2
86.         pop esi
87.         pop ecx
88.         pop ebx
89.         add esi,ni*N*4;
90.         add edi,ni*N*4;
91.         dec ecx
92.         ja l1
93.     };
94. }

Результаты для матриц 512x512 - 30 мс, 1024x1024 - 230 мс, 2048x2048 - 1.9 с, 4096x4096 - 18 c, 8192x8192 - 150 c.

Данная версия алгорита быстрее исходного примерно в 4 раза. Правда если считать в MAC'ах/такт(умножение с накоплением), то скорость падает с почти 2 до 1.5 при увеличении размера матрицы. Вообще даже 2 MAC/такт как-то совсем не впечатляют. У кого есть идеи по оптимизации? Понятно что можно попробовать на несколько ядер разбить, только вот даже два потока перемножающих свои матрицы взаимно тормозят друг друга процентов на 20. Поэтому пока интересуют идеи применимые к одному ядру. Да и вообще любопытно понять сколько максимум MAC можно выжать из SSE.

 

о кучи
http://http.developer.nvidia.com/GPUGems/gpugems_part01.html
http://http.developer.nvidia.com/GPUGems2/gpugems2_part01.html
http://http.developer.nvidia.com/GPUGems3/gpugems3_part01.html
3 бесплатных вводных туториальных учебника по гпу/куда от невидии. с картинками.
можно скачать оффлайново метров по 15-18 будут. есть и так в сети

эхх, жалко, что на моей гф го 7150 попробовать не выйдет..

 

Black_mirror, у тебя вторая матрица типа уже транспонирована?

 

Кстати если начинать сегодня, то OpenCL пока ещё слишком сырой. Он очень близок к CUDA, но CUDA уже отлажена, есть сообщество, куча доков, отличный SDK, а по OpenCL ещё голяк. Поэтому лучше начинать с CUDA. Перейти потом на OpenCL труда не составит.

 

W4FhLF
У меня все матрицы хранятся по строкам. Вариант, когда одна из матриц транспонированная, пробовал перемножать таким алгоритмом: for(i) for(j) {s=0 for(k) {s+=a[k]*b[j][k] } c[j]=s}, но он оказался в полтора раза менее эффективным(даже если разбивать внешние циклы) и я на него забил. А сейчас попробовал перемножать алгоритмом Штрассена (всё что меньше 1024x1024 умножаю модифицированным matrix_mul3, с тем отличием что суммирую сразу по 4 строки на SSE, разбиение внешних циклов не делаю - и так в кеш лезет). Результаты для Штрассена: 1024x1024 - 200 мс, 2048x2048 - 1.4c, 4096x4096 - 10c. Вообще в алгоритме нужно сделать оптимизацию по использованию памяти, потому что для матриц 8192x8192 памяти мне уже не хватило, но если кому интересно, то это безумие вот:

Код (Text):

1. float *alloc(float**heap,int sz){
2.     float *p=*heap;
3.     *heap+=sz;
4.     return p;
5. }
6.  
7. void line_mul4(int sz,float* d,float* s,int l,float* k)
8. {
9.     __asm{
10.         mov eax,[k]
11.         movdqa xmm0,[eax]
12.         pshufd xmm3,xmm0,255
13.         pshufd xmm2,xmm0,170
14.         pshufd xmm1,xmm0,85
15.         pshufd xmm0,xmm0,0
16.         mov edi,[d]
17.         mov esi,[s]
18.         mov ebx,[l]
19.         mov ecx,[sz]
20.         imul ecx,-4
21.         sub edi,ecx
22.         sub esi,ecx
23.         lea edx,[esi+ebx]
24.         lea eax,[edx+ebx*2]
25.         add ebx,edx
26.     }
27. lab:
28.     __asm{
29.         movdqa xmm4,[esi+ecx]
30.         movdqa xmm5,[edx+ecx]
31.         mulps xmm4,xmm0
32.         movdqa xmm6,[ebx+ecx]
33.         mulps xmm5,xmm1
34.         movdqa xmm7,[eax+ecx]
35.         mulps xmm6,xmm2
36.         mulps xmm7,xmm3
37.         addps xmm5,xmm4
38.         addps xmm7,xmm6
39.         movdqa xmm4,[edi+ecx]
40.         addps xmm5,xmm7
41.         addps xmm4,xmm5
42.         movdqa [edi+ecx],xmm4
43.         add ecx,16
44.         js lab
45.     }
46. }
47.  
48. void matrix_mul4(int sz,float* m3,int d3,float* m1,int d1,float* m2,int d2){
49.     for(int i=0;i<sz;i++)
50.         memset(m3+i*d3,0,sz*sizeof(*m3));
51.     for(int i=0;i<sz;i++){
52.         float* l1=m1;
53.         float* l2=m2;
54.         for(int k=0;k<sz;k+=4){
55.             line_mul4(sz,m3,l2,d2,l1);
56.             l1+=4;
57.             l2+=4*d2;
58.         }
59.         m1+=d1;
60.         m3+=d3;
61.     }
62. }
63.  
64. void line_add(int sz,float* m3,float* m1,float* m2)
65. {
66.     __asm mov ecx,[sz];
67.     __asm mov ebx,[m1]
68.     __asm mov esi,[m2]
69.     __asm mov edi,[m3]
70.     __asm imul ecx,-4;
71.     __asm sub ebx,ecx;
72.     __asm sub esi,ecx;
73.     __asm sub edi,ecx
74. lab:
75.     __asm movdqa xmm0,[ebx+ecx];
76.     __asm movdqa xmm2,[ebx+ecx+16];
77.     __asm movdqa xmm1,[esi+ecx];
78.     __asm movdqa xmm3,[esi+ecx+16];
79.     __asm addps xmm0,xmm1;
80.     __asm addps xmm2,xmm3;
81.     __asm movdqa [edi+ecx],xmm0;
82.     __asm movdqa [edi+ecx+16],xmm2;
83.     __asm add ecx,32;
84.     __asm js lab;
85. }
86.  
87. void line_sub(int sz,float* m3,float* m1,float* m2)
88. {
89.     __asm mov ecx,[sz];
90.     __asm mov ebx,[m1]
91.     __asm mov esi,[m2]
92.     __asm mov edi,[m3]
93.     __asm imul ecx,-4;
94.     __asm sub ebx,ecx;
95.     __asm sub esi,ecx;
96.     __asm sub edi,ecx
97. lab:
98.     __asm movdqa xmm0,[ebx+ecx];
99.     __asm movdqa xmm2,[ebx+ecx+16];
100.     __asm movdqa xmm1,[esi+ecx];
101.     __asm movdqa xmm3,[esi+ecx+16];
102.     __asm subps xmm0,xmm1;
103.     __asm subps xmm2,xmm3;
104.     __asm movdqa [edi+ecx],xmm0;
105.     __asm movdqa [edi+ecx+16],xmm2;
106.     __asm add ecx,32;
107.     __asm js lab;
108. }
109.  
110. void matrix_add(int sz,float* m3,int d3,float* m1,int d1,float* m2,int d2){
111.     for(int i=0;i<sz;i++){
112.         line_add(sz,m3,m1,m2);
113.         m1+=d1;
114.         m2+=d2;
115.         m3+=d3;
116.     }
117. }
118.  
119. void matrix_sub(int sz,float* m3,int d3,float* m1,int d1,float* m2,int d2){
120.     for(int i=0;i<sz;i++){
121.         line_sub(sz,m3,m1,m2);
122.         m1+=d1;
123.         m2+=d2;
124.         m3+=d3;
125.     }
126. }
127.  
128. void matrix_mul8(int sz,float* m3,int d3,float* m1,int d1,float* m2,int d2,float*heap){
129.     if(sz<512){
130.         matrix_mul4(sz,m3,d3,m1,d1,m2,d2);
131.         return;
132.     }
133.     int h=sz/2;
134.     int half=h*h;
135.     float* a11=m1;
136.     float* a12=a11+h;
137.     float* a21=m1+d1*h;
138.     float* a22=a21+h;
139.     float* b11=m2;
140.     float* b12=b11+h;
141.     float* b21=m2+d2*h;
142.     float* b22=b21+h;
143.     float* c11=m3;
144.     float* c12=c11+h;
145.     float* c21=m3+d3*h;
146.     float* c22=c21+h;
147.     float* s1=alloc(&heap,half);
148.     float* s2=alloc(&heap,half);
149.     float* s3=alloc(&heap,half);
150.     float* s4=alloc(&heap,half);
151.     float* s5=alloc(&heap,half);
152.     float* s6=alloc(&heap,half);
153.     float* s7=alloc(&heap,half);
154.     float* s8=alloc(&heap,half);
155.     float* p1=alloc(&heap,half);
156.     float* p2=alloc(&heap,half);
157.     float* p3=alloc(&heap,half);
158.     float* p4=alloc(&heap,half);
159.     float* p5=alloc(&heap,half);
160.     float* p6=alloc(&heap,half);
161.     float* p7=alloc(&heap,half);
162.     float* t1=alloc(&heap,half);
163.     float* t2=alloc(&heap,half);
164.     float* t3=alloc(&heap,half);
165.     matrix_add(h,s1,h,a21,sz,a22,sz);
166.     matrix_sub(h,s2,h,s1,h,a11,sz);
167.     matrix_sub(h,s3,h,a11,sz,a21,sz);
168.     matrix_sub(h,s4,h,a12,sz,s2,h);
169.     matrix_sub(h,s5,h,b12,sz,b11,sz);
170.     matrix_sub(h,s6,h,b22,sz,s5,h);
171.     matrix_sub(h,s7,h,b22,sz,b12,sz);
172.     matrix_sub(h,s8,h,s6,h,b21,sz);
173.     matrix_mul8(h,p1,h,s2,h,s6,h,heap);
174.     matrix_mul8(h,p2,h,a11,sz,b11,sz,heap);
175.     matrix_mul8(h,p3,h,a12,sz,b21,sz,heap);
176.     matrix_mul8(h,p4,h,s3,h,s7,h,heap);
177.     matrix_mul8(h,p5,h,s1,h,s5,h,heap);
178.     matrix_mul8(h,p6,h,s4,h,b22,sz,heap);
179.     matrix_mul8(h,p7,h,a22,sz,s8,h,heap);
180.     matrix_add(h,t1,h,p1,h,p2,h);
181.     matrix_add(h,t2,h,t1,h,p4,h);
182.     matrix_add(h,t3,h,p5,h,p6,h);
183.     matrix_add(h,c11,sz,p2,h,p3,h);
184.     matrix_add(h,c12,sz,t1,h,t3,h);
185.     matrix_sub(h,c21,sz,t2,h,p7,h);
186.     matrix_add(h,c22,sz,t2,h,p5,h);
187. }

 

GPU уже в сторону DB шагает.

A reimplementation of portions of the SQLite database to execute on a GPU, part of the GPGPU-3 workshop. Speed Up 70x.

CyberManiac, привет