1*ad6ad4e1SHansol Suh const char help[] = "Test VecCreateMatDense()\n\n"; 2*ad6ad4e1SHansol Suh 3*ad6ad4e1SHansol Suh #include <petscdevice_cuda.h> 4*ad6ad4e1SHansol Suh #include <petscmat.h> 5*ad6ad4e1SHansol Suh #include <petscconf.h> 6*ad6ad4e1SHansol Suh #include <assert.h> 7*ad6ad4e1SHansol Suh 8*ad6ad4e1SHansol Suh int main(int argc, char **args) 9*ad6ad4e1SHansol Suh { 10*ad6ad4e1SHansol Suh Mat A; 11*ad6ad4e1SHansol Suh Vec X; 12*ad6ad4e1SHansol Suh PetscInt N = 20; 13*ad6ad4e1SHansol Suh 14*ad6ad4e1SHansol Suh PetscFunctionBeginUser; 15*ad6ad4e1SHansol Suh PetscCall(PetscInitialize(&argc, &args, NULL, help)); 16*ad6ad4e1SHansol Suh 17*ad6ad4e1SHansol Suh PetscOptionsBegin(PETSC_COMM_WORLD, NULL, "Creating Mat from Vec example", NULL); 18*ad6ad4e1SHansol Suh PetscCall(PetscOptionsGetInt(NULL, NULL, "-n", &N, NULL)); 19*ad6ad4e1SHansol Suh PetscOptionsEnd(); 20*ad6ad4e1SHansol Suh 21*ad6ad4e1SHansol Suh PetscCall(VecCreate(PETSC_COMM_WORLD, &X)); 22*ad6ad4e1SHansol Suh PetscCall(VecSetSizes(X, PETSC_DECIDE, N)); 23*ad6ad4e1SHansol Suh PetscCall(VecSetFromOptions(X)); 24*ad6ad4e1SHansol Suh PetscCall(VecSetUp(X)); 25*ad6ad4e1SHansol Suh 26*ad6ad4e1SHansol Suh PetscCall(VecCreateMatDense(X, PETSC_DECIDE, PETSC_DECIDE, N, N, NULL, &A)); 27*ad6ad4e1SHansol Suh PetscCall(MatSetFromOptions(A)); 28*ad6ad4e1SHansol Suh PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY)); 29*ad6ad4e1SHansol Suh PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY)); 30*ad6ad4e1SHansol Suh 31*ad6ad4e1SHansol Suh MPI_Comm X_comm = PetscObjectComm((PetscObject)X); 32*ad6ad4e1SHansol Suh MPI_Comm A_comm = PetscObjectComm((PetscObject)X); 33*ad6ad4e1SHansol Suh PetscMPIInt comp; 34*ad6ad4e1SHansol Suh PetscCallMPI(MPI_Comm_compare(X_comm, A_comm, &comp)); 35*ad6ad4e1SHansol Suh PetscAssert(comp == MPI_IDENT || comp == MPI_CONGRUENT, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Failed communicator guarantee in MatCreateDenseMatchingVec()"); 36*ad6ad4e1SHansol Suh 37*ad6ad4e1SHansol Suh PetscMemType X_memtype, A_memtype; 38*ad6ad4e1SHansol Suh const PetscScalar *array; 39*ad6ad4e1SHansol Suh PetscCall(VecGetArrayReadAndMemType(X, &array, &X_memtype)); 40*ad6ad4e1SHansol Suh PetscCall(VecRestoreArrayReadAndMemType(X, &array)); 41*ad6ad4e1SHansol Suh PetscCall(MatDenseGetArrayReadAndMemType(A, &array, &A_memtype)); 42*ad6ad4e1SHansol Suh PetscCall(MatDenseRestoreArrayReadAndMemType(A, &array)); 43*ad6ad4e1SHansol Suh PetscAssert(A_memtype == X_memtype, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Failed memtype guarantee in MatCreateDenseMatchingVec()"); 44*ad6ad4e1SHansol Suh 45*ad6ad4e1SHansol Suh /* test */ 46*ad6ad4e1SHansol Suh PetscCall(MatViewFromOptions(A, NULL, "-ex19_mat_view")); 47*ad6ad4e1SHansol Suh PetscCall(MatDestroy(&A)); 48*ad6ad4e1SHansol Suh PetscCall(VecDestroy(&X)); 49*ad6ad4e1SHansol Suh PetscCall(PetscFinalize()); 50*ad6ad4e1SHansol Suh return 0; 51*ad6ad4e1SHansol Suh } 52*ad6ad4e1SHansol Suh 53*ad6ad4e1SHansol Suh /*TEST 54*ad6ad4e1SHansol Suh 55*ad6ad4e1SHansol Suh test: 56*ad6ad4e1SHansol Suh suffix: cuda 57*ad6ad4e1SHansol Suh requires: cuda 58*ad6ad4e1SHansol Suh args: -vec_type cuda -ex19_mat_view 59*ad6ad4e1SHansol Suh 60*ad6ad4e1SHansol Suh test: 61*ad6ad4e1SHansol Suh suffix: mpicuda 62*ad6ad4e1SHansol Suh requires: cuda 63*ad6ad4e1SHansol Suh args: -vec_type mpicuda -ex19_mat_view 64*ad6ad4e1SHansol Suh 65*ad6ad4e1SHansol Suh test: 66*ad6ad4e1SHansol Suh suffix: hip 67*ad6ad4e1SHansol Suh requires: hip 68*ad6ad4e1SHansol Suh args: -vec_type hip -ex19_mat_view 69*ad6ad4e1SHansol Suh 70*ad6ad4e1SHansol Suh test: 71*ad6ad4e1SHansol Suh suffix: standard 72*ad6ad4e1SHansol Suh args: -vec_type standard -ex19_mat_view 73*ad6ad4e1SHansol Suh 74*ad6ad4e1SHansol Suh test: 75*ad6ad4e1SHansol Suh suffix: kokkos_cuda 76*ad6ad4e1SHansol Suh requires: kokkos kokkos_kernels cuda 77*ad6ad4e1SHansol Suh args: -vec_type kokkos -ex19_mat_view 78*ad6ad4e1SHansol Suh 79*ad6ad4e1SHansol Suh test: 80*ad6ad4e1SHansol Suh suffix: kokkos_hip 81*ad6ad4e1SHansol Suh requires: kokkos kokkos_kernels hip 82*ad6ad4e1SHansol Suh args: -vec_type kokkos -ex19_mat_view 83*ad6ad4e1SHansol Suh 84*ad6ad4e1SHansol Suh test: 85*ad6ad4e1SHansol Suh suffix: kokkos 86*ad6ad4e1SHansol Suh requires: kokkos kokkos_kernels !cuda !hip 87*ad6ad4e1SHansol Suh args: -vec_type kokkos -ex19_mat_view 88*ad6ad4e1SHansol Suh TEST*/ 89