<div dir="ltr">Rohan,<br> You could try the petsc/kokkos backend. I have not tested it, but I guess it should handle 64 bit CUDA index types.<br> I guess the petsc/cuda 32-bit limit came from old CUDA versions where only 32-bit indices were supported such that the original developers hardwired the type to THRUSTINTARRAY32. We try to support generations of cuda toolkits and thus have the current code.<br><div><br></div><div> Anyway, this should be fixed.</div><div><div><div dir="ltr" class="gmail_signature" data-smartmail="gmail_signature"><div dir="ltr">--Junchao Zhang</div></div></div><br></div></div><br><div class="gmail_quote"><div dir="ltr" class="gmail_attr">On Fri, Aug 11, 2023 at 1:07 PM Barry Smith <<a href="mailto:bsmith@petsc.dev">bsmith@petsc.dev</a>> wrote:<br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div><div><br></div> We do not currently have any code for using 64 bit integer sizes on the GPUs. <div><br></div><div> Given the current memory available on GPUs is 64 bit integer support needed? I think even a single vector of length 2^31 will use up most of the GPU's memory? Are the practical, not synthetic, situations that require 64 bit integer support on GPUs immediately? For example, is the vector length of the entire parallel vector across all GPUs limited to 32 bits? </div><div><br></div><div> We will certainly add such support, but it is a question of priorities; there are many things we need to do to improve PETSc GPU support, and they take time. Unless we have practical use cases, 64 bit integer support for integer sizes on the GPU is not at the top of the list. Of course, we would be very happy with a merge request that would provide this support at any time.</div><div><br></div><div> Barry</div><div><br></div><div><br><div><br><blockquote type="cite"><div>On Aug 11, 2023, at 1:23 PM, Rohan Yadav <<a href="mailto:rohany@alumni.cmu.edu" target="_blank">rohany@alumni.cmu.edu</a>> wrote:</div><br><div><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><div 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 --------------------------------------------------------------</div><div 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</div><div 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</div><div 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.</div><div 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</div><div 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</div><div 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</div><div 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</div><div 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</div><div 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</div><div 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</div><div 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</div><div 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</div><div 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"">--------------------------------------------------------------------------</div><div 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</div><div 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.</div><div 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></div><div 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.</div><div 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</div><div 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.</div><div 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"">--------------------------------------------------------------------------</div><div 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"">```</div><div 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></div><div 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</div><div 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"">```</div><div style="margin:0px;font-style:normal;font-variant:normal;font-weight:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue";font-kerning:auto;font-feature-settings:normal"><span> </span>** On entry to cusparseCreateCsr(): dimension mismatch for CUSPARSE_INDEX_32I, cols (4294967294) + base (0) > INT32_MAX (2147483647)</div><div style="margin:0px;font-style:normal;font-variant:normal;font-weight:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue";font-kerning:auto;font-feature-settings:normal;min-height:15px"><br></div><div style="margin:0px;font-style:normal;font-variant:normal;font-weight:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue";font-kerning:auto;font-feature-settings:normal">[0]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------</div><div style="margin:0px;font-style:normal;font-variant:normal;font-weight:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue";font-kerning:auto;font-feature-settings:normal">[0]PETSC ERROR: GPU error</div><div style="margin:0px;font-style:normal;font-variant:normal;font-weight:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue";font-kerning:auto;font-feature-settings:normal">[0]PETSC ERROR: cuSPARSE errorcode 3 (CUSPARSE_STATUS_INVALID_VALUE) : invalid value</div><div style="margin:0px;font-style:normal;font-variant:normal;font-weight:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue";font-kerning:auto;font-feature-settings:normal">[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.</div><div style="margin:0px;font-style:normal;font-variant:normal;font-weight:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue";font-kerning:auto;font-feature-settings:normal">[0]PETSC ERROR: Petsc Release Version 3.19.4, unknown</div><div style="margin:0px;font-style:normal;font-variant:normal;font-weight:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue";font-kerning:auto;font-feature-settings:normal">[0]PETSC ERROR: ./gpu-bug on a<span> </span>named sean-dgx2 by rohany Fri Aug 11 09:43:07 2023</div><div style="margin:0px;font-style:normal;font-variant:normal;font-weight:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue";font-kerning:auto;font-feature-settings:normal">[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</div><div style="margin:0px;font-style:normal;font-variant:normal;font-weight:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue";font-kerning:auto;font-feature-settings:normal">[0]PETSC ERROR: #1 MatSeqAIJCUSPARSECopyToGPU() at /local/home/rohany/petsc/src/mat/impls/aij/seq/seqcusparse/<a href="http://aijcusparse.cu:2503/" target="_blank">aijcusparse.cu:2503</a></div><div style="margin:0px;font-style:normal;font-variant:normal;font-weight:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue";font-kerning:auto;font-feature-settings:normal">[0]PETSC ERROR: #2 MatMultAddKernel_SeqAIJCUSPARSE() at /local/home/rohany/petsc/src/mat/impls/aij/seq/seqcusparse/<a href="http://aijcusparse.cu:3544/" target="_blank">aijcusparse.cu:3544</a></div><div style="margin:0px;font-style:normal;font-variant:normal;font-weight:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue";font-kerning:auto;font-feature-settings:normal">[0]PETSC ERROR: #3 MatMult_SeqAIJCUSPARSE() at /local/home/rohany/petsc/src/mat/impls/aij/seq/seqcusparse/<a href="http://aijcusparse.cu:3485/" target="_blank">aijcusparse.cu:3485</a></div><div style="margin:0px;font-style:normal;font-variant:normal;font-weight:normal;font-stretch:normal;font-size:13px;line-height:normal;font-family:"Helvetica Neue";font-kerning:auto;font-feature-settings:normal">[0]PETSC ERROR: #4 MatMult_MPIAIJCUSPARSE() at /local/home/rohany/petsc/src/mat/impls/aij/mpi/mpicusparse/<a href="http://mpiaijcusparse.cu:452/" target="_blank">mpiaijcusparse.cu:452</a></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"">
</p><div 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</div><div 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"">```</div><div 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></div><div 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,</div><div 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></div><div 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</div></div><div><br></div></div>
</div></blockquote></div><br></div></div></blockquote></div>