Actual source code: ex243.c
1: static char help[] = "Test conversion of ScaLAPACK matrices.\n\n";
3: #include <petscmat.h>
5: int main(int argc, char **argv)
6: {
7: Mat A, A_scalapack;
8: PetscInt i, j, M = 10, N = 5, nloc, mloc, nrows, ncols;
9: PetscMPIInt rank, size;
10: IS isrows, iscols;
11: const PetscInt *rows, *cols;
12: PetscScalar *v;
13: MatType type;
14: PetscBool isDense, isAIJ, flg;
17: PetscInitialize(&argc, &argv, (char *)0, help);
18: MPI_Comm_size(PETSC_COMM_WORLD, &size);
19: MPI_Comm_rank(PETSC_COMM_WORLD, &rank);
20: PetscOptionsGetInt(NULL, NULL, "-M", &M, NULL);
21: PetscOptionsGetInt(NULL, NULL, "-N", &N, NULL);
23: /* Create a matrix */
24: MatCreate(PETSC_COMM_WORLD, &A);
25: mloc = PETSC_DECIDE;
26: PetscSplitOwnershipEqual(PETSC_COMM_WORLD, &mloc, &M);
27: nloc = PETSC_DECIDE;
28: PetscSplitOwnershipEqual(PETSC_COMM_WORLD, &nloc, &N);
29: MatSetSizes(A, mloc, nloc, M, N);
30: MatSetType(A, MATDENSE);
31: MatSetFromOptions(A);
32: MatSetUp(A);
34: /* Set local matrix entries */
35: MatGetOwnershipIS(A, &isrows, &iscols);
36: ISGetLocalSize(isrows, &nrows);
37: ISGetIndices(isrows, &rows);
38: ISGetLocalSize(iscols, &ncols);
39: ISGetIndices(iscols, &cols);
40: PetscMalloc1(nrows * ncols, &v);
42: for (i = 0; i < nrows; i++) {
43: for (j = 0; j < ncols; j++) {
44: if (size == 1) {
45: v[i * ncols + j] = (PetscScalar)(i + j);
46: } else {
47: v[i * ncols + j] = (PetscScalar)rank + j * 0.1;
48: }
49: }
50: }
51: MatSetValues(A, nrows, rows, ncols, cols, v, INSERT_VALUES);
52: MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
53: MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);
55: /* Test MatSetValues() by converting A to A_scalapack */
56: MatGetType(A, &type);
57: if (size == 1) {
58: PetscObjectTypeCompare((PetscObject)A, MATSEQDENSE, &isDense);
59: PetscObjectTypeCompare((PetscObject)A, MATSEQAIJ, &isAIJ);
60: } else {
61: PetscObjectTypeCompare((PetscObject)A, MATMPIDENSE, &isDense);
62: PetscObjectTypeCompare((PetscObject)A, MATMPIAIJ, &isAIJ);
63: }
65: if (isDense || isAIJ) {
66: Mat Aexplicit;
67: MatConvert(A, MATSCALAPACK, MAT_INITIAL_MATRIX, &A_scalapack);
68: MatComputeOperator(A_scalapack, isAIJ ? MATAIJ : MATDENSE, &Aexplicit);
69: MatMultEqual(Aexplicit, A_scalapack, 5, &flg);
71: MatDestroy(&Aexplicit);
73: /* Test MAT_REUSE_MATRIX which is only supported for inplace conversion */
74: MatConvert(A, MATSCALAPACK, MAT_INPLACE_MATRIX, &A);
75: MatMultEqual(A_scalapack, A, 5, &flg);
77: MatDestroy(&A_scalapack);
78: }
80: ISRestoreIndices(isrows, &rows);
81: ISRestoreIndices(iscols, &cols);
82: ISDestroy(&isrows);
83: ISDestroy(&iscols);
84: PetscFree(v);
85: MatDestroy(&A);
86: PetscFinalize();
87: return 0;
88: }
90: /*TEST
92: build:
93: requires: scalapack
95: test:
96: nsize: 6
98: test:
99: suffix: 2
100: nsize: 6
101: args: -mat_type aij
102: output_file: output/ex243_1.out
104: test:
105: suffix: 3
106: nsize: 6
107: args: -mat_type scalapack
108: output_file: output/ex243_1.out
110: TEST*/