Современные информационные технологии/Вычислительная техника и программирование

Мясищев А.А.

Хмельницкий национальный университет, Украина

ОЦЕНКА ЦЕЛЕСООБРАЗНОСТИ ИСПОЛЬЗОВАНИЯ ВИДЕОКАРТЫ NVIDIA ПО СРАВНЕНИЮ С 6-И ЯДЕРНЫМ ПРОЦЕССОРОМ AMD ДЛЯ МАТРИЧНОГО УМНОЖЕНИЯ

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

Сопоставим эффективность использования для расчета произведения квадратных матриц для одинарной точности вычислительной системой, в которой установлены 6-и ядерный процессор AMD Phenom II X6 1090T и видеокарта NVIDIA GeForce GTX 480. Причем в первом случае расчет будет выполняться распараллеливание по 6-ядрам процессора с использованием библиотек ScaLAPACK и библиотек ATLAS (автоматически настраиваемое программное обеспечение для решения задач линейной алгебры). Библиотеки ATLAS при компиляции настраиваются под конкретную архитектуру процессора вычислительной системы. Во втором случае  расчет будет выполняться на видеокарте NVIDIA GeForce GTX 480 с использованием технологии CUDA[1]. Вычислительная система работает под управлением операционной системы Linux Ubuntu ver. 10.04 desktop. Здесь установлены компилятор фортран  F77, библиотеки MPI [2], ScaLAPACK[3] и ATLAS[4] для 6-и ядерного процессора. Для программирования на NVIDIA GeForce GTX 480 инсталлируется видеодрайвер nvidia и программное обеспечение  с сайта http://developer.nvidia.com/cuda-toolkit-archive. Последовательность установки программного обеспечения представлена в источнике [5].

    Рассмотрим установку библиотеки ATLAS. Предполагается, что библиотека ScaLAPACK установлена, как представлено в работе [6]. С сайта http://sourceforge.net/projects/math-atlas/files/ копируем последнюю версию файла atlas3.9.51.tar.bz2. Разархивируем его командами:

bzip2 –d atlas3.9.51.tar.bz2

tar xf atlas3.9.51.tar

Далее переходим в каталог ATLAS

cd ATLAS

Создаем рабочий каталог my:

mkdir my

Переходим в него

cd my

Запускаем конфигурационный скрипт

../configure -Si cputhrchk 0 -b 64 -Si archdef 1

Здесь

-Si cputhrchk 0 - игнорирование изучения работы процессора в режиме дросселирования,

-b 64 -  компиляция для 64-х разрядной операционной системы,

-Si archdef 1  -  включение режима обнаружения архитектуры процессора по умолчанию.

Запускаем компиляцию библиотеки

make

переходим в каталог lib

cd lib

и копируем необходимые скомпилированные библиотеки в  рабочий каталог

cp libatlas.a  /home/alex/

cp libf77blas.a  /home/alex/

В рабочем каталоге должны находиться библиотеки scalapack, blacs и blas. Компиляция программы выглядит следующим образом:

f77  O2  I/usr/local/ch_shmem/include  f  o  lab_sca_s  lab_sca_s.f scalapack_LINUX.a blacsF77init_MPI-LINUX-0.a  blacs_MPI-LINUX-0.a blacsF77init_MPI-LINUX-0.a  libf77blas.a libatlas.a /usr/local/ch_shmem/lib/libmpich.a

Здесь  lab_sca_s.f  - текст программы на фортране.

Запуск на исполнение на 6-и ядрах:

/usr/local/ch_shmem/bin/mpirun –np 6 -machinefile ../machine lab_sca_s

Файл machine имеет вид:

78.152.183.39:6

