1 static char help[] = "Test MatDenseGetSubMatrix() on a CUDA matrix.\n"; 2 3 #include <petscmat.h> 4 5 int main(int argc,char **argv) 6 { 7 Mat A,B; 8 PetscScalar *b; 9 PetscInt n=4,lda=5,i,k; 10 PetscBool cuda=PETSC_FALSE; 11 12 PetscCall(PetscInitialize(&argc,&argv,0,help)); 13 PetscCall(PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL)); 14 PetscCall(PetscOptionsGetInt(NULL,NULL,"-lda",&lda,NULL)); 15 PetscCall(PetscOptionsGetBool(NULL,NULL,"-cuda",&cuda,NULL)); 16 PetscCheck(lda>=n,PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"lda %" PetscInt_FMT " < n %" PetscInt_FMT,lda,n); 17 18 #if defined(PETSC_HAVE_CUDA) 19 if (cuda) PetscCall(MatCreateSeqDenseCUDA(PETSC_COMM_SELF,lda,n,NULL,&A)); 20 else 21 #endif 22 PetscCall(MatCreateSeqDense(PETSC_COMM_SELF,lda,n,NULL,&A)); 23 24 for (k=0;k<3;k++) { 25 PetscCall(MatDenseGetSubMatrix(A,0,n,0,n,&B)); 26 PetscCall(MatDenseGetArray(B,&b)); 27 for (i=0; i<n; i++) { 28 b[i+i*lda] = 2.0*(i+1); 29 if (i>0) b[i+(i-1)*lda] = (PetscReal)(k+1); 30 } 31 PetscCall(MatDenseRestoreArray(B,&b)); 32 PetscCall(MatDenseRestoreSubMatrix(A,&B)); 33 PetscCall(MatView(A,NULL)); 34 } 35 36 PetscCall(MatDestroy(&A)); 37 PetscCall(PetscFinalize()); 38 return 0; 39 } 40 41 /*TEST 42 43 testset: 44 output_file: output/ex257_1.out 45 diff_args: -j 46 test: 47 suffix: 1 48 test: 49 suffix: 1_cuda 50 args: -cuda 51 requires: cuda 52 filter: sed -e "s/seqdensecuda/seqdense/" 53 54 TEST*/ 55