1c4762a1bSJed Brown static char help[] = "Tests MatTranspose() and MatEqual() for MPIAIJ matrices.\n\n";
2c4762a1bSJed Brown
3c4762a1bSJed Brown #include <petscmat.h>
4c4762a1bSJed Brown
main(int argc,char ** argv)5d71ae5a4SJacob Faibussowitsch int main(int argc, char **argv)
6d71ae5a4SJacob Faibussowitsch {
7c4762a1bSJed Brown Mat A, B;
8c4762a1bSJed Brown PetscInt m = 7, n, i, rstart, rend, cols[3];
9c4762a1bSJed Brown PetscScalar v[3];
10c4762a1bSJed Brown PetscBool equal;
11c4762a1bSJed Brown const char *eq[2];
12c4762a1bSJed Brown
13327415f7SBarry Smith PetscFunctionBeginUser;
14*c8025a54SPierre Jolivet PetscCall(PetscInitialize(&argc, &argv, NULL, help));
159566063dSJacob Faibussowitsch PetscCall(PetscOptionsGetInt(NULL, NULL, "-m", &m, NULL));
16c4762a1bSJed Brown n = m;
17c4762a1bSJed Brown
18c4762a1bSJed Brown /* ------- Assemble matrix, --------- */
19c4762a1bSJed Brown
209566063dSJacob Faibussowitsch PetscCall(MatCreateAIJ(PETSC_COMM_WORLD, PETSC_DECIDE, PETSC_DECIDE, m, n, 0, 0, 0, 0, &A));
219566063dSJacob Faibussowitsch PetscCall(MatSetOption(A, MAT_NEW_NONZERO_LOCATION_ERR, PETSC_FALSE));
229566063dSJacob Faibussowitsch PetscCall(MatGetOwnershipRange(A, &rstart, &rend));
23c4762a1bSJed Brown if (!rstart) {
24c4762a1bSJed Brown cols[0] = 0;
25c4762a1bSJed Brown cols[1] = 1;
269371c9d4SSatish Balay v[0] = 2.0;
279371c9d4SSatish Balay v[1] = -1.0;
289566063dSJacob Faibussowitsch PetscCall(MatSetValues(A, 1, &rstart, 2, cols, v, INSERT_VALUES));
29c4762a1bSJed Brown rstart++;
30c4762a1bSJed Brown }
31c4762a1bSJed Brown if (rend == m) {
32c4762a1bSJed Brown rend--;
33c4762a1bSJed Brown cols[0] = rend - 1;
34c4762a1bSJed Brown cols[1] = rend;
359371c9d4SSatish Balay v[0] = -1.0;
369371c9d4SSatish Balay v[1] = 2.0;
379566063dSJacob Faibussowitsch PetscCall(MatSetValues(A, 1, &rend, 2, cols, v, INSERT_VALUES));
38c4762a1bSJed Brown }
399371c9d4SSatish Balay v[0] = -1.0;
409371c9d4SSatish Balay v[1] = 2.0;
419371c9d4SSatish Balay v[2] = -1.0;
42c4762a1bSJed Brown for (i = rstart; i < rend; i++) {
43c4762a1bSJed Brown cols[0] = i - 1;
44c4762a1bSJed Brown cols[1] = i;
45c4762a1bSJed Brown cols[2] = i + 1;
469566063dSJacob Faibussowitsch PetscCall(MatSetValues(A, 1, &i, 3, cols, v, INSERT_VALUES));
47c4762a1bSJed Brown }
489566063dSJacob Faibussowitsch PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
499566063dSJacob Faibussowitsch PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
50c4762a1bSJed Brown
519566063dSJacob Faibussowitsch PetscCall(MatTranspose(A, MAT_INITIAL_MATRIX, &B));
52c4762a1bSJed Brown
539566063dSJacob Faibussowitsch PetscCall(MatEqual(A, B, &equal));
54c4762a1bSJed Brown
55c4762a1bSJed Brown eq[0] = "not equal";
56c4762a1bSJed Brown eq[1] = "equal";
579566063dSJacob Faibussowitsch PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Matrices are %s\n", eq[equal]));
58c4762a1bSJed Brown
599566063dSJacob Faibussowitsch PetscCall(MatTranspose(A, MAT_REUSE_MATRIX, &B));
609566063dSJacob Faibussowitsch PetscCall(MatEqual(A, B, &equal));
619566063dSJacob Faibussowitsch if (!equal) PetscCall(PetscPrintf(PETSC_COMM_WORLD, "MatTranspose with MAT_REUSE_MATRIX failed"));
62c4762a1bSJed Brown
63c4762a1bSJed Brown /* Free data structures */
649566063dSJacob Faibussowitsch PetscCall(MatDestroy(&A));
659566063dSJacob Faibussowitsch PetscCall(MatDestroy(&B));
66c4762a1bSJed Brown
679566063dSJacob Faibussowitsch PetscCall(PetscFinalize());
68b122ec5aSJacob Faibussowitsch return 0;
69c4762a1bSJed Brown }
70c4762a1bSJed Brown
71c4762a1bSJed Brown /*TEST
72c4762a1bSJed Brown
73c4762a1bSJed Brown test:
74c4762a1bSJed Brown
75c4762a1bSJed Brown test:
76c4762a1bSJed Brown suffix: 2
77c4762a1bSJed Brown nsize: 2
78c4762a1bSJed Brown output_file: output/ex58_1.out
79c4762a1bSJed Brown
80c4762a1bSJed Brown TEST*/
81