1: static char help[] = "Two SuperLU_DIST factorizations alive at once on one communicator.\n"
  2:                      "Exercises the 'communicator attribute already in use' path of\n"
  3:                      "MatLUFactorSymbolic_SuperLU_DIST and the later attribute delete callback.\n\n";

  5: #include <petscksp.h>

  7: static PetscErrorCode MakeMat(PetscInt n, Mat *A)
  8: {
  9:   PetscInt rs, re;

 11:   PetscFunctionBeginUser;
 12:   PetscCall(MatCreateAIJ(PETSC_COMM_WORLD, PETSC_DECIDE, PETSC_DECIDE, n, n, 3, NULL, 2, NULL, A));
 13:   PetscCall(MatGetOwnershipRange(*A, &rs, &re));
 14:   for (PetscInt i = rs; i < re; i++) {
 15:     PetscScalar d = 4.0, o = -1.0;
 16:     PetscInt    im = i - 1, ip = i + 1;

 18:     PetscCall(MatSetValues(*A, 1, &i, 1, &i, &d, INSERT_VALUES));
 19:     if (i > 0) PetscCall(MatSetValues(*A, 1, &i, 1, &im, &o, INSERT_VALUES));
 20:     if (i < n - 1) PetscCall(MatSetValues(*A, 1, &i, 1, &ip, &o, INSERT_VALUES));
 21:   }
 22:   PetscCall(MatAssemblyBegin(*A, MAT_FINAL_ASSEMBLY));
 23:   PetscCall(MatAssemblyEnd(*A, MAT_FINAL_ASSEMBLY));
 24:   PetscFunctionReturn(PETSC_SUCCESS);
 25: }

 27: static PetscErrorCode SolveOnce(Mat A, KSP *ksp)
 28: {
 29:   Vec x, b;

 31:   PetscFunctionBeginUser;
 32:   PetscCall(KSPCreate(PETSC_COMM_WORLD, ksp));
 33:   PetscCall(KSPSetFromOptions(*ksp));
 34:   PetscCall(KSPSetOperators(*ksp, A, A));
 35:   PetscCall(MatCreateVecs(A, &x, &b));
 36:   PetscCall(VecSet(b, 1.0));
 37:   PetscCall(KSPSolve(*ksp, b, x));
 38:   PetscCall(VecDestroy(&x));
 39:   PetscCall(VecDestroy(&b));
 40:   PetscFunctionReturn(PETSC_SUCCESS);
 41: }

 43: int main(int argc, char **argv)
 44: {
 45:   Mat A1, A2, B;
 46:   KSP k1, k2, k3, kb;

 48:   PetscFunctionBeginUser;
 49:   PetscCall(PetscInitialize(&argc, &argv, NULL, help));

 51:   /* first factorization: grid stored in the communicator attribute */
 52:   PetscCall(MakeMat(18, &B));
 53:   PetscCall(SolveOnce(B, &kb));
 54:   PetscCall(KSPDestroy(&kb));
 55:   PetscCall(MatDestroy(&B));

 57:   /* three factorizations kept alive at the same time */
 58:   PetscCall(MakeMat(18, &A1));
 59:   PetscCall(MakeMat(18, &A2));
 60:   PetscCall(SolveOnce(A1, &k1));
 61:   PetscCall(SolveOnce(A2, &k2));
 62:   PetscCall(SolveOnce(A1, &k3));

 64:   /* a fourth one created and destroyed while the others are alive */
 65:   PetscCall(MakeMat(18, &B));
 66:   PetscCall(SolveOnce(B, &kb));
 67:   PetscCall(KSPDestroy(&kb));
 68:   PetscCall(MatDestroy(&B));

 70:   /* destroying the survivors frees the inner communicator and runs the attribute delete callback */
 71:   PetscCall(KSPDestroy(&k1));
 72:   PetscCall(KSPDestroy(&k2));
 73:   PetscCall(KSPDestroy(&k3));
 74:   PetscCall(MatDestroy(&A1));
 75:   PetscCall(MatDestroy(&A2));
 76:   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "done\n"));
 77:   PetscCall(PetscFinalize());
 78:   return 0;
 79: }

 81: /*TEST

 83:   build:
 84:     requires: superlu_dist

 86:   test:
 87:     args: -ksp_type richardson -ksp_max_it 1 -pc_type lu -pc_factor_mat_solver_type superlu_dist

 89:   test:
 90:     suffix: 2
 91:     nsize: 2
 92:     args: -ksp_type richardson -ksp_max_it 1 -pc_type lu -pc_factor_mat_solver_type superlu_dist
 93:     output_file: output/ex91_1.out

 95: TEST*/