78.152.183.39 - ip адрес компьютера, 6 - количество ядер.

         Ниже представлена программа lab_sca_s.f перемножения квадратных матриц для чисел с одинарной точностью, которая была написана в соответствии с материалами работы [7], работающая  на многоядерном процессоре.

      program mult_matrix

      include 'mpif.h'

      parameter (nm = 6144, nxn = nm*nm)

      real a(nxn), b(nxn), c(nxn), mem(nm)

      real aa(nm,nm),bb(nm,nm)

      double precision time(6), ops, total, t1

      INTEGER DESCA(9), DESCB(9), DESCC(9)

      CALL BLACS_PINFO( IAM, NPROCS )

      NPROW = INT(SQRT(REAL(NPROCS)))

      NPCOL = NPROCS/NPROW

      IF( IAM.EQ.0 ) THEN

      print *,'Input N and NB: '

      read *, N, NB

      END IF

      call MPI_BCAST(N, 1, MPI_INTEGER, 0, MPI_COMM_WORLD, ierr)

      call MPI_BCAST(NB,1, MPI_INTEGER, 0, MPI_COMM_WORLD, ierr)

      ops = (2.0d0*dfloat(n)-1)*dfloat(n)*dfloat(n)

      CALL BLACS_GET( -1, 0, ICTXT )

      CALL BLACS_GRIDINIT( ICTXT, 'Row-major', NPROW, NPCOL )

      CALL BLACS_GRIDINFO( ICTXT, NPROW, NPCOL, MYROW, MYCOL)

      IF( MYROW.GE.NPROW .OR. MYCOL.GE.NPCOL ) GOTO 500

      NP = NUMROC( N, NB, MYROW, 0, NPROW )

      NQ = NUMROC( N, NB, MYCOL, 0, NPCOL )

      CALL DESCINIT( DESCA, N, N, NB,NB,0,0,ICTXT, NP, INFO )

      CALL DESCINIT( DESCB, N, N, NB,NB,0,0,ICTXT, NP, INFO )

      CALL DESCINIT( DESCC, N, N, NB,NB,0,0,ICTXT, NP, INFO )

      lda = DESCA(9)

c initial data

      do 10 i=1,n

      do 10 j=1,n

      aa(i,j)=1.0*(i+2*j)

      bb(i,j)=1.0/aa(i,j)

10   continue

c initial matrix in processors

      call pmatgen(a,aa,DESCA,np,nq,b,bb,DESCB,nprow,npcol, myrow,mycol,nm)

      t1 = MPI_Wtime()

      CALL PSGEMM('N','N', N, N, N, 1.0, A, 1, 1, DESCA,

     $  B, 1, 1, DESCB, 0.0, C, 1, 1, DESCC)

      time(2) = MPI_Wtime() - t1

      if (IAM.eq.0) then

      end if

      CALL PSLAPRNT( 1, 1, C, n/256+1, n/128+1,DESCC, 0, 0, 'C', 6, MEM )

      CALL PSLAPRNT( 1, 1, C, 3*n/4+1, 5*n/16+1,DESCC, 0, 0, 'C', 6, MEM )

      CALL PSLAPRNT( 1, 1, C, n-4+1, n-2+1, DESCC, 0, 0, 'C', 6, MEM )

      total=time(2)

      time(4)=ops/(1.0d6*total)

      if (IAM.EQ.0) then

      write(6,110) time(2), time(4)

110   format(1x,'Time = ',f12.4, ' sec.\n',

     $         ' Mflops = ',f14.4)

      end if

      CALL BLACS_GRIDEXIT( ICTXT )

      CALL BLACS_EXIT(0)

500   continue

      stop

      end

      subroutine pmatgen(a,aa,DESCA,np,nq,b,bb,DESCB,nprow,

    &npcol,myrow,mycol,nm)

      integer i, j, DESCA(*), DESCB(*), nprow, npcol,myrow,mycol

      real a(*), b(*),aa(6144,6144),bb(6144,6144)

      nb = DESCA(5)

      k = 0

      do 250 j = 1,nq

      jc = (j-1)/nb

      jb = mod(j-1,nb)

      jg = mycol*nb + jc*npcol*nb + jb + 1

      do 240 i = 1,np

      ic = (i-1)/nb

      ib = mod(i-1,nb)

      ig = myrow*nb + ic*nprow*nb + ib + 1

      k = k + 1

      a(k) = aa(ig,jg)

      b(k) = bb(ig,jg)

240   continue

250   continue

      return

      end

Рассмотрим программы расчета произведения квадратных заполненных матриц для  NVIDIA GeForce GTX 480(GPU). В первом случае для простоты примем, что элементы матриц хранятся в глобальной памяти GPU. Для составления даже этой простой программы необходимо подробнее рассмотреть программную модель GPU для компилятора Си. Верхний уровень ядра, называемый сеткой(grid) состоит из блоков. В свою очередь блоки представляют собой либо одноразмерную, либо двухразмерную сеть нитей.  Вышесказанное может быть  проиллюстрировано рис.1. Номер нити в блоке или номер блока в grid синтаксически специфицируется оператором <<<…>>> и может быть переменной типа int или dim3.

