<div dir="ltr">Hi,<br><div><br></div><div>I was wondering what the official status of 64-bit integer support in the PETSc GPU backend is (specifically CUDA). This question comes from the result of benchmarking some PETSc code and looking at some sources. In particular, I found that PETSc's call to cuSPARSE SpMV seems to always be using the 32-bit integer call, even if I compile PETSc with `--with-64-bit-indices`. After digging around more, I see that PETSc always only creates 32-bit cuSPARSE matrices as well: <a href="https://gitlab.com/petsc/petsc/-/blob/v3.19.4/src/mat/impls/aij/seq/seqcusparse/aijcusparse.cu?ref_type=tags#L2501" target="_blank">https://gitlab.com/petsc/petsc/-/blob/v3.19.4/src/mat/impls/aij/seq/seqcusparse/aijcusparse.cu?ref_type=tags#L2501</a>. I was looking around for a switch somewhere to 64 bit integers inside this code, but everything seems to be pretty hardcoded with `<code>THRUSTINTARRAY32</code>`.</div><div><br></div><div>As expected, this all works when the range of coordinates in each sparse matrix partition is less than INT_MAX, but PETSc GPU code breaks in different ways (calling cuBLAS and cuSPARSE) when trying a (synthetic) problem that needs 64 bit integers:</div><div><br></div><div>```</div><div>#include "petscmat.h"<br>#include "petscvec.h"<br>#include "petsc.h"<br><br>int main(int argc, char** argv) {<br>  PetscInt ierr;<br>  PetscInitialize(&argc, &argv, (char *)0, "GPU bug");<br><br>  PetscInt numRows = 1;<br>  PetscInt numCols = PetscInt(INT_MAX) * 2;<br><br>  Mat A;<br>  PetscInt rowStart, rowEnd;<br>  ierr = MatCreate(PETSC_COMM_WORLD, &A); CHKERRQ(ierr);<br>  MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, numRows, numCols);<br>  MatSetType(A, MATMPIAIJ);<br>  MatSetFromOptions(A);<br><br>  MatSetValue(A, 0, 0, 1.0, INSERT_VALUES);<br>  MatSetValue(A, 0, numCols - 1, 1.0, INSERT_VALUES);<br>  MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);<br>  MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);<br><br>  Vec b;<br>  ierr = VecCreate(PETSC_COMM_WORLD, &b); CHKERRQ(ierr);<br>  VecSetSizes(b, PETSC_DECIDE, numCols);<br>  VecSetFromOptions(b);<br>  VecSet(b, 0.0);<br>  VecSetValue(b, 0, 42.0, INSERT_VALUES);<br>  VecSetValue(b, numCols - 1, 58.0, INSERT_VALUES);<br>  VecAssemblyBegin(b);<br>  VecAssemblyEnd(b);<br><br>  Vec x;<br>  ierr = VecCreate(PETSC_COMM_WORLD, &x); CHKERRQ(ierr);<br>  VecSetSizes(x, PETSC_DECIDE, numRows);<br>  VecSetFromOptions(x);<br>  VecSet(x, 0.0);<br><br>  MatMult(A, b, x);<br>  PetscScalar result;<br>  VecSum(x, &result);<br>  PetscPrintf(PETSC_COMM_WORLD, "Result of mult: %f\n", result);<br>  PetscFinalize();<br>}<br></div><div>```</div><div><br></div><div>When this program is run on CPUs, it outputs 100.0, as expected.</div><div><br></div><div>When run on a single GPU with `-vec_type cuda -mat_type aijcusparse -use_gpu_aware_mpi 0` it fails with</div><div>```</div><div>





