2

更新: タイトルは誤解を招くものです。もともと、block以下のコードでループを展開することでエラーを消すことができました。今では、単純なコードの変更でさえ、それを消すことができます. (以下のコード例を参照してください)。

バックグラウンド：

12x12 行列のコレスキー分解の CUDA C カーネル実装は、非常に大きな CUDA カーネル (280 コード行、多数のループ) になります。

セットアップを減らしてエラーを再現しました（以下のコード）。NVCC (CUDA 4.2) で呼び出される

nvcc -arch sm_20 -o main main.cu

Linux 上の Fermi アーキテクチャで実行:

kernel call: unspecified launch failure

カーネル本体には、条件付きプリプロセッサ ブロック#if 1、#else、が含まれています#endif。これを挿入して、動作中のバージョンと動作していないバージョンを簡単に切り替えました。最初の代替案をコンパイルすると、unspecified launch failure. 一方、2 番目の選択肢は正常に動作します。

注意が必要なのは、実際に実行されるコードはどちらの場合も同じでなければならないということです。(hasOrderedRep は true です!!)

#if 1ステートメントをそのままにしておいても、エラーを消すことができます。そのため、blockループを展開する必要があります。これがタイトルの由来です。

#include "cuda.h"
#include "stdio.h"
#include <iostream>
#include <string>

using namespace std;

/////////////////////////////////////
//
// First some basic types I need
// Implementation of a templated
// scalar and complex number

template<class T> class RScalar
{
public:
  __device__  RScalar() {}

  __device__  ~RScalar() {}

  template<class T1>  __device__   
  RScalar(const RScalar<T1>& rhs) : F(rhs.elem()) {}

  template<class T1>  __device__
  RScalar(const T1& rhs) : F(rhs) {}

  template<class T1> __device__ inline
  RScalar& operator=(const RScalar<T1>& rhs)  {
    elem() = rhs.elem();
    return *this;
  }

public:
  __device__ T& elem() {return F;}
  __device__ const T& elem() const {return F;}

private:
  T F;
};



template<class T> class RComplex
{
public:
  __device__  RComplex() {}

  __device__  ~RComplex() {}

  template<class T1, class T2>  __device__
  RComplex(const RScalar<T1>& _re, const RScalar<T2>& _im): 
    re(_re.elem()), im(_im.elem()) {}

  template<class T1, class T2>  __device__
  RComplex(const T1& _re, const T2& _im): re(_re), im(_im) {}

  template<class T1>  __device__
  RComplex(const T1& _re): re(_re), im() {}

  template<class T1>
  __device__ inline
  RComplex& operator*=(const RScalar<T1>& rhs) 
    {
      real() *= rhs.elem();
      imag() *= rhs.elem();
      return *this;
    }

  template<class T1>
  __device__ inline
  RComplex& operator-=(const RComplex<T1>& rhs) 
    {
      real() -= rhs.real();
      imag() -= rhs.imag();
      return *this;
    }



  template<class T1>  __device__ inline
  RComplex& operator/=(const RComplex<T1>& rhs) 
    {
      RComplex<T> d;
      d = *this / rhs;
      real() = d.real();
      imag() = d.imag();
      return *this;
    }


public:
  __device__ T& real() {return re;}
  __device__ const T& real() const {return re;}

  __device__ T& imag() {return im;}
  __device__ const T& imag() const {return im;}

private:
  T re;
  T im;
};


template<class T> __device__ RComplex<T>
operator*(const RComplex<T>& __restrict__ l, 
      const RComplex<T>& __restrict__ r) 
{
  return RComplex<T>(l.real()*r.real() - l.imag()*r.imag(),
             l.real()*r.imag() + l.imag()*r.real());
}


template<class T> __device__ RComplex<T>
operator/(const RComplex<T>& l, const RComplex<T>& r)
{
  T tmp = T(1.0) / (r.real()*r.real() + r.imag()*r.imag());

  return RComplex<T>((l.real()*r.real() + l.imag()*r.imag()) * tmp,
             (l.imag()*r.real() - l.real()*r.imag()) * tmp);
}


template<class T> __device__ RComplex<T>
operator*(const RComplex<T>& l, const RScalar<T>& r)
{
  return RComplex<T>(l.real()*r.elem(), 
             l.imag()*r.elem());
}

//
//
//////////////////////////////////////////////



#define REALT float
#define Nc 3

struct PrimitiveClovTriang
{
  RScalar<REALT>   diag[2][2*Nc];
  RComplex<REALT>  offd[2][2*Nc*Nc-Nc];
};

