Open Source & Linux Lab

It's better when it's simple

User Tools

Site Tools


etc:common_activities:intel_students_cup:tour2

Differences

This shows you the differences between two versions of the page.

Link to this comparison view

Both sides previous revisionPrevious revision
Next revision
Previous revision
etc:common_activities:intel_students_cup:tour2 [2007/10/31 23:22] jcmvbkbcetc:common_activities:intel_students_cup:tour2 [2008/01/03 02:32] (current) – external edit 127.0.0.1
Line 1: Line 1:
 ====== Исследование вопроса ====== ====== Исследование вопроса ======
-пока я выложу свои изыскания. если появятся альтернативы -- сделаем разные странички, или секции -- как будет удобнее.+ 
 +  * [[tour2_tasks|Текущие задачи]] 
 +  * [[algo_details|Детали алгоритма]]
  
 ===== План ===== ===== План =====
Line 150: Line 152:
  LIBS_ALL   = $(LIBS)  LIBS_ALL   = $(LIBS)
 </file> </file>
 +
  
 ==== Менее очевидные шаги ==== ==== Менее очевидные шаги ====
 +=== Использование Math Kernel Library ===
 +[[http://www.intel.com/cd/software/products/asmo-na/eng/307757.htm|Math Kernel Library]]([[http://www.google.ru/url?sa=t&ct=res&cd=1&url=http%3A%2F%2Fwww.intel.com%2Fsoftware%2Fproducts%2Fmkl%2Ftechtopics%2Fmklman52.pdf&ei=txIpR7fmGIfs0AS7oMj7Aw&usg=AFQjCNGAKhvQ7sI3I9hp7BI6QU5mxa-RUg&sig2=ywkFBVb2LBI5chtUSOhN3Q|Описание]])
 +  * библиотека в основном хорошо реализует матричные операции;
 +Для исследования Intel Kernel Math Library написал следующее:
 +
 +Код "глупой" программы, которая вызывает не думая стандартный синус.
 +<code cpp> 
 +#include <stdlib.h>
 +#include <math.h>
 +
 +int main()
 +{
 +  double X[32][32];
 +  double F[32][32];
 +
 +  for(int i=0; i<32; ++i)
 +    for(int j=0; j<32; ++j)
 +      X[i][j] = rand()%1024;
 +
 +  for(int x=0;x<800;++x)
 +    for(int y=0;y<600;++y)
 +      for(int i=0;i<32;++i)
 +        for(int j=0;j<32;++j)
 +           F[i][j] = sinf(X[i][j]);
 +  return 0;
 +}
 +</code>
 +
 +Код программы, с использованием MKL:
 +В нем применил функцию, вычисляющую синус элементов вектора.
 +<code cpp>
 +#include "mkl.h"
 +#include <stdlib.h>
 +
 +int main()
 +{
 +  double X[32][32];
 +  double F[32][32];
 +
 +  for(int i=0; i<32; ++i)
 +    for(int j=0; j<32; ++j)
 +      X[i][j] = rand()%1024;
 +
 +  for(int x=0;x<800;++x)
 +    for(int y=0;y<600;++y)
 +      vdSin(32*32,(const double *)X,(double *)F);
 +  return 0;
 +}
 +</code>
 +
 +MakeFile:
 +<code>
 +default: stupid fast
 +
 +main.o: main.cpp
 +        g++ main.cpp -c -o main.o
 +
 +stupid: main.o
 +        g++ main.o -o stupid
 +clean:
 +        rm -f main.o stupid
 +fast: imkl_main.o
 +        g++ -L/opt/intel/mkl/9.1.023/lib/32 imkl_main.o -lguide -lmkl_p4m -lmkl_ia32 -lm -lirc -o fast
 +
 +imkl_main.o: imkl_main.cpp
 +        g++ -I/opt/intel/mkl/9.1.023/include -c imkl_main.cpp -o imkl_main.o
 +</code>
 +Третьим шагом было изменение в "быстрой" программе всех double на float. и вызов функции vsSin
 +Результат запуска:
 +<code bash>
 +make && time ./stupid && time ./fast && time ./floatfast
 +
 +real    0m41.288s
 +user    0m39.579s
 +sys     0m0.143s
 +
 +real    0m18.878s
 +user    0m18.158s
 +sys     0m0.086s
 +
 +real    0m7.799s
 +user    0m7.412s
 +sys     0m0.049s
 +</code>
 +
 +помимо всего прочего MKL имеет реализацию __одновременного__ вычисления синуса и косинуса в одной функции
 +
 +Очевидно, что если в "тупую" программу добавить рядом с вызовом синуса вызов косинуса, то врямя возрастет в два раза, что и произошло при опытной проверке.
 +
 +В случае же с MKL, интересней. Далее сравнительные времена выполнения двух программ с применением MKL:
 +  - вызываются vsSin и vsCos 
 +  - вызывается vsSinCos
 +<code>
 +real    0m16.126s
 +user    0m15.261s
 +sys     0m0.026s
 +
 +real    0m13.289s
 +user    0m12.670s
 +sys     0m0.028s
 +</code>
 +
 +=== Использование Intel C Compiler ===
 +[[http://www.intel.com/cd/software/products/asmo-na/eng/284132.htm|Компилятор]]([[http://www.intel.com/cd/software/products/asmo-na/eng/compilers/clin/277618.htm|Описание]])
 +  * должен позволить автоматически ипользовать SIMD команды (SSE, SSE2...) для оптимизации вычислений в основном в циклах;
 +  * возможно подскажет где что можно ещё распараллелить;
 +
 +[[icc|Результаты]]
 +
 +=== Использование Integrated Performance Primitives ===
 +[[http://www.intel.com/cd/software/products/asmo-na/eng/302910.htm|Intel Performance Primitives]]({{intel:sw5r.pdf|Описание}})
 +  * как заявляется производительность растёт в том числе и за счёт оптимизации библиотеки под различные модели процессоров;
 +  * можно попробовать исопльзовать оттуда не только тригонометрию но и функции работы с изображениями 2D;
 +P.S. Менеджер проекта этой библиотеки из Нижнего Новгорода откуда и сам sunset :)
 +
 +тестовая программка, а-ля zps:
 +<code cpp>
 +#include <stdlib.h>
 +#ifdef USE_IPP
 +#include <ippvm.h>
 +#endif
 +#include <math.h>
 +#include <stdio.h>
 +
 +#define COUNT(a) (sizeof(a)/sizeof(*(a)))
 +
 +float a[1024];
 +float r1[1024];
 +float r2[1024];
 +
 +int main()
 +{
 +  for(size_t i=0;i<COUNT(a);++i)
 +    a[i]=(drand48()-.5)*20;
 +
 +#ifdef USE_IPP
 +  for(int i=0;i<800*600;++i)
 +    ippsSin_32f_A21(a,r1,COUNT(r1));
 +#else
 +  for(int i=0;i<800*600;++i)
 +    for(size_t j=0;j<COUNT(r1);++j)
 +      r2[j]=sinf(a[j]);
 +#endif
 +
 +#if 0
 +  double s=0;
 +  for(size_t i=0;i<COUNT(r1);++i)
 +  {
 +    s+=fabs(r1[i]-r2[i]);
 +  }
 +  printf("%lg\n",s);
 +#endif
 +}
 +</code>
 +<code>
 +default: stupid fast
 +
 +stupid: main.o
 +  g++ main.o -o stupid
 +clean:
 +  rm -f main.o stupid
 +fast: ipp_main.o
 +  g++ -L/opt/intel/ipp/5.2/ia32/sharedlib ipp_main.o -lippcore -lippvm -o fast
 +
 +ipp_main.o: main.cpp
 +  g++ -I/opt/intel/ipp/5.2/ia32/include -DUSE_IPP -c main.cpp -o ipp_main.o
 +</code>
 +
 +результаты:
 +<code>
 +$ time ./fast && time ./stupid
 +
 +real    0m3.666s
 +user    0m3.661s
 +sys     0m0.002s
 +
 +real    0m37.244s
 +user    0m37.095s
 +sys     0m0.046s
 +</code>
 +итого -- ускорение в 10 раз. замена A21 на A11 дает
 +<code>
 +$ time ./fast && time ./stupid
 +
 +real    0m2.975s
 +user    0m2.964s
 +sys     0m0.007s
 +
 +real    0m36.754s
 +user    0m36.675s
 +sys     0m0.033s
 +</code>
  
-  - ускорить счет синуса. например линейной интерполяцией по таблице:+Для тригонометрии в hotspot 2 использовал ippsSin.
 <file> <file>
-diff -ruN src-org/sunset.cpp src/sunset.cpp+diff -ruN src-org/Makefile src1/Makefile 
 +--- src-org/Makefile    2007-09-17 17:43:08.000000000 +0400 
 ++++ src1/Makefile       2007-11-02 00:49:53.000000000 +0300 
 +@@ -56,19 +56,19 @@ 
 + CC         = gcc 
 + CPLAT      = 
 + CPROC      = 
 +-CINC       = -I$(SRC_DIR) 
 ++CINC       = -I$(SRC_DIR) -I/opt/intel/ipp/5.2/ia32/include 
 + CDEFS      = 
 + COBJ       = -c -o$(OBJ_DIR)/$@ 
 +-CDEFOPT    = -O2 
 ++CDEFOPT    = -g -pg -fprofile-arcs -ftest-coverage 
 + COPT       = 
 +-CFLAGS     = 
 ++CFLAGS     = -O3 -ffast-math -ffinite-math-only -fno-math-errno -funsafe-math-optimizations -fno-trapping-math -march=prescott 
 + CFLAGS_ALL = $(CFLAGS) $(CINC) $(CDEFS) $(CDEFOPT) $(CPROC) $(CPLAT) 
 + 
 + LD         = g++ 
 + LDPLAT     = 
 +-LDFLAGS    = 
 ++LDFLAGS    = -g -pg -ax -fprofile-arcs -ftest-coverage -L/opt/intel/ipp/5.2/ia32/sharedlib 
 + LDOUTOPT   = -o "$(OUT_DIR)/$(BENCHMARK)" 
 +-LIBS       = -lm -lc 
 ++LIBS       = -lm -lc -lippcore -lippvm 
 + LIBS_ALL   = $(LIBS) 
 + 
 + endif 
 +diff -ruN src-org/sunset.cpp src1/sunset.cpp
 --- src-org/sunset.cpp  2007-09-16 12:04:44.000000000 +0400 --- src-org/sunset.cpp  2007-09-16 12:04:44.000000000 +0400
-+++ src/sunset.cpp  2007-10-24 19:39:04.000000000 +0400 ++++ src1/sunset.cpp     2007-11-02 00:34:38.000000000 +0300 
-@@ -84,6 +84,27 @@ +@@ -45,6 +45,@@ 
-     return a * (float)k; + #include <omp.h> 
- }+ #endif 
 + #include "sunset.h" 
 ++#include <ippvm.h>
  
-+#define SIN_TAB_SZ 1024 + #define MIN(x,y   (((x(y)) ? (x) : (y)) 
-+#define PI2 (2*3.141592653589793f) + #define MAX(x,y   (((x(y)) ? (y) : (x)
-+float g_sinTab[2*SIN_TAB_SZ+1]; +@@ -747,+748,10 @@
-+
-+void fillSinTab() +
-+{ +
-+  for(int i=0;i<2*SIN_TAB_SZ+1;++i) +
-+  { +
-+    g_sinTab[i]=sinf((i-SIN_TAB_SZ)*PI2/SIN_TAB_SZ); +
-+  } +
-+} +
-+
-+inline float tab_sinf(float v) +
-+{ +
-+  float i=fmodf(v,PI2)*SIN_TAB_SZ/PI2+SIN_TAB_SZ; +
-+  size_t idx=i; +
-+  float d=i-idx; +
-+
-+  return g_sinTab[idx]*(1-d)+g_sinTab[idx+1]*d; +
-+} +
-+
- /* memory allocation function for two-dimensional arrays */ +
- void** alloc2D(int M, int N, size_t size) +
- { +
-@@ -320,+341,@@ +
- /*---------------------------------------------------------------------------*/ +
-     if(iCurFrame == 1) +
-     { +
-+        fillSinTab(); +
-+
-         iWaveMeshSize = iWaveHarmNum * iAngleHarmNum;+
  
-         if(iAllocated == 1) +                 pFlTmp flArgSin[currentthread].aptr;
-@@ -749,7 +772,7 @@+
  
-                 #pragma ivdep +-                #pragma ivdep 
-                 for(t=0; t<iWaveMeshSize; t++)+-                for(t=0; t<iWaveMeshSize; t++)
 -                    pFlTmp[t] = (float)sinf(pFlTmp[t]); -                    pFlTmp[t] = (float)sinf(pFlTmp[t]);
-                   pFlTmp[t] = (float)tab_sinf(pFlTmp[t]);+               ippsSin_32f_A21(pFlTmp,pFlTmp,iWaveMeshSize); 
 ++                //#pragma ivdep 
 ++                //for(t=0; t<iWaveMeshSize; t++) 
 ++                //    pFlTmp[t] = (float)tab_sinf(pFlTmp[t]);
  
                  /* initialize the values of derivation */                  /* initialize the values of derivation */
                  flDerivX = 0.0f;                  flDerivX = 0.0f;
 </file> </file>
 +Результат -- 10.044/кадр, 0.7% отличий. Понижение точности до 11 бит дает 9.139/кадр, 1% отличий. Однако, по-простецки с openmp оно дружить не захотело -- segfault.
 +
 +
 +zps.
 +Максимальный результат:
 +diff -ruN src-org src-ipp
 +<code>
 +diff src-org/Makefile src-ipp/Makefile
 +59c59
 +< CINC       = -I$(SRC_DIR)
 +---
 +> CINC       = -I$(SRC_DIR) -I/opt/intel/ipp/5.2/ia32/include
 +62c62,63
 +< CDEFOPT    = -O2
 +---
 +> #CDEFOPT    = -O2
 +> CDEFOPT    = -g
 +64c65
 +< CFLAGS     =
 +---
 +> CFLAGS     = -O3 -ffast-math -ffinite-math-only -fno-math-errno -funsafe-math-optimizations -fno-trapping-math -march=pentium4m 
 +69c70
 +< LDFLAGS    = 
 +---
 +> LDFLAGS    = -g -pg -ax -fprofile-arcs -ftest-coverage -L/opt/intel/ipp/5.2/ia32/sharedlib
 +71c72
 +< LIBS       = -lm -lc
 +---
 +> LIBS       = -lm -lc -lippcore -lippvm -lguide -lipps -lippm
 +diff src-org/sunset.cpp src-ipp/sunset.cpp
 +48c48,51
 +
 +---
 +> #include <ippvm.h>
 +> #include <ipps.h>
 +> #include <ippm.h>
 +> #include <iostream>
 +307c310,311
 +<                 float   OT, KX1, KY1;
 +---
 +>                 float   *OT, *KX1, *KY1;
 +>                 float   kx1, ky1, ot; // zps
 +339a344,346
 +>             free(OT);
 +>             free(KX1);
 +>             free(KY1); 
 +357c364,366
 +
 +---
 +>         OT           = (float*)malloc(iWaveHarmNum * sizeof(float));     
 +>         KX1          = (float*)malloc(iWaveHarmNum * sizeof(float));     
 +>         KY1          = (float*)malloc(iWaveHarmNum * sizeof(float));     
 +362a372
 +>         
 +657a668,669
 +>     ippsMulC_32f(flOmega, flTime, OT, iWaveHarmNum);
 +>     float xxx[iAngleHarmNum];
 +663c675
 +<         private(currentthread, OT, KX1, KY1) \
 +---
 +>         private(currentthread, OT, KX1, KY1, kx1, ky1, ot) \
 +679a692
 +>     
 +732c745,749
 +< */
 +---
 +> *///              std::cerr << "NKMAX: "<< NKMAX << std::endl;
 +>                  
 +>                 ippsMulC_32f(flK, flDecartX[i][j], KX1, NKMAX);
 +>                 ippsMulC_32f(flK, flDecartY[i][j], KY1, NKMAX);
 +
 +735,738c752,755
 +<                     OT  = flOmega[t] * flTime;
 +<                     KX1 = flK[t] * flDecartX[i][j];
 +<                     KY1 = flK[t] * flDecartY[i][j];
 +
 +---
 +>                     kx1 = KX1[t];
 +>                     ky1 = KY1[t];
 +>                     ot =  OT[t];
 +>                     int len = t * iAngleHarmNum;
 +741,745c758,761
 +<                         iSinIndex1 = t * iAngleHarmNum + l;
 +<                         flArgSin[currentthread].aptr[iSinIndex1] = OT - 
 +<                             KX1 * flAzimuthCosFi[l] - KY1 * flAzimuthSinFi[l] +
 +<                             flRandomPhase[t*iAngleHarmNum + l];
 +<                     } /* end for l */
 +---
 +>                         iSinIndex1 = len + l;
 +>                         flArgSin[currentthread].aptr[iSinIndex1] = ot -
 +>                            kx1*flAzimuthCosFi[l] - ky1*flAzimuthSinFi[l] + flRandomPhase[iSinIndex1];
 +>                    } /* end for l */
 +747c763,776
 +
 +---
 +> #if 0
 +>                 float * dest;
 +>                 for(t = 0; t < NKMAX; t++)
 +>                 
 +>                     dest = &flArgSin[currentthread].aptr[t*iAngleHarmNum];
 +>                     kx1 = KX1[t];
 +>                     ky1 = KY1[t];
 +>                     ippmLComb_vv_32f(flAzimuthCosFi, 4,kx1, flAzimuthSinFi, 4,ky1, xxx, 4, iAngleHarmNum);
 +>                     ippsSubCRev_32f(xxx,OT[t],dest, iAngleHarmNum);
 +>                     
 +>                 } /* end for t */
 +>                 dest = NULL;
 +> #endif
 +> //                ippsAdd_32f_I(flRandomPhase, flArgSin[currentthread].aptr, iWaveMeshSize);
 +750,752c779,783
 +<                 #pragma ivdep
 +<                 for(t=0; t<iWaveMeshSize; t++)
 +<                     pFlTmp[t] = (float)sinf(pFlTmp[t]);
 +---
 +> ippsSin_32f_A21(pFlTmp,pFlTmp,iWaveMeshSize);
 +
 +>                 //#pragma ivdep
 +>                 //for(t=0; t<iWaveMeshSize; t++)
 +>                 //    pFlTmp[t] = (float)sinf(pFlTmp[t]);
 +758a790,801
 +> #if 0
 +>                           
 +>                   const float *pr[]={flAmplitudeX,flAmplitudeY};
 +>                   float __r[2];
 +
 +>                   ippmDotProduct_vav_32f_L(pr,0,4,pFlTmp,4,__r,iWaveMeshSize,2);
 +
 +>                   flDerivX = __r[0];
 +>                   flDerivY = __r[1];
 +>                 }
 +> #endif
 +> #if 1
 +767a811
 +> #endif
 +Только в src-ipp: .sunset.cpp.swp
 +</code>
 +
 +результат:
 +<code bash>
 +Frame 1 of 16 ... frame time 8.534
 +Frame 2 of 16 ... frame time 8.508
 +Frame 3 of 16 ... frame time 8.167
 +Frame 4 of 16 ... frame time 8.538
 +Frame 5 of 16 ... frame time 8.467
 +Frame 6 of 16 ... frame time 8.302
 +Frame 7 of 16 ... frame time 8.565
 +Frame 8 of 16 ... frame time 8.425
 +Frame 9 of 16 ... frame time 8.549
 +Frame 10 of 16 ... frame time 8.398
 +Frame 11 of 16 ... frame time 8.425
 +Frame 12 of 16 ... frame time 8.595
 +Frame 13 of 16 ... frame time 8.419
 +Frame 14 of 16 ... frame time 8.312
 +Frame 15 of 16 ... frame time 8.471
 +Frame 16 of 16 ... frame time 8.417
 +=================================
 +Timing:
 +        Total time is 135.098 sec., average frame time is 8.444 sec.
 +
 +Correctness check:
 +        Max RGB difference is 11.
 +        Number of different color pixels is 3564 (0.7%).
 +
 +</code>
 +
 +{{:intel:result_kel.txt|Результаты работы на конфигурации P4(2,6)Windows}}:\\
 +Лучший из полученных результатов:
 +<code bash>
 +Frame 1 of 16 ... frame time 3.591
 +Frame 2 of 16 ... frame time 3.492
 +Frame 3 of 16 ... frame time 3.497
 +Frame 4 of 16 ... frame time 3.499
 +Frame 5 of 16 ... frame time 3.485
 +Frame 6 of 16 ... frame time 3.492
 +Frame 7 of 16 ... frame time 3.494
 +Frame 8 of 16 ... frame time 3.492
 +Frame 9 of 16 ... frame time 3.523
 +Frame 10 of 16 ... frame time 3.491
 +Frame 11 of 16 ... frame time 3.496
 +Frame 12 of 16 ... frame time 3.495
 +Frame 13 of 16 ... frame time 3.494
 +Frame 14 of 16 ... frame time 3.492
 +Frame 15 of 16 ... frame time 3.503
 +Frame 16 of 16 ... frame time 3.495
 +=================================
 +Timing:
 +        Total time is 56.036 sec., average frame time is 3.502 sec.
  
-Даже такая простецкая реализация дает время 23.724/кадрили ускорение в 1.5 раза. При этом количество отличающихся пикселов растет по сравнению с оригинальным алгоритмом незначительно. Если интерполяцию заменить просто поиском ближайшего узла таблицы, точность падает значительно (условие 5% отличных пикселов не выполняется). +Correctness check: 
-Полагаю, более тонкая (ассемблерная?) реализация синуса даст еще небольшой прирост скорости.+        Max RGB difference is 20. 
 +        Number of different color pixels is 5437 (1.1%). 
 +</code> 
 +Достигнут:\\ 
 +  указанием ключа fast при компилляции 
 +  заменой дублирующих вычислений в циклах в очевидных местах 
 +  выносом вычисления синусов и косинусов (ipp) полярных координат
  
-Откомпилировать другими компиляторамиесли времена для одного и того же кода будут заметно отличаться -- смотреть в ассемблер. или просто использовать другой компилятор и дальше (.+==== Совсем неочевидные/исследовательские шаги ==== 
 +  * [[tab_sin|Табличный синус]] (плавающие числа) -- тупиковая ветвь. 
 +  * Целочисленная реализация алгоритма расчета вектора нормали -- выигрыша по скорости нет, по точности -- проигрыш. Тупиковая ветвь
 +  * [[harm_reduction|Уменьшение количества гармоник]] -- тупиковая ветвь.
etc/common_activities/intel_students_cup/tour2.1193862176.txt.gz · Last modified: 2008/01/03 02:32 (external edit)