Рис.1

Число нитей, входящих в блок определяется встроенной переменной blockDim. Индекс нити внутри блока определяется переменной threadIdx, а индекс блока внутри grid – переменной blockIdx. Нити в блоке являются непосредственными исполнителями вычислений.  Нить является 3-х компонентным вектором, т.е. может идентифицироваться, используя одно размерный, двух размерный и трех размерный индекс нити, образуя в свою очередь одноразмерный, двух размерный и трех размерный блок нитей. Существует ограничение количества нитей на один блок, которое не может превышать 1024 нити.

         Расширение Си  CUDA позволяет вызвать функцию, называемую ядром так, что она будет параллельно выполнять N разных нитей CUDA. Такие функция декларируется спецификатором __global__  .  В качестве иллюстрации ниже рассмотрен пример сложения двух векторов размером N с сохранением результата в векторе C.

// Kernel definition

__global__ void VecAdd(float* A, float* B, float* C)

{

int i = threadIdx.x;

C[i] = A[i] + B[i];

}

int main()

{

...

// Kernel invocation with N threads

VecAdd<<<1, N>>>(A, B, C);

}

Здесь каждая из N нитей, которая выполняется функцией VecAdd(), производит только одну пару сложения. Задействован только один блок. Если используется множество блоков, каждый из которых  состоит из 16x16 нитей и матрицы являются двухмерными, то код программы будет выглядеть несколько сложнее:

// Kernel definition

__global__ void  MatAdd(float A[N][N], float B[N][N],

float C[N][N])

{

int i = blockIdx.x * blockDim.x + threadIdx.x;

int j = blockIdx.y * blockDim.y + threadIdx.y;

if (i < N && j < N)

C[i][j] = A[i][j] + B[i][j];

}

int main()

{

...

// Kernel invocation

dim3 threadsPerBlock(16, 16);

dim3 numBlocks(N / threadsPerBlock.x, N / threadsPerBlock.y);

MatAdd<<<numBlocks, threadsPerBlock>>>(A, B, C);

}

Ниже представлена программа  по расчету произведения матрицы  в случае размещения элементов матрицы в глобальной памяти. Предполагается, что массивы a, b и c в которых размещаются исходные матрицы A, B и матрица произведения C являются линейными.

#include <stdio.h>

#define BLOCK 16 //Устанавливаем размер блока

// Функция умножения двух матриц

__global__ void mulMatr(float* a, float* b, float* c, int n)

{

  //Получаем id текущей нити.

  int i = threadIdx.y+blockIdx.y*blockDim.y;

  int j = threadIdx.x+blockIdx.x*blockDim.x;

  //Расчитываем результат.

     float sum=0.0f;

     for(int p = 0; p < n; p++){

     sum+=  a[i*n + p] * b[p*n + j];}

  c[i*n+j] = sum;

}

__host__ int main()

