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