1: #include <petscksp.h>

  3: // Keep this driver and the companion Python driver ex100.py in sync.

  5: /* ------------------------------------------------------- */

  7: PetscErrorCode RunTest(void)
  8: {
  9:   PetscInt  N = 100, its = 0;
 10:   PetscBool draw = PETSC_FALSE, test = PETSC_FALSE;
 11:   PetscReal rnorm;
 12:   Mat       A;
 13:   Vec       b, x, r;
 14:   KSP       ksp;
 15:   PC        pc;

 17:   PetscFunctionBegin;
 18:   PetscCall(PetscOptionsGetInt(NULL, NULL, "-N", &N, NULL));
 19:   PetscCall(PetscOptionsGetBool(NULL, NULL, "-test", &test, NULL));
 20:   PetscCall(PetscOptionsGetBool(NULL, NULL, "-draw", &draw, NULL));

 22:   PetscCall(MatCreate(PETSC_COMM_WORLD, &A));
 23:   PetscCall(MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, N, N));
 24:   PetscCall(MatSetType(A, MATPYTHON));
 25:   PetscCall(MatPythonSetType(A, "example100.py:Laplace1D"));
 26:   PetscCall(MatSetUp(A));

 28:   PetscCall(MatCreateVecs(A, &x, &b));
 29:   PetscCall(VecSet(b, 1));

 31:   PetscCall(KSPCreate(PETSC_COMM_WORLD, &ksp));
 32:   PetscCall(KSPSetType(ksp, KSPPYTHON));
 33:   PetscCall(KSPPythonSetType(ksp, "example100.py:ConjGrad"));

 35:   PetscCall(KSPGetPC(ksp, &pc));
 36:   PetscCall(PCSetType(pc, PCPYTHON));
 37:   PetscCall(PCPythonSetType(pc, "example100.py:Jacobi"));

 39:   PetscCall(KSPSetOperators(ksp, A, A));
 40:   PetscCall(KSPSetFromOptions(ksp));
 41:   PetscCall(KSPSolve(ksp, b, x));

 43:   if (test) {
 44:     PetscCall(KSPGetTotalIterations(ksp, &its));
 45:     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Number of KSP iterations = %" PetscInt_FMT "\n", its));
 46:   } else {
 47:     PetscCall(VecDuplicate(b, &r));
 48:     PetscCall(MatMult(A, x, r));
 49:     PetscCall(VecAYPX(r, -1, b));
 50:     PetscCall(VecNorm(r, NORM_2, &rnorm));
 51:     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "error norm = %g\n", (double)rnorm));
 52:     PetscCall(VecDestroy(&r));
 53:   }

 55:   if (draw) {
 56:     PetscCall(VecView(x, PETSC_VIEWER_DRAW_WORLD));
 57:     PetscCall(PetscSleep(2));
 58:   }

 60:   PetscCall(VecDestroy(&x));
 61:   PetscCall(VecDestroy(&b));
 62:   PetscCall(MatDestroy(&A));
 63:   PetscCall(KSPDestroy(&ksp));
 64:   PetscFunctionReturn(PETSC_SUCCESS);
 65: }

 67: /* ------------------------------------------------------- */

 69: static char help[] = "Python-implemented Mat/KSP/PC.\n"
 70:                      "The companion Python driver ex100.py uses the same implementations in example100.py.\n"
 71:                      "Run it with python ex100.py and the same PETSc options as this example.\n\n";

 73: #if !defined(PYTHON_EXE)
 74:   #define PYTHON_EXE 0
 75: #endif
 76: #if !defined(PYTHON_LIB)
 77:   #define PYTHON_LIB 0
 78: #endif

 80: int main(int argc, char *argv[])
 81: {
 82:   PetscFunctionBeginUser;
 83:   PetscCall(PetscInitialize(&argc, &argv, 0, help));
 84:   PetscCall(PetscPythonInitialize(PYTHON_EXE, PYTHON_LIB));
 85:   PetscCall(RunTest());
 86:   PetscCall(PetscPythonPrintError());
 87:   PetscCall(PetscFinalize());
 88:   return 0;
 89: }

 91: /*TEST

 93:     test:
 94:       args: -ksp_monitor
 95:       requires: petsc4py
 96:       localrunfiles: example100.py

 98: TEST*/