{

int N;

int M;

float mf=0.0f;

 printf ( "Input N->");

 scanf ( "%d",&N);

 printf ( "Matrix = %dx%d elements\n", N,N );

 M=N*N;

//Выделяем память под вектора

  float* a = new float[M];

  float* b = new float[M];

  float* c = new float[M];

//Инициализируем значения векторов

  for (int i = 0; i < N; i++){

   for (int j = 0; j < N; j++)

    {

      a[i*N+j] = 1.0f*((i+1)+2*(j+1));

      b[i*N+j] = 1.0f/a[i*N+j];

    }

  }

//Указатели на память в видеокарте

  float* deva;

  float* devb;

  float* devc;

//Выделяем память для векторов на видеокарте

  cudaMalloc((void**)&deva, sizeof(float) * M);

  cudaMalloc((void**)&devb, sizeof(float) * M);

  cudaMalloc((void**)&devc, sizeof(float) * M);

// create cuda event handle

cudaEvent_t start, stop;

float gpuTime=0.0f;

cudaEventCreate ( &start );

cudaEventCreate ( &stop );

//Копируем данные в память видеокарты

  cudaMemcpy(deva, a, sizeof(float) * M, cudaMemcpyHostToDevice);

  cudaMemcpy(devb, b, sizeof(float) * M, cudaMemcpyHostToDevice);

//Выполняем вызов функции ядра

dim3 threads = dim3(BLOCK,BLOCK); // Количество нитей в блоке

dim3 blocks  = dim3(N/BLOCK,N/BLOCK); // Количество блоков в grid-е

cudaEventRecord(start, 0); // привязываем событие к началу выполнения ядра

mulMatr<<<blocks, threads>>>(deva, devb, devc,N);

cudaEventRecord(stop, 0); // привязываем событие к концу выполнения ядра

//Получаем результат расчета

 cudaMemcpy(c, devc, sizeof(float) * M, cudaMemcpyDeviceToHost);

cudaEventSynchronize(stop); //Дожидаемся выполнение ядра, синхронизируя

// по событию stop

cudaEventElapsedTime (&gpuTime,start,stop);//Запрашиваем время между start,

// stop

  //Результаты расчета

 gpuTime=gpuTime/1000;

  mf=((2.0*N-1)*N*N)/(gpuTime*1000000.0);

 printf("time=%.4fsec\nspeed=%0.2fMFlops\n",gpuTime,mf);

 printf("i=%d\tj=%d\tC=%.5f\n",N/256,N/128,c[(N/256)*N+(N/128)]);

 printf("i=%d\tj=%d\tC=%.5f\n",3*N/4,5*N/16,c[(3*N/4)*N+(5*N/16)]);

 printf("i=%d\tj=%d\tC=%.5f\n",N-4,N-2,c[(N-4)*N+(N-2)]);

  // Высвобождаем ресурсы

  cudaEventDestroy(start);

  cudaEventDestroy(stop);

  cudaFree(deva);

  cudaFree(devb);

  cudaFree(devc);

  delete[] a; a = 0;

  delete[] b; b = 0;

  delete[] c; c = 0;

}

         Программа по расчету произведения двух квадратных матриц в разделяемой памяти имеет вид:

#include <stdio.h>

#define BLOCK_SIZE 16

#define AS(i, j) As[i][j]

#define BS(i, j) Bs[i][j]

__global__ void

matrixMul( float* C, float* A, float* B, int wA, int wB)

{

    // Block index

    int bx = blockIdx.x;

    int by = blockIdx.y;

    // Thread index

    int tx = threadIdx.x;

    int ty = threadIdx.y;

    // Index of the first sub-matrix of A processed by the block

    int aBegin = wA * BLOCK_SIZE * by;

    // Index of the last sub-matrix of A processed by the block

    int aEnd   = aBegin + wA - 1;

    // Step size used to iterate through the sub-matrices of A

    int aStep  = BLOCK_SIZE;

    // Index of the first sub-matrix of B processed by the block

    int bBegin = BLOCK_SIZE * bx;

    // Step size used to iterate through the sub-matrices of B

    int bStep  = BLOCK_SIZE * wB;

    // Csub is used to store the element of the block sub-matrix

    // that is computed by the thread

    float Csub = 0;

    // Loop over all the sub-matrices of A and B

    // required to compute the block sub-matrix

    for (int a = aBegin, b = bBegin;

             a <= aEnd;

             a += aStep, b += bStep) {

        // Declaration of the shared memory array As used to

        // store the sub-matrix of A

        __shared__ float As[BLOCK_SIZE][BLOCK_SIZE];

        // Declaration of the shared memory array Bs used to

        // store the sub-matrix of B

        __shared__ float Bs[BLOCK_SIZE][BLOCK_SIZE];

        // Load the matrices from device memory

        // to shared memory; each thread loads

        // one element of each matrix

        AS(ty, tx) = A[a + wA * ty + tx];

        BS(ty, tx) = B[b + wB * ty + tx];

        // Synchronize to make sure the matrices are loaded

        __syncthreads();

        // Multiply the two matrices together;

        // each thread computes one element

        // of the block sub-matrix

        for (int k = 0; k < BLOCK_SIZE; ++k)

            Csub += AS(ty, k) * BS(k, tx);

        // Synchronize to make sure that the preceding

        // computation is done before loading two new

        // sub-matrices of A and B in the next iteration

        __syncthreads();

    }

    // Write the block sub-matrix to device memory;

    // each thread writes one element

    int c = wB * BLOCK_SIZE * by + BLOCK_SIZE * bx;

    C[c + wB * ty + tx] = Csub;

}

