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