xref: /petsc/src/mat/tests/ex56.c (revision 732aec7a18f2199fb53bb9a2f3aef439a834ce31)
1c4762a1bSJed Brown static char help[] = "Test the use of MatSetValuesBlocked(), MatZeroRows() for rectangular MatBAIJ matrix, test MatSetValuesBlocked() for MatSBAIJ matrix (-test_mat_sbaij).";
2c4762a1bSJed Brown 
3c4762a1bSJed Brown #include <petscmat.h>
4c4762a1bSJed Brown 
main(int argc,char ** args)5d71ae5a4SJacob Faibussowitsch int main(int argc, char **args)
6d71ae5a4SJacob Faibussowitsch {
7c4762a1bSJed Brown   Mat         A;
8c4762a1bSJed Brown   PetscInt    bs = 3, m = 4, n = 6, i, j, val = 10, row[2], col[3], eval, rstart;
9c4762a1bSJed Brown   PetscMPIInt size, rank;
10c4762a1bSJed Brown   PetscScalar x[6][9], y[3][3], one = 1.0;
11c4762a1bSJed Brown   PetscBool   flg, testsbaij = PETSC_FALSE;
12c4762a1bSJed Brown 
13327415f7SBarry Smith   PetscFunctionBeginUser;
14*c8025a54SPierre Jolivet   PetscCall(PetscInitialize(&argc, &args, NULL, help));
159566063dSJacob Faibussowitsch   PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
169566063dSJacob Faibussowitsch   PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
17c4762a1bSJed Brown 
189566063dSJacob Faibussowitsch   PetscCall(PetscOptionsHasName(NULL, NULL, "-test_mat_sbaij", &testsbaij));
19c4762a1bSJed Brown 
20c4762a1bSJed Brown   if (testsbaij) {
219566063dSJacob Faibussowitsch     PetscCall(MatCreateSBAIJ(PETSC_COMM_WORLD, bs, m * bs, n * bs, PETSC_DECIDE, PETSC_DECIDE, 1, NULL, 1, NULL, &A));
22c4762a1bSJed Brown   } else {
239566063dSJacob Faibussowitsch     PetscCall(MatCreateBAIJ(PETSC_COMM_WORLD, bs, m * bs, n * bs, PETSC_DECIDE, PETSC_DECIDE, 1, NULL, 1, NULL, &A));
24c4762a1bSJed Brown   }
259566063dSJacob Faibussowitsch   PetscCall(MatSetOption(A, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_FALSE));
26c4762a1bSJed Brown   eval = 9;
27c4762a1bSJed Brown 
289566063dSJacob Faibussowitsch   PetscCall(PetscOptionsHasName(NULL, NULL, "-ass_extern", &flg));
29c4762a1bSJed Brown   if (flg && (size != 1)) rstart = m * ((rank + 1) % size);
30c4762a1bSJed Brown   else rstart = m * (rank);
31c4762a1bSJed Brown 
329371c9d4SSatish Balay   row[0] = rstart + 0;
339371c9d4SSatish Balay   row[1] = rstart + 2;
349371c9d4SSatish Balay   col[0] = rstart + 0;
359371c9d4SSatish Balay   col[1] = rstart + 1;
369371c9d4SSatish Balay   col[2] = rstart + 3;
37c4762a1bSJed Brown   for (i = 0; i < 6; i++) {
38c4762a1bSJed Brown     for (j = 0; j < 9; j++) x[i][j] = (PetscScalar)val++;
39c4762a1bSJed Brown   }
40c4762a1bSJed Brown 
419566063dSJacob Faibussowitsch   PetscCall(MatSetValuesBlocked(A, 2, row, 3, col, &x[0][0], INSERT_VALUES));
42c4762a1bSJed Brown 
439566063dSJacob Faibussowitsch   PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
449566063dSJacob Faibussowitsch   PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
45c4762a1bSJed Brown 
46c4762a1bSJed Brown   /*
47c4762a1bSJed Brown   This option does not work for rectangular matrices
489566063dSJacob Faibussowitsch   PetscCall(MatSetOption(A,MAT_NEW_NONZERO_LOCATION_ERR,PETSC_TRUE));
49c4762a1bSJed Brown   */
50c4762a1bSJed Brown 
519566063dSJacob Faibussowitsch   PetscCall(MatSetValuesBlocked(A, 2, row, 3, col, &x[0][0], INSERT_VALUES));
52c4762a1bSJed Brown 
53c4762a1bSJed Brown   /* Do another MatSetValues to test the case when only one local block is specified */
54c4762a1bSJed Brown   for (i = 0; i < 3; i++) {
55c4762a1bSJed Brown     for (j = 0; j < 3; j++) y[i][j] = (PetscScalar)(10 + i * eval + j);
56c4762a1bSJed Brown   }
579566063dSJacob Faibussowitsch   PetscCall(MatSetValuesBlocked(A, 1, row, 1, col, &y[0][0], INSERT_VALUES));
589566063dSJacob Faibussowitsch   PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
599566063dSJacob Faibussowitsch   PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
60c4762a1bSJed Brown 
619566063dSJacob Faibussowitsch   PetscCall(PetscOptionsHasName(NULL, NULL, "-zero_rows", &flg));
62c4762a1bSJed Brown   if (flg) {
63c4762a1bSJed Brown     col[0] = rstart * bs + 0;
64c4762a1bSJed Brown     col[1] = rstart * bs + 1;
65c4762a1bSJed Brown     col[2] = rstart * bs + 2;
669566063dSJacob Faibussowitsch     PetscCall(MatZeroRows(A, 3, col, one, 0, 0));
67c4762a1bSJed Brown   }
68c4762a1bSJed Brown 
699566063dSJacob Faibussowitsch   PetscCall(MatView(A, PETSC_VIEWER_STDOUT_WORLD));
70c4762a1bSJed Brown 
719566063dSJacob Faibussowitsch   PetscCall(MatDestroy(&A));
729566063dSJacob Faibussowitsch   PetscCall(PetscFinalize());
73b122ec5aSJacob Faibussowitsch   return 0;
74c4762a1bSJed Brown }
75c4762a1bSJed Brown 
76c4762a1bSJed Brown /*TEST
77c4762a1bSJed Brown 
78c4762a1bSJed Brown    test:
798cc725e6SPierre Jolivet       filter: grep -v " MPI process"
80c4762a1bSJed Brown 
81c4762a1bSJed Brown    test:
82c4762a1bSJed Brown       suffix: 4
83c4762a1bSJed Brown       nsize: 3
84c4762a1bSJed Brown       args: -ass_extern
858cc725e6SPierre Jolivet       filter: grep -v " MPI process"
86c4762a1bSJed Brown 
87c4762a1bSJed Brown    test:
88c4762a1bSJed Brown       suffix: 5
89c4762a1bSJed Brown       nsize: 3
90c4762a1bSJed Brown       args: -ass_extern -zero_rows
918cc725e6SPierre Jolivet       filter: grep -v " MPI process"
92c4762a1bSJed Brown 
93c4762a1bSJed Brown TEST*/
94