__host__ int main()

{

int N;

int M;

float mf=0.0f;

 

 printf ( "Input N->");

 scanf ( "%d",&N);

 printf ( "Matrix = %dx%d elements\n", N,N );

 M=N*N;

//Выделяем память под вектора

  float* a = new float[M];

  float* b = new float[M];

  float* c = new float[M];

//Инициализируем значения векторов

  for (int i = 0; i < N; i++){

   for (int j = 0; j < N; j++)

    {

      a[i*N+j] = 1.0f*((i+1)+2*(j+1));

      b[i*N+j] = 1.0f/a[i*N+j];

    }

  }

//Указатели на память в видеокарте

  float* deva;

  float* devb;

  float* devc;

//Выделяем память для векторов на видеокарте

  cudaMalloc((void**)&deva, sizeof(float) * M);

  cudaMalloc((void**)&devb, sizeof(float) * M);

  cudaMalloc((void**)&devc, sizeof(float) * M);

// create cuda event handle

cudaEvent_t start, stop;

float gpuTime=0.0f;

cudaEventCreate ( &start );

cudaEventCreate ( &stop );

//Копируем данные в память видеокарты

  cudaMemcpy(deva, a, sizeof(float) * M, cudaMemcpyHostToDevice);

  cudaMemcpy(devb, b, sizeof(float) * M, cudaMemcpyHostToDevice);

//Выполняем вызов функции ядра

dim3 threads = dim3(BLOCK_SIZE,BLOCK_SIZE);

dim3 blocks  = dim3(N/BLOCK_SIZE,N/BLOCK_SIZE);

cudaEventRecord(start, 0); // привязываем событие к началу выполнения ядра

matrixMul<<<blocks, threads>>>(devc, deva, devb,N,N);

cudaEventRecord(stop, 0); // привязываем событие к концу выполнения ядра

//Получаем результат расчета

 cudaMemcpy(c, devc, sizeof(float) * M, cudaMemcpyDeviceToHost);

cudaEventSynchronize(stop); //Дожидаемся выполнение ядра, синхронизируя //по событию stop

cudaEventElapsedTime (&gpuTime,start,stop);//Запрашиваем время между start, //stop

  //Результаты расчета

 gpuTime=gpuTime/1000;

  mf=((2.0*N-1)*N*N)/(gpuTime*1000000.0);

 printf("time=%.4fsec\nspeed=%0.2fMFlops\n",gpuTime,mf);

 printf("i=%d\tj=%d\tC=%.5f\n",N/256,N/128,c[(N/256)*N+(N/128)]);

 printf("i=%d\tj=%d\tC=%.5f\n",3*N/4,5*N/16,c[(3*N/4)*N+(5*N/16)]);

 printf("i=%d\tj=%d\tC=%.5f\n",N-4,N-2,c[(N-4)*N+(N-2)]);

  // Высвобождаем ресурсы

  cudaEventDestroy(start);

  cudaEventDestroy(stop);

  cudaFree(deva);

  cudaFree(devb);

  cudaFree(devc);

  delete[] a; a = 0;

  delete[] b; b = 0;

  delete[] c; c = 0;

}

Здесь принимается, что из-за малого размера разделяемой памяти, исходные квадратные матрицы разбиваются на блочные матрицы размером 16x16 элементов  (BLOCK_SIZE 16), как в работе [1] .  Текст представленной программы, которая выполняется на GPU (__global__  void  matrixMul), соответствует тексту программы-примера, которая поставляется с сайта http://developer.nvidia.com/cuda-toolkit-32-downloads   под именем matrixMul_kernel.cu.

         Ниже представлена программа умножения матриц с использованием библиотеки CUBLAS.

#include <stdio.h>

#include <cublas.h>

int main ( int argc, char** argv )

