xref: /petsc/src/mat/tests/ex108.c (revision a916153043bedfa1a5cd23fc2b8cd7320ba6b2e5)
1c4762a1bSJed Brown static char help[] = "Testing MatCreateSeqBAIJWithArrays() and MatCreateSeqSBAIJWithArrays().\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, As;
8c4762a1bSJed Brown   const PetscInt *ai, *aj;
9c4762a1bSJed Brown   PetscInt        i, j, k, nz, n, asi[] = {0, 2, 3, 4, 6, 7};
10c4762a1bSJed Brown   PetscInt        asj[] = {0, 4, 1, 2, 3, 4, 4};
11c4762a1bSJed Brown   PetscScalar     asa[7], *aa;
12c4762a1bSJed Brown   PetscRandom     rctx;
13c4762a1bSJed Brown   PetscMPIInt     size;
14c4762a1bSJed Brown   PetscBool       flg;
15c4762a1bSJed Brown 
16327415f7SBarry Smith   PetscFunctionBeginUser;
17c8025a54SPierre Jolivet   PetscCall(PetscInitialize(&argc, &argv, NULL, help));
189566063dSJacob Faibussowitsch   PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
19be096a46SBarry Smith   PetscCheck(size == 1, PETSC_COMM_WORLD, PETSC_ERR_WRONG_MPI_SIZE, "This is a uniprocessor example only!");
20c4762a1bSJed Brown 
21c4762a1bSJed Brown   /* Create a aij matrix for checking */
229566063dSJacob Faibussowitsch   PetscCall(MatCreateSeqAIJ(PETSC_COMM_SELF, 5, 5, 2, NULL, &A));
23*2fc555ebSPierre Jolivet   PetscCall(MatSetOption(A, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_FALSE));
249566063dSJacob Faibussowitsch   PetscCall(PetscRandomCreate(PETSC_COMM_WORLD, &rctx));
259566063dSJacob Faibussowitsch   PetscCall(PetscRandomSetFromOptions(rctx));
26c4762a1bSJed Brown 
27c4762a1bSJed Brown   k = 0;
28c4762a1bSJed Brown   for (i = 0; i < 5; i++) {
29c4762a1bSJed Brown     nz = asi[i + 1] - asi[i]; /* length of i_th row of A */
30c4762a1bSJed Brown     for (j = 0; j < nz; j++) {
319566063dSJacob Faibussowitsch       PetscCall(PetscRandomGetValue(rctx, &asa[k]));
329566063dSJacob Faibussowitsch       PetscCall(MatSetValues(A, 1, &i, 1, &asj[k], &asa[k], INSERT_VALUES));
33c4762a1bSJed Brown       if (i != asj[k]) { /* insert symmetric entry */
349566063dSJacob Faibussowitsch         PetscCall(MatSetValues(A, 1, &asj[k], 1, &i, &asa[k], INSERT_VALUES));
35c4762a1bSJed Brown       }
36c4762a1bSJed Brown       k++;
37c4762a1bSJed Brown     }
38c4762a1bSJed Brown   }
399566063dSJacob Faibussowitsch   PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
409566063dSJacob Faibussowitsch   PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
41c4762a1bSJed Brown 
42c4762a1bSJed Brown   /* Create a baij matrix using MatCreateSeqBAIJWithArrays() */
439566063dSJacob Faibussowitsch   PetscCall(MatGetRowIJ(A, 0, PETSC_FALSE, PETSC_FALSE, &n, &ai, &aj, &flg));
449566063dSJacob Faibussowitsch   PetscCall(MatSeqAIJGetArray(A, &aa));
45c4762a1bSJed Brown   /* WARNING: This sharing is dangerous if either A or B is later assembled */
469566063dSJacob Faibussowitsch   PetscCall(MatCreateSeqBAIJWithArrays(PETSC_COMM_SELF, 1, 5, 5, (PetscInt *)ai, (PetscInt *)aj, aa, &B));
479566063dSJacob Faibussowitsch   PetscCall(MatSeqAIJRestoreArray(A, &aa));
489566063dSJacob Faibussowitsch   PetscCall(MatRestoreRowIJ(A, 0, PETSC_FALSE, PETSC_FALSE, &n, &ai, &aj, &flg));
499566063dSJacob Faibussowitsch   PetscCall(MatMultEqual(A, B, 10, &flg));
5028b400f6SJacob Faibussowitsch   PetscCheck(flg, PETSC_COMM_SELF, PETSC_ERR_ARG_NOTSAMETYPE, "MatMult(A,B) are NOT equal");
51c4762a1bSJed Brown 
52c4762a1bSJed Brown   /* Create a sbaij matrix using MatCreateSeqSBAIJWithArrays() */
539566063dSJacob Faibussowitsch   PetscCall(MatCreateSeqSBAIJWithArrays(PETSC_COMM_SELF, 1, 5, 5, asi, asj, asa, &As));
549566063dSJacob Faibussowitsch   PetscCall(MatMultEqual(A, As, 10, &flg));
5528b400f6SJacob Faibussowitsch   PetscCheck(flg, PETSC_COMM_SELF, PETSC_ERR_ARG_NOTSAMETYPE, "MatMult(A,As) are NOT equal");
56c4762a1bSJed Brown 
57c4762a1bSJed Brown   /* Free spaces */
589566063dSJacob Faibussowitsch   PetscCall(PetscRandomDestroy(&rctx));
599566063dSJacob Faibussowitsch   PetscCall(MatDestroy(&A));
609566063dSJacob Faibussowitsch   PetscCall(MatDestroy(&B));
619566063dSJacob Faibussowitsch   PetscCall(MatDestroy(&As));
629566063dSJacob Faibussowitsch   PetscCall(PetscFinalize());
63b122ec5aSJacob Faibussowitsch   return 0;
64c4762a1bSJed Brown }
65ceb5bf51SJacob Faibussowitsch 
66ceb5bf51SJacob Faibussowitsch /*TEST
67ceb5bf51SJacob Faibussowitsch 
68ceb5bf51SJacob Faibussowitsch   test:
69*2fc555ebSPierre Jolivet     output_file: output/empty.out
70ceb5bf51SJacob Faibussowitsch 
71ceb5bf51SJacob Faibussowitsch TEST*/
72