Implemeting multidimensional stencil in Cuda C/C++












0












$begingroup$


I've constructed a system that contains multiple 2D lattice and a 1D cable all in a single vector of type double. Each 2D lattice is coupled to a site along the cable. For instance, there are 800 nodes on the cable which stores the state of each 800 50-by-50 lattice that it is coupled to. The size of the vector is then 800*50*50.
Using CUDA C, I implement a generalized "stencil" kernel function that computes the average of n - nearest neighbors of a given index. In this example, the stencil is generalized to be used for both the 1D cable and the 2D lattice. The code runs fine, except it takes longer than the code without CUDA! Obviously, I am not implementing it correctly. The following snippet demonstrates CUDA protocol. The OFFSETS values take the iterator to the start and end points of one of the 2D lattice that is being operated on.



first = v.begin() + OFFSET0;    
last = v.begin() + OFFSET1;
std::copy(first,last,host_vector);
cudaMemcpy(d_in,host_vector,size,cudaMemcpyHostToDevice);
_1Dstencil<<<M,THREADS_PER_BLOCK>>>(d_in,d_out,nx,ny,nz,X_,Y_,Z_);
cudaMemcpy(host_vector,d_out,size,cudaMemcpyDeviceToHost);


The following is my implementation of the generalized kernel stencil



__global__ void _1Dstencil(double* d_in,double* d_out,int nx,int ny,int nz,int X_,int Y_,int Z_){
int i = threadIdx.x+ blockIdx.x* blockDim.x ;
int i_ = (i/(Z_*Y_))%X_ ,
j_ = (i/Z_)%Y_,
k_ = i%Z_ ;
double sum = 0;
//stencil around i
for(int nnx = 2*nx+1; nnx> 0; nnx--){//X
int start = i_-(nnx-nx-1);
int ii = cp(start,X_,false);
for(int nny = 2*ny+1; nny> 0; nny--){//Y
start = j_-(nny-ny-1);
int jj = cp(start,Y_,false);
for(int nnz = 2*nz+1; nnz> 0; nnz--){//Z
start = k_-(nnz-nz-1);
int kk = cp(start,Z_,false);
sum += *(d_in + _ND_1D(ii, jj, kk, X_, Y_, Z_));
}
}
}
d_out[i] = sum/((2*nx+1)*(2*ny+1)*(2*nz+1));
}


The functions cp(...) and _ND_1D(...) are declared as device functions which return the correct position based on flux boundaries and convert N-Dimensional coordinates to 1D coordinates, respectively. Could someone help me clarify the correct way of distributing the task on the deice? Thanks in advance!










share|improve this question







New contributor




Lrom is a new contributor to this site. Take care in asking for clarification, commenting, and answering.
Check out our Code of Conduct.