{

float time_seconds=0.0f;

float mf=0.0f;

int N;

cudaEvent_t start, stop;

cudaEventCreate ( &start );

cudaEventCreate ( &stop );

printf ( "Input N->");

scanf ( "%d",&N);

printf ( "Matrix = %dx%d elements\n", N,N );

int M=N*N;

float *d_A, *d_B, *d_C;

float* A = new float[M];

float* B = new float[M];

float* C = new float[M];

  for (int j = 0; j < N; j++){

   for (int i = 0; i < N; i++)

    {

      A[i+j*N] = 1.0f*((i+1)+2*(j+1));

      B[i+j*N] = 1.0f/A[i+j*N];

    }

  }

cublasInit();

cublasAlloc ( N * N, sizeof(float), (void**)&d_A);

cublasAlloc ( N * N, sizeof(float), (void**)&d_B);

cublasAlloc ( N * N, sizeof(float), (void**)&d_C);

cublasSetMatrix ( N, N, sizeof(float), (void *) A, N, (void *) d_A, N);

cublasSetMatrix ( N, N, sizeof(float), (void *) B, N, (void *) d_B, N);

cudaEventRecord(start, 0);

cublasSgemm ( 'n', 'n', N, N, N, 1.0f, d_A, N, d_B, N, 0.0f, d_C, N );

cudaEventRecord(stop, 0);

cudaEventSynchronize(stop);

cudaEventElapsedTime (&time_seconds,start,stop);

cublasGetMatrix ( N, N, sizeof(float), (void *) d_C, N, (void *) C, N );

cublasFree (d_A);

cublasFree (d_B);

cublasFree (d_C);

cublasShutdown();

time_seconds=time_seconds/1000;

mf=((2.0*N-1)*N*N)/(time_seconds*1000000.0);

printf("time=%.4fsec\nspeed=%0.2fMFlops\n",time_seconds,mf);

printf("i=%d\tj=%d\tC=%.5f\n",N/256,N/128,C[(N/256)+(N/128)*N]);

printf("i=%d\tj=%d\tC=%.5f\n",3*N/4,5*N/16,C[(3*N/4)+(5*N/16)*N]);

printf("i=%d\tj=%d\tC=%.5f\n",N-4,N-2,C[(N-4)+(N-2)*N]);

}

Результат расчета для одинарной точности NVIDIA (float):

alex@alex-desktop:~/cuda/cuda_matr$ ./mulmatr_blas

Input N->5120

Matrix = 5120x5120 elements

time=0.3254sec

speed=824912.06MFlops

i=20    j=40    C=9647.39844

i=3840  j=1600  C=7792.25586

i=5116  j=5118  C=4011.71313

Результат расчета для одинарной точности на процессоре AMD Phenom II X6 1090T:

alex@alex-desktop:~/fortran/SCALAPACK$ /usr/local/ch_shmem/bin/mpirun -np 6 -machinefile ../machines lab_sca_s

Input N and NB:

5120 128

C(    21,    41)=   .96474023E+04

C(  3841,  1601)=   .77922549E+04

C(  5117,  5119)=   .40117131E+04

 Time =       2.8354 sec.

 Mflops =     94663.1914

Для сопоставления представлен расчет, выполненный на AMD Phenom II X6 1090T с двойной точностью:

alex@alex-desktop:~/fortran/SCALAPACK$ /usr/local/ch_shmem/bin/mpirun -np 6 -machinefile ../machines lab_sca

  Input N and NB:

5120 128

 The following parameter values will be used:

 N= 5120  NB= 128

 nprow= 2  npcol= 3

 Matrix C...

C(    21,    41)=       .964739510929363314E+04

C(  3841,  1601)=       .779225572072732757E+04

C(  5117,  5119)=       .401171136522860297E+04

 times for array with leading dimension of2560

  Time calculation:       5.6869 sec. Mflops =     47197.9957

Видно, что по точности расчетов предпочтение можно отдать технологии CUDA:

NVIDIA(C=9647.39844 – по сравнению с .964739510929363314E+04),

AMD Phenom II X6 (.96474023E+04 – по сравнению с .964739510929363314E+04).

         В таблице 1 сведены результаты расчетов по представленным выше программам для процессора и видеокарты при одинарной точности вычислений. Для процессора даются два варианта: расчет выполняется одним ядром и 6-ю ядрами. Представлена величина ускорения, которая нигде не достигает шестикратной величины. Результаты представлены в виде дроби: числитель – время счета в секундах, знаменатель – производительность в Гигафлопсах в секунду. Производительность определялась как отношение числа операций с плавающей точкой при матричном произведении к  затраченному времени. Число операций в программах определяется по выражению: op=N*N(2*N-1), где N-число строк или столбцов квадратной матрицы.