<p style="margin:0px;font-variant-numeric:normal;font-variant-east-asian:normal;font-variant-alternates:normal;font-kerning:auto;font-feature-settings:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue"">[0]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------</p>
<p style="margin:0px;font-variant-numeric:normal;font-variant-east-asian:normal;font-variant-alternates:normal;font-kerning:auto;font-feature-settings:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue"">[0]PETSC ERROR: Argument out of range</p>
<p style="margin:0px;font-variant-numeric:normal;font-variant-east-asian:normal;font-variant-alternates:normal;font-kerning:auto;font-feature-settings:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue"">[0]PETSC ERROR: 4294967294 is too big for cuBLAS, which may be restricted to 32-bit integers</p>
<p style="margin:0px;font-variant-numeric:normal;font-variant-east-asian:normal;font-variant-alternates:normal;font-kerning:auto;font-feature-settings:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue"">[0]PETSC ERROR: See <a href="https://petsc.org/release/faq/" target="_blank"><span style="color:rgb(220,161,13)">https://petsc.org/release/faq/</span></a> for trouble shooting.</p>
<p style="margin:0px;font-variant-numeric:normal;font-variant-east-asian:normal;font-variant-alternates:normal;font-kerning:auto;font-feature-settings:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue"">[0]PETSC ERROR: Petsc Release Version 3.19.4, unknown</p>
<p style="margin:0px;font-variant-numeric:normal;font-variant-east-asian:normal;font-variant-alternates:normal;font-kerning:auto;font-feature-settings:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue"">[0]PETSC ERROR: ./gpu-bug on a<span>  </span>named sean-dgx2 by rohany Fri Aug 11 09:34:10 2023</p>
<p style="margin:0px;font-variant-numeric:normal;font-variant-east-asian:normal;font-variant-alternates:normal;font-kerning:auto;font-feature-settings:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue"">[0]PETSC ERROR: Configure options --with-cuda=1 --prefix=/local/home/rohany/petsc/petsc-install/ --with-cuda-dir=/usr/local/cuda-11.7/ CXXFLAGS=-O3 COPTFLAGS=-O3 CXXOPTFLAGS=-O3 FOPTFLAGS=-O3 --download-fblaslapack=1 --with-debugging=0 --with-64-bit-indices</p>
<p style="margin:0px;font-variant-numeric:normal;font-variant-east-asian:normal;font-variant-alternates:normal;font-kerning:auto;font-feature-settings:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue"">[0]PETSC ERROR: #1 checkCupmBlasIntCast() at /local/home/rohany/petsc/include/petsc/private/cupmblasinterface.hpp:435</p>
<p style="margin:0px;font-variant-numeric:normal;font-variant-east-asian:normal;font-variant-alternates:normal;font-kerning:auto;font-feature-settings:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue"">[0]PETSC ERROR: #2 VecAllocateCheck_() at /local/home/rohany/petsc/include/petsc/private/veccupmimpl.h:335</p>
<p style="margin:0px;font-variant-numeric:normal;font-variant-east-asian:normal;font-variant-alternates:normal;font-kerning:auto;font-feature-settings:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue"">[0]PETSC ERROR: #3 VecCUPMAllocateCheck_() at /local/home/rohany/petsc/include/petsc/private/veccupmimpl.h:360</p>
<p style="margin:0px;font-variant-numeric:normal;font-variant-east-asian:normal;font-variant-alternates:normal;font-kerning:auto;font-feature-settings:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue"">[0]PETSC ERROR: #4 DeviceAllocateCheck_() at /local/home/rohany/petsc/include/petsc/private/veccupmimpl.h:389</p>
<p style="margin:0px;font-variant-numeric:normal;font-variant-east-asian:normal;font-variant-alternates:normal;font-kerning:auto;font-feature-settings:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue"">[0]PETSC ERROR: #5 GetArray() at /local/home/rohany/petsc/include/petsc/private/veccupmimpl.h:545</p>
<p style="margin:0px;font-variant-numeric:normal;font-variant-east-asian:normal;font-variant-alternates:normal;font-kerning:auto;font-feature-settings:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue"">[0]PETSC ERROR: #6 VectorArray() at /local/home/rohany/petsc/include/petsc/private/veccupmimpl.h:273</p>
<p style="margin:0px;font-variant-numeric:normal;font-variant-east-asian:normal;font-variant-alternates:normal;font-kerning:auto;font-feature-settings:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue"">--------------------------------------------------------------------------</p>
<p style="margin:0px;font-variant-numeric:normal;font-variant-east-asian:normal;font-variant-alternates:normal;font-kerning:auto;font-feature-settings:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue"">MPI_ABORT was invoked on rank 0 in communicator MPI_COMM_SELF</p>
<p style="margin:0px;font-variant-numeric:normal;font-variant-east-asian:normal;font-variant-alternates:normal;font-kerning:auto;font-feature-settings:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue"">with errorcode 63.</p>
<p style="margin:0px;font-variant-numeric:normal;font-variant-east-asian:normal;font-variant-alternates:normal;font-kerning:auto;font-feature-settings:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue";min-height:15px"><br></p>
<p style="margin:0px;font-variant-numeric:normal;font-variant-east-asian:normal;font-variant-alternates:normal;font-kerning:auto;font-feature-settings:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue"">NOTE: invoking MPI_ABORT causes Open MPI to kill all MPI processes.</p>
<p style="margin:0px;font-variant-numeric:normal;font-variant-east-asian:normal;font-variant-alternates:normal;font-kerning:auto;font-feature-settings:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue"">You may or may not see output from other processes, depending on</p>
<p style="margin:0px;font-variant-numeric:normal;font-variant-east-asian:normal;font-variant-alternates:normal;font-kerning:auto;font-feature-settings:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue"">exactly when Open MPI kills them.</p>
<p style="margin:0px;font-variant-numeric:normal;font-variant-east-asian:normal;font-variant-alternates:normal;font-kerning:auto;font-feature-settings:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue"">--------------------------------------------------------------------------</p><p style="margin:0px;font-variant-numeric:normal;font-variant-east-asian:normal;font-variant-alternates:normal;font-kerning:auto;font-feature-settings:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue"">```</p><p style="margin:0px;font-variant-numeric:normal;font-variant-east-asian:normal;font-variant-alternates:normal;font-kerning:auto;font-feature-settings:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue""><br></p><p style="margin:0px;font-variant-numeric:normal;font-variant-east-asian:normal;font-variant-alternates:normal;font-kerning:auto;font-feature-settings:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue"">and when run with just `-mat_type aijcusparse -use_gpu_aware_mpi 0` it fails with</p><p style="margin:0px;font-variant-numeric:normal;font-variant-east-asian:normal;font-variant-alternates:normal;font-kerning:auto;font-feature-settings:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue"">```</p><p class="gmail-p1" style="margin:0px;font:13px "Helvetica Neue""><span class="gmail-Apple-converted-space"> </span>** On entry to cusparseCreateCsr(): dimension mismatch for CUSPARSE_INDEX_32I, cols (4294967294) + base (0) > INT32_MAX (2147483647)</p><p class="gmail-p2" style="margin:0px;font:13px "Helvetica Neue";min-height:15px"><br></p><p class="gmail-p1" style="margin:0px;font:13px "Helvetica Neue"">[0]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------</p><p class="gmail-p1" style="margin:0px;font:13px "Helvetica Neue"">[0]PETSC ERROR: GPU error</p><p class="gmail-p1" style="margin:0px;font:13px "Helvetica Neue"">[0]PETSC ERROR: cuSPARSE errorcode 3 (CUSPARSE_STATUS_INVALID_VALUE) : invalid value</p><p class="gmail-p1" style="margin:0px;font:13px "Helvetica Neue"">[0]PETSC ERROR: See <a href="https://petsc.org/release/faq/"><span class="gmail-s1" style="color:rgb(220,161,13)">https://petsc.org/release/faq/</span></a> for trouble shooting.</p><p class="gmail-p1" style="margin:0px;font:13px "Helvetica Neue"">[0]PETSC ERROR: Petsc Release Version 3.19.4, unknown</p><p class="gmail-p1" style="margin:0px;font:13px "Helvetica Neue"">[0]PETSC ERROR: ./gpu-bug on a<span class="gmail-Apple-converted-space">  </span>named sean-dgx2 by rohany Fri Aug 11 09:43:07 2023</p><p class="gmail-p1" style="margin:0px;font:13px "Helvetica Neue"">[0]PETSC ERROR: Configure options --with-cuda=1 --prefix=/local/home/rohany/petsc/petsc-install/ --with-cuda-dir=/usr/local/cuda-11.7/ CXXFLAGS=-O3 COPTFLAGS=-O3 CXXOPTFLAGS=-O3 FOPTFLAGS=-O3 --download-fblaslapack=1 --with-debugging=0 --with-64-bit-indices</p><p class="gmail-p1" style="margin:0px;font:13px "Helvetica Neue"">[0]PETSC ERROR: #1 MatSeqAIJCUSPARSECopyToGPU() at /local/home/rohany/petsc/src/mat/impls/aij/seq/seqcusparse/<a href="http://aijcusparse.cu:2503">aijcusparse.cu:2503</a></p><p class="gmail-p1" style="margin:0px;font:13px "Helvetica Neue"">[0]PETSC ERROR: #2 MatMultAddKernel_SeqAIJCUSPARSE() at /local/home/rohany/petsc/src/mat/impls/aij/seq/seqcusparse/<a href="http://aijcusparse.cu:3544">aijcusparse.cu:3544</a></p><p class="gmail-p1" style="margin:0px;font:13px "Helvetica Neue"">[0]PETSC ERROR: #3 MatMult_SeqAIJCUSPARSE() at /local/home/rohany/petsc/src/mat/impls/aij/seq/seqcusparse/<a href="http://aijcusparse.cu:3485">aijcusparse.cu:3485</a></p><p class="gmail-p1" style="margin:0px;font:13px "Helvetica Neue"">[0]PETSC ERROR: #4 MatMult_MPIAIJCUSPARSE() at /local/home/rohany/petsc/src/mat/impls/aij/mpi/mpicusparse/<a href="http://mpiaijcusparse.cu:452">mpiaijcusparse.cu:452</a></p><p style="margin:0px;font-variant-numeric:normal;font-variant-east-asian:normal;font-variant-alternates:normal;font-kerning:auto;font-feature-settings:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue"">


