__global__ void kernel(bool hasOrderedRep, int * siteTable, 
               PrimitiveClovTriang* tri)
{
  RScalar<REALT> zip=0;
  int N = 2*Nc;
  int site;


  //
  // First if-block results in an error,
  // second, runs fine! Since hasOrderedRep
  // is true, the code blocks should be
  // identical.
  //
#if 1
  if (hasOrderedRep) {
    site = blockDim.x * blockIdx.x + 
      blockDim.x * gridDim.x * blockIdx.y + 
      threadIdx.x;
  } else {
    int idx0 = blockDim.x * blockIdx.x + 
      blockDim.x * gridDim.x * blockIdx.y + 
      threadIdx.x;
    site = ((int*)(siteTable))[idx0];
  }
#else
  site = blockDim.x * blockIdx.x + blockDim.x * gridDim.x * blockIdx.y + threadIdx.x;
#endif

  int site_neg_logdet=0;
  for(int block=0; block < 2; block++) { 
    RScalar<REALT> inv_d[6];
    RComplex<REALT> inv_offd[15];
    RComplex<REALT> v[6];
    RScalar<REALT>  diag_g[6];
    for(int i=0; i < N; i++) { 
      inv_d[i] = tri[site].diag[block][i];
    }
    for(int i=0; i < 15; i++) { 
      inv_offd[i]  =tri[site].offd[block][i];
    }
    for(int j=0; j < N; ++j) { 
      for(int i=0; i < j; i++) { 
    int elem_ji = j*(j-1)/2 + i;
    RComplex<REALT> A_ii = RComplex<REALT>( inv_d[i], zip );
    v[i] = A_ii*RComplex<REALT>(inv_offd[elem_ji].real(),-inv_offd[elem_ji].imag());
      }
      v[j] = RComplex<REALT>(inv_d[j],zip);
      for(int k=0; k < j; k++) { 
    int elem_jk = j*(j-1)/2 + k;
    v[j] -= inv_offd[elem_jk]*v[k];
      }
      inv_d[j].elem() = v[j].real();
      for(int k=j+1; k < N; k++) { 
    int elem_kj = k*(k-1)/2 + j;
    for(int l=0; l < j; l++) { 
      int elem_kl = k*(k-1)/2 + l;
      inv_offd[elem_kj] -= inv_offd[elem_kl] * v[l];
    }
    inv_offd[elem_kj] /= v[j];
      }
    }
    RScalar<REALT> one;
    one.elem() = (REALT)1;
    for(int i=0; i < N; i++) { 
      diag_g[i].elem() = one.elem()/inv_d[i].elem();
      //      ((PScalar<PScalar<RScalar<float> > > *)(args->dev_ptr[ 1 ] ))[site] .elem().elem().elem() += log(fabs(inv_d[i].elem()));
      if( inv_d[i].elem() < 0 ) { 
    site_neg_logdet++;
      }
    }
    RComplex<REALT> sum;
    for(int k = 0; k < N; ++k) {
      for(int i = 0; i < k; ++i) {
    v[i].real()=v[i].imag()=0;
      }
      v[k] = RComplex<REALT>(diag_g[k],zip);
      for(int i = k+1; i < N; ++i) {
    v[i].real()=v[i].imag()=0;
    for(int j = k; j < i; ++j) {
      int elem_ij = i*(i-1)/2+j;    
      v[i] -= inv_offd[elem_ij] *inv_d[j]*v[j];
    }
    v[i] *= diag_g[i];
      }
      for(int i = N-2; (int)i >= (int)k; --i) {
    for(int j = i+1; j < N; ++j) {
      int elem_ji = j*(j-1)/2 + i;
      v[i] -= RComplex<REALT>(inv_offd[elem_ji].real(),-inv_offd[elem_ji].imag()) * v[j];
    }
      }
      inv_d[k].elem() = v[k].real();
      for(int i = k+1; i < N; ++i) {
    int elem_ik = i*(i-1)/2+k;
    inv_offd[elem_ik] = v[i];
      }
    }
    for(int i=0; i < N; i++) { 
      tri[site].diag[block][i] = inv_d[i];
    }
    for(int i=0; i < 15; i++) { 
      tri[site].offd[block][i] = inv_offd[i];
    }
  }
  if( site_neg_logdet != 0 ) { 
  }

}




int main()
{
  int sites=1;

  dim3  blocksPerGrid( 1 , 1 , 1 );
  dim3  threadsPerBlock( sites , 1, 1);

  PrimitiveClovTriang* tri_dev;
  int * siteTable;
  cudaMalloc( (void**)&tri_dev , sizeof(PrimitiveClovTriang) * sites );
  cudaMalloc( (void**)&siteTable , sizeof(int) * sites );
  bool ord=true;

  kernel<<< blocksPerGrid , threadsPerBlock , 0 >>>( ord , siteTable , tri_dev );

  cudaDeviceSynchronize();
  cudaError_t kernel_call = cudaGetLastError();

  cout << "kernel call: " << string(cudaGetErrorString(kernel_call)) << endl;
  cudaFree(tri_dev);
  cudaFree(siteTable);

  return(0);
}

siteTable には初期化されていないデータが含まれていることに注意してください。未使用なので大丈夫です。エラーを表示するために必要です。

アップデート：

CUDA 4.0 がインストールされている別のマシンで試してみました。エラーは表示されません (同じ Fermi カード モデル)。CUDA 4.2 の NVCC バグかもしれません。彼らは CUDA 4.1 から LLVM に切り替えたので、これはおそらくバグです。

4

1 に答える 1

1

この問題は、LLVM ベースの CUDA ツールチェーンの初期リリースにおける微妙なコンパイラのバグによって引き起こされたようです。このエラーは、ツールチェーンの最新バージョン (CUDA 6.x および 7.x) では再現されません。

[この回答は、CUDA タグの回答リストから質問を削除するために、コメントと編集からコミュニティ wiki エントリとして追加されました]

于 2016-03-23T09:55:25.623 に答える