1 static char help[] = "Testing MatCreateSeqBAIJWithArrays() and MatCreateSeqSBAIJWithArrays().\n\n"; 2 3 #include <petscmat.h> 4 5 int main(int argc, char **argv) 6 { 7 Mat A, B, As; 8 const PetscInt *ai, *aj; 9 PetscInt i, j, k, nz, n, asi[] = {0, 2, 3, 4, 6, 7}; 10 PetscInt asj[] = {0, 4, 1, 2, 3, 4, 4}; 11 PetscScalar asa[7], *aa; 12 PetscRandom rctx; 13 PetscMPIInt size; 14 PetscBool flg; 15 16 PetscFunctionBeginUser; 17 PetscCall(PetscInitialize(&argc, &argv, (char *)0, help)); 18 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size)); 19 PetscCheck(size == 1, PETSC_COMM_WORLD, PETSC_ERR_WRONG_MPI_SIZE, "This is a uniprocessor example only!"); 20 21 /* Create a aij matrix for checking */ 22 PetscCall(MatCreateSeqAIJ(PETSC_COMM_SELF, 5, 5, 2, NULL, &A)); 23 PetscCall(PetscRandomCreate(PETSC_COMM_WORLD, &rctx)); 24 PetscCall(PetscRandomSetFromOptions(rctx)); 25 26 k = 0; 27 for (i = 0; i < 5; i++) { 28 nz = asi[i + 1] - asi[i]; /* length of i_th row of A */ 29 for (j = 0; j < nz; j++) { 30 PetscCall(PetscRandomGetValue(rctx, &asa[k])); 31 PetscCall(MatSetValues(A, 1, &i, 1, &asj[k], &asa[k], INSERT_VALUES)); 32 PetscCall(MatSetValues(A, 1, &i, 1, &asj[k], &asa[k], INSERT_VALUES)); 33 if (i != asj[k]) { /* insert symmetric entry */ 34 PetscCall(MatSetValues(A, 1, &asj[k], 1, &i, &asa[k], INSERT_VALUES)); 35 } 36 k++; 37 } 38 } 39 PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY)); 40 PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY)); 41 42 /* Create a baij matrix using MatCreateSeqBAIJWithArrays() */ 43 PetscCall(MatGetRowIJ(A, 0, PETSC_FALSE, PETSC_FALSE, &n, &ai, &aj, &flg)); 44 PetscCall(MatSeqAIJGetArray(A, &aa)); 45 /* WARNING: This sharing is dangerous if either A or B is later assembled */ 46 PetscCall(MatCreateSeqBAIJWithArrays(PETSC_COMM_SELF, 1, 5, 5, (PetscInt *)ai, (PetscInt *)aj, aa, &B)); 47 PetscCall(MatSeqAIJRestoreArray(A, &aa)); 48 PetscCall(MatRestoreRowIJ(A, 0, PETSC_FALSE, PETSC_FALSE, &n, &ai, &aj, &flg)); 49 PetscCall(MatMultEqual(A, B, 10, &flg)); 50 PetscCheck(flg, PETSC_COMM_SELF, PETSC_ERR_ARG_NOTSAMETYPE, "MatMult(A,B) are NOT equal"); 51 52 /* Create a sbaij matrix using MatCreateSeqSBAIJWithArrays() */ 53 PetscCall(MatCreateSeqSBAIJWithArrays(PETSC_COMM_SELF, 1, 5, 5, asi, asj, asa, &As)); 54 PetscCall(MatMultEqual(A, As, 10, &flg)); 55 PetscCheck(flg, PETSC_COMM_SELF, PETSC_ERR_ARG_NOTSAMETYPE, "MatMult(A,As) are NOT equal"); 56 57 /* Free spaces */ 58 PetscCall(PetscRandomDestroy(&rctx)); 59 PetscCall(MatDestroy(&A)); 60 PetscCall(MatDestroy(&B)); 61 PetscCall(MatDestroy(&As)); 62 PetscCall(PetscFinalize()); 63 return 0; 64 } 65 66 /*TEST 67 68 test: 69 TODO: MatCreateSeqBAIJWithArrays() is broken, it leaks imax and ilen arrays 70 71 TEST*/ 72