</p><p class="gmail-p1" style="margin:0px;font-variant-numeric:normal;font-variant-east-asian:normal;font-variant-alternates:normal;font-kerning:auto;font-feature-settings:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue"">[0]PETSC ERROR: #5 MatMult() at /local/home/rohany/petsc/src/mat/interface/matrix.c:2599</p><p style="margin:0px;font-variant-numeric:normal;font-variant-east-asian:normal;font-variant-alternates:normal;font-kerning:auto;font-feature-settings:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue"">```</p><p style="margin:0px;font-variant-numeric:normal;font-variant-east-asian:normal;font-variant-alternates:normal;font-kerning:auto;font-feature-settings:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue""><br></p><p style="margin:0px;font-variant-numeric:normal;font-variant-east-asian:normal;font-variant-alternates:normal;font-kerning:auto;font-feature-settings:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue"">Thanks,</p><p style="margin:0px;font-variant-numeric:normal;font-variant-east-asian:normal;font-variant-alternates:normal;font-kerning:auto;font-feature-settings:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue""><br></p><p style="margin:0px;font-variant-numeric:normal;font-variant-east-asian:normal;font-variant-alternates:normal;font-kerning:auto;font-feature-settings:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue"">Rohan Yadav</p></div><div><br></div></div>