Таблица 1

Матрица

NxN

Процессор AMD Phenom

Видеокарта GeForce GTX 480

1 ядро

6 ядер

Уск.

Глоб.память

Разд.память

CuBLAS

1024x1024

0.110/

19.50

0.032/

67.07

3.44

0.031/

70.37

0.009/

239.24

0.003/

649.67

2048x2048

0.842/

20.39

0.218/

78.66

3.86

0.243/

70.80

0.071/

240.64

0.022/

769.96

3072x3072

2.806/

20.66

0.624/

92.98

4.50

0.818/

70.89

0.240/

241.74

0.070/

833.12

4096x4096

6.574/

20.90

1.462/

94.02

4.50

2.176/

63.16

0.570/

241.05

0.171/

805.17

5120x5120

12.50/

21.47

2.842/

94.44

4.40

4.148/

64.71

1.110/

241.74

0.325/

824.91

6144x6144

21.60/

21.48

4.635/

100.07

4.66

7.186/

64.55

1.911/

242.67

0.551/

841.45

 

Выводы

1.Производительность видеокарты GeForce GTX 480 с использованием библиотеки CuBLAS (841.45Гфлопс/с) в 8.4 раза выше максимальной производительности процессора AMD Phenom (100.07Гфлопс/с) с использованием библиотек ScaLAPACK и ATLAS.

2.Несмотря на использование библиотеки CuBLAS не удается достичь пиковой производительности видеокарты GeForce GTX 480, которая согласно источнику [8] равна 1344.9Гфлопс/с.

3.Производительность видеокарты GeForce GTX 480 при написании программ традиционным способом без учета тонкостей архитектуры ее вычислительной системы приводит к существенно меньшим значениям производительности (70.80 и 242.67Гфлопс/с) по сравнению с использованием  CuBLAS. Однако если не использовать библиотеки ScaLAPACK и ATLAS для процессора, а распределять матрицы по его ядрам традиционным способом, как это описано в работе [9], наибольшая производительность  процессора не будет превышать 4-х Гфлопс/с, что больше чем в десять раз меньше самой низкой производительности видеокарты при использовании глобальной памяти.

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

5.Результаты расчетов для многоядерного процессора показали, что не удается достичь быстродействия, кратное числу ядер (максимальное ускорение при расчете 6-ю ядрами по сравнению с 1-м ядром не превышает величины 4.66)

 

Литература

1. Боресков А.В., Харламов А.А. Основы работы с технологией CUDA. - М.: ДМК Пресс, 2010. -232 с.: ил.

2. Антонов А.С. Параллельное программирование с использованием технологии MPI: Учебное пособие. – М.: Изд-во МГУ, 2004. – 71 с.

3.The ScaLAPACK Project. http://www.netlib.org

4.Automatically Tuned Linear Algebra Software (ATLAS). http://math-atlas.sourceforge.net/

5. Установка CUDA Toolkit и драйвера NVIDIA для разработчиков. http://forum.ubuntu.ru/index.php?topic=114802.0

6. Мясищев А.А. Достижение наибольшей производительности перемножения матриц на системах с многоядерными процессорами. http://ism.tup.km.ua

7. Мясищев А.А. Оптимизация матричного умножения за счет использования библиотек ScaLAPACK для систем с многоядерными процессорами. Materialy V Miedzynarodowej naukowi-praktycznej konferencji "Nauka i inowacja-2009" Volume 12. Nowoczesne informacyjne technologie. Matematyka.: Przemysl. Nauka i studia - 11-30

8. CUDA. Материал из Википедии — свободной энциклопедии http://ru.wikipedia.org/wiki/CUDA

9.Мясщев А.А. Об эффективности использования 4-х и 6-и ядерных процессоров  Intel и AMD для задач матричного умножения. Материали за 7-а международна научна практична конференция – 2011. Том 28. София, «Бял ГРАД-БГ» ООД -96 стр.