$endgroup$

















    0












    $begingroup$


    I've constructed a system that contains multiple 2D lattice and a 1D cable all in a single vector of type double. Each 2D lattice is coupled to a site along the cable. For instance, there are 800 nodes on the cable which stores the state of each 800 50-by-50 lattice that it is coupled to. The size of the vector is then 800*50*50.
    Using CUDA C, I implement a generalized "stencil" kernel function that computes the average of n - nearest neighbors of a given index. In this example, the stencil is generalized to be used for both the 1D cable and the 2D lattice. The code runs fine, except it takes longer than the code without CUDA! Obviously, I am not implementing it correctly. The following snippet demonstrates CUDA protocol. The OFFSETS values take the iterator to the start and end points of one of the 2D lattice that is being operated on.



    first = v.begin() + OFFSET0;    
    last = v.begin() + OFFSET1;
    std::copy(first,last,host_vector);
    cudaMemcpy(d_in,host_vector,size,cudaMemcpyHostToDevice);
    _1Dstencil<<<M,THREADS_PER_BLOCK>>>(d_in,d_out,nx,ny,nz,X_,Y_,Z_);
    cudaMemcpy(host_vector,d_out,size,cudaMemcpyDeviceToHost);


    The following is my implementation of the generalized kernel stencil



    __global__ void _1Dstencil(double* d_in,double* d_out,int nx,int ny,int nz,int X_,int Y_,int Z_){
    int i = threadIdx.x+ blockIdx.x* blockDim.x ;
    int i_ = (i/(Z_*Y_))%X_ ,
    j_ = (i/Z_)%Y_,
    k_ = i%Z_ ;
    double sum = 0;
    //stencil around i
    for(int nnx = 2*nx+1; nnx> 0; nnx--){//X
    int start = i_-(nnx-nx-1);
    int ii = cp(start,X_,false);
    for(int nny = 2*ny+1; nny> 0; nny--){//Y
    start = j_-(nny-ny-1);
    int jj = cp(start,Y_,false);
    for(int nnz = 2*nz+1; nnz> 0; nnz--){//Z
    start = k_-(nnz-nz-1);
    int kk = cp(start,Z_,false);
    sum += *(d_in + _ND_1D(ii, jj, kk, X_, Y_, Z_));
    }
    }
    }
    d_out[i] = sum/((2*nx+1)*(2*ny+1)*(2*nz+1));
    }


    The functions cp(...) and _ND_1D(...) are declared as device functions which return the correct position based on flux boundaries and convert N-Dimensional coordinates to 1D coordinates, respectively. Could someone help me clarify the correct way of distributing the task on the deice? Thanks in advance!










    share|improve this question







    New contributor




    Lrom is a new contributor to this site. Take care in asking for clarification, commenting, and answering.
    Check out our Code of Conduct.







    $endgroup$















      0












      0








      0





      $begingroup$


      I've constructed a system that contains multiple 2D lattice and a 1D cable all in a single vector of type double. Each 2D lattice is coupled to a site along the cable. For instance, there are 800 nodes on the cable which stores the state of each 800 50-by-50 lattice that it is coupled to. The size of the vector is then 800*50*50.
      Using CUDA C, I implement a generalized "stencil" kernel function that computes the average of n - nearest neighbors of a given index. In this example, the stencil is generalized to be used for both the 1D cable and the 2D lattice. The code runs fine, except it takes longer than the code without CUDA! Obviously, I am not implementing it correctly. The following snippet demonstrates CUDA protocol. The OFFSETS values take the iterator to the start and end points of one of the 2D lattice that is being operated on.



      first = v.begin() + OFFSET0;    
      last = v.begin() + OFFSET1;
      std::copy(first,last,host_vector);
      cudaMemcpy(d_in,host_vector,size,cudaMemcpyHostToDevice);
      _1Dstencil<<<M,THREADS_PER_BLOCK>>>(d_in,d_out,nx,ny,nz,X_,Y_,Z_);
      cudaMemcpy(host_vector,d_out,size,cudaMemcpyDeviceToHost);


      The following is my implementation of the generalized kernel stencil



      __global__ void _1Dstencil(double* d_in,double* d_out,int nx,int ny,int nz,int X_,int Y_,int Z_){
      int i = threadIdx.x+ blockIdx.x* blockDim.x ;
      int i_ = (i/(Z_*Y_))%X_ ,
      j_ = (i/Z_)%Y_,
      k_ = i%Z_ ;
      double sum = 0;
      //stencil around i
      for(int nnx = 2*nx+1; nnx> 0; nnx--){//X
      int start = i_-(nnx-nx-1);
      int ii = cp(start,X_,false);
      for(int nny = 2*ny+1; nny> 0; nny--){//Y
      start = j_-(nny-ny-1);
      int jj = cp(start,Y_,false);
      for(int nnz = 2*nz+1; nnz> 0; nnz--){//Z
      start = k_-(nnz-nz-1);
      int kk = cp(start,Z_,false);
      sum += *(d_in + _ND_1D(ii, jj, kk, X_, Y_, Z_));
      }
      }
      }
      d_out[i] = sum/((2*nx+1)*(2*ny+1)*(2*nz+1));
      }


      The functions cp(...) and _ND_1D(...) are declared as device functions which return the correct position based on flux boundaries and convert N-Dimensional coordinates to 1D coordinates, respectively. Could someone help me clarify the correct way of distributing the task on the deice? Thanks in advance!










      share|improve this question







      New contributor




      Lrom is a new contributor to this site. Take care in asking for clarification, commenting, and answering.
      Check out our Code of Conduct.







      $endgroup$




      I've constructed a system that contains multiple 2D lattice and a 1D cable all in a single vector of type double. Each 2D lattice is coupled to a site along the cable. For instance, there are 800 nodes on the cable which stores the state of each 800 50-by-50 lattice that it is coupled to. The size of the vector is then 800*50*50.
      Using CUDA C, I implement a generalized "stencil" kernel function that computes the average of n - nearest neighbors of a given index. In this example, the stencil is generalized to be used for both the 1D cable and the 2D lattice. The code runs fine, except it takes longer than the code without CUDA! Obviously, I am not implementing it correctly. The following snippet demonstrates CUDA protocol. The OFFSETS values take the iterator to the start and end points of one of the 2D lattice that is being operated on.



      first = v.begin() + OFFSET0;    
      last = v.begin() + OFFSET1;
      std::copy(first,last,host_vector);
      cudaMemcpy(d_in,host_vector,size,cudaMemcpyHostToDevice);
      _1Dstencil<<<M,THREADS_PER_BLOCK>>>(d_in,d_out,nx,ny,nz,X_,Y_,Z_);
      cudaMemcpy(host_vector,d_out,size,cudaMemcpyDeviceToHost);


      The following is my implementation of the generalized kernel stencil



      __global__ void _1Dstencil(double* d_in,double* d_out,int nx,int ny,int nz,int X_,int Y_,int Z_){
      int i = threadIdx.x+ blockIdx.x* blockDim.x ;
      int i_ = (i/(Z_*Y_))%X_ ,
      j_ = (i/Z_)%Y_,
      k_ = i%Z_ ;
      double sum = 0;
      //stencil around i
      for(int nnx = 2*nx+1; nnx> 0; nnx--){//X
      int start = i_-(nnx-nx-1);
      int ii = cp(start,X_,false);
      for(int nny = 2*ny+1; nny> 0; nny--){//Y
      start = j_-(nny-ny-1);
      int jj = cp(start,Y_,false);
      for(int nnz = 2*nz+1; nnz> 0; nnz--){//Z
      start = k_-(nnz-nz-1);
      int kk = cp(start,Z_,false);
      sum += *(d_in + _ND_1D(ii, jj, kk, X_, Y_, Z_));
      }
      }
      }
      d_out[i] = sum/((2*nx+1)*(2*ny+1)*(2*nz+1));
      }


      The functions cp(...) and _ND_1D(...) are declared as device functions which return the correct position based on flux boundaries and convert N-Dimensional coordinates to 1D coordinates, respectively. Could someone help me clarify the correct way of distributing the task on the deice? Thanks in advance!







      c++ cuda






      share|improve this question







      New contributor




      Lrom is a new contributor to this site. Take care in asking for clarification, commenting, and answering.
      Check out our Code of Conduct.











      share|improve this question







      New contributor




      Lrom is a new contributor to this site. Take care in asking for clarification, commenting, and answering.
      Check out our Code of Conduct.









      share|improve this question




      share|improve this question






      New contributor




      Lrom is a new contributor to this site. Take care in asking for clarification, commenting, and answering.
      Check out our Code of Conduct.









      asked 11 mins ago









      LromLrom

      11




      11




      New contributor




      Lrom is a new contributor to this site. Take care in asking for clarification, commenting, and answering.
      Check out our Code of Conduct.





      New contributor





      Lrom is a new contributor to this site. Take care in asking for clarification, commenting, and answering.
      Check out our Code of Conduct.






      Lrom is a new contributor to this site. Take care in asking for clarification, commenting, and answering.
      Check out our Code of Conduct.






















          0






          active

          oldest

          votes











          Your Answer





          StackExchange.ifUsing("editor", function () {
          return StackExchange.using("mathjaxEditing", function () {
          StackExchange.MarkdownEditor.creationCallbacks.add(function (editor, postfix) {
          StackExchange.mathjaxEditing.prepareWmdForMathJax(editor, postfix, [["\$", "\$"]]);
          });
          });
          }, "mathjax-editing");

          StackExchange.ifUsing("editor", function () {
          StackExchange.using("externalEditor", function () {
          StackExchange.using("snippets", function () {
          StackExchange.snippets.init();
          });
          });
          }, "code-snippets");

          StackExchange.ready(function() {
          var channelOptions = {
          tags: "".split(" "),
          id: "196"
          };
          initTagRenderer("".split(" "), "".split(" "), channelOptions);

          StackExchange.using("externalEditor", function() {
          // Have to fire editor after snippets, if snippets enabled
          if (StackExchange.settings.snippets.snippetsEnabled) {
          StackExchange.using("snippets", function() {
          createEditor();
          });
          }
          else {
          createEditor();
          }
          });

          function createEditor() {
          StackExchange.prepareEditor({
          heartbeatType: 'answer',
          autoActivateHeartbeat: false,
          convertImagesToLinks: false,
          noModals: true,
          showLowRepImageUploadWarning: true,
          reputationToPostImages: null,
          bindNavPrevention: true,
          postfix: "",
          imageUploader: {
          brandingHtml: "Powered by u003ca class="icon-imgur-white" href="https://imgur.com/"u003eu003c/au003e",
          contentPolicyHtml: "User contributions licensed under u003ca href="https://creativecommons.org/licenses/by-sa/3.0/"u003ecc by-sa 3.0 with attribution requiredu003c/au003e u003ca href="https://stackoverflow.com/legal/content-policy"u003e(content policy)u003c/au003e",
          allowUrls: true
          },
          onDemand: true,
          discardSelector: ".discard-answer"
          ,immediatelyShowMarkdownHelp:true
          });


          }
          });






          Lrom is a new contributor. Be nice, and check out our Code of Conduct.










          draft saved

          draft discarded


















          StackExchange.ready(
          function () {
          StackExchange.openid.initPostLogin('.new-post-login', 'https%3a%2f%2fcodereview.stackexchange.com%2fquestions%2f214005%2fimplemeting-multidimensional-stencil-in-cuda-c-c%23new-answer', 'question_page');
          }
          );

          Post as a guest















          Required, but never shown

























          0






          active

          oldest

          votes








          0






          active

          oldest

          votes









          active

          oldest

          votes






          active

          oldest

          votes








          Lrom is a new contributor. Be nice, and check out our Code of Conduct.










          draft saved

          draft discarded


















          Lrom is a new contributor. Be nice, and check out our Code of Conduct.













          Lrom is a new contributor. Be nice, and check out our Code of Conduct.












          Lrom is a new contributor. Be nice, and check out our Code of Conduct.
















          Thanks for contributing an answer to Code Review Stack Exchange!


          • Please be sure to answer the question. Provide details and share your research!

          But avoid



          • Asking for help, clarification, or responding to other answers.

          • Making statements based on opinion; back them up with references or personal experience.


          Use MathJax to format equations. MathJax reference.


          To learn more, see our tips on writing great answers.




          draft saved


          draft discarded














          StackExchange.ready(
          function () {
          StackExchange.openid.initPostLogin('.new-post-login', 'https%3a%2f%2fcodereview.stackexchange.com%2fquestions%2f214005%2fimplemeting-multidimensional-stencil-in-cuda-c-c%23new-answer', 'question_page');
          }
          );

          Post as a guest















          Required, but never shown





















































          Required, but never shown














          Required, but never shown












          Required, but never shown







          Required, but never shown

































          Required, but never shown














          Required, but never shown












          Required, but never shown







          Required, but never shown







          Popular posts from this blog

          Сан-Квентин

          8-я гвардейская общевойсковая армия

          Алькесар