1 static char help[] = "Test device/host memory allocation in MatDenseSeqCUDA()\n\n"; 2 3 /* Contributed by: Victor Eijkhout <eijkhout@tacc.utexas.edu> */ 4 5 #include <petscmat.h> 6 int main(int argc, char** argv) 7 { 8 PetscErrorCode ierr; 9 PetscInt global_size=100; 10 Mat cuda_matrix; 11 Vec input,output; 12 MPI_Comm comm = PETSC_COMM_SELF; 13 PetscReal nrm = 1; 14 15 ierr = PetscInitialize(&argc,&argv,NULL,help);if (ierr) return ierr; 16 CHKERRQ(MatCreateDenseCUDA(comm,global_size,global_size,global_size,global_size,NULL,&cuda_matrix)); 17 CHKERRQ(MatAssemblyBegin(cuda_matrix,MAT_FINAL_ASSEMBLY)); 18 CHKERRQ(MatAssemblyEnd(cuda_matrix,MAT_FINAL_ASSEMBLY)); 19 20 CHKERRQ(VecCreateSeqCUDA(comm,global_size,&input)); 21 CHKERRQ(VecDuplicate(input,&output)); 22 CHKERRQ(VecSet(input,1.)); 23 CHKERRQ(VecSet(output,2.)); 24 CHKERRQ(MatMult(cuda_matrix,input,output)); 25 CHKERRQ(VecNorm(output,NORM_2,&nrm)); 26 PetscCheckFalse(nrm > PETSC_SMALL,PETSC_COMM_SELF,PETSC_ERR_PLIB,"PETSc generated wrong result. Should be 0, but is %g",(double)nrm); 27 CHKERRQ(VecDestroy(&input)); 28 CHKERRQ(VecDestroy(&output)); 29 CHKERRQ(MatDestroy(&cuda_matrix)); 30 ierr = PetscFinalize(); 31 return ierr; 32 } 33 34 /*TEST 35 build: 36 requires: cuda 37 38 test: 39 nsize: 1 40 41 TEST*/ 42