1 const char help[] = "Test MatCreateDenseFromVecType()\n\n";
2
3 #include <petscdevice_cuda.h>
4 #include <petscmat.h>
5 #include <petscconf.h>
6 #include <assert.h>
7
main(int argc,char ** args)8 int main(int argc, char **args)
9 {
10 Mat A;
11 Vec X;
12 VecType vtype;
13 PetscInt n = 20, lda = PETSC_DECIDE;
14
15 PetscFunctionBeginUser;
16 PetscCall(PetscInitialize(&argc, &args, NULL, help));
17
18 PetscOptionsBegin(PETSC_COMM_WORLD, NULL, "Creating Mat from Vec type example", NULL);
19 PetscCall(PetscOptionsGetInt(NULL, NULL, "-lda", &lda, NULL));
20 PetscOptionsEnd();
21 if (lda > 0) lda += n;
22
23 PetscCall(VecCreate(PETSC_COMM_WORLD, &X));
24 PetscCall(VecSetSizes(X, n, PETSC_DECIDE));
25 PetscCall(VecSetFromOptions(X));
26 PetscCall(VecSetUp(X));
27 PetscCall(VecGetType(X, &vtype));
28
29 PetscCall(MatCreateDenseFromVecType(PETSC_COMM_WORLD, vtype, n, n, PETSC_DECIDE, PETSC_DECIDE, lda, NULL, &A));
30 PetscCall(MatSetFromOptions(A));
31 PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
32 PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
33
34 PetscMemType X_memtype, A_memtype;
35 const PetscScalar *array;
36 PetscCall(VecGetArrayReadAndMemType(X, &array, &X_memtype));
37 PetscCall(VecRestoreArrayReadAndMemType(X, &array));
38 PetscCall(MatDenseGetArrayReadAndMemType(A, &array, &A_memtype));
39 PetscCall(MatDenseRestoreArrayReadAndMemType(A, &array));
40 PetscAssert(A_memtype == X_memtype, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Failed memtype guarantee in MatCreateDenseFromVecType");
41
42 /* test */
43 PetscCall(MatViewFromOptions(A, NULL, "-ex19_mat_view"));
44 PetscCall(MatDestroy(&A));
45 PetscCall(VecDestroy(&X));
46 PetscCall(PetscFinalize());
47 return 0;
48 }
49
50 /*TEST
51 testset:
52 args: -lda {{0 1}} -ex19_mat_view
53 filter: grep -v -i type
54 output_file: output/ex19.out
55
56 test:
57 suffix: cuda
58 requires: cuda
59 args: -vec_type {{cuda mpicuda}}
60
61 test:
62 suffix: hip
63 requires: hip
64 args: -vec_type hip
65
66 test:
67 suffix: standard
68 args: -vec_type standard
69
70 test:
71 suffix: kokkos
72 # we don't have MATDENSESYCL yet
73 requires: kokkos_kernels !sycl
74 args: -vec_type kokkos
75 TEST*/
76