Add a silly residual monitor

pism · Oct 14, 2020 · 3ef667e · 3ef667e
1 parent 66b0c06
commit 3ef667e
Showing 1 changed file with 34 additions and 0 deletions.
diff --git a/src/stressbalance/blatter/Poisson3.cc b/src/stressbalance/blatter/Poisson3.cc
@@ -521,6 +521,38 @@ PetscErrorCode p3_coarsening_hook(DM dm_fine, DM dm_coarse, void *ctx) {
   return 0;
 }
 
+PetscErrorCode p3_monitor(SNES snes, PetscInt its, PetscReal norm, void *mctx) {
+  Vec F;
+  DM da;
+  double ***a;
+
+  SNESGetDM(snes, &da);
+  SNESGetFunction(snes, &F, NULL, NULL);
+
+  DMDALocalInfo info;
+  DMDAGetLocalInfo(da, &info);
+
+  DMDAVecGetArray(da, F, &a);
+  {
+    PetscPrintf(PETSC_COMM_SELF, "residual (i = %d)\n", its);
+
+    for (int k = info.zs; k < info.zs + info.zm; k++) {
+      PetscPrintf(PETSC_COMM_SELF, "k = %d\n", k);
+      for (int j = info.ys; j < info.ys + info.ym; j++) {
+        for (int i = info.xs; i < info.xs + info.xm; i++) {
+          PetscPrintf(PETSC_COMM_SELF, "%g ", a[k][j][i]);
+        }
+        PetscPrintf(PETSC_COMM_SELF, "\n");
+      }
+      PetscPrintf(PETSC_COMM_SELF, "\n");
+    }
+    PetscPrintf(PETSC_COMM_SELF, "\n");
+  }
+  DMDAVecRestoreArray(da, F, &a);
+
+  return 0;
+}
+
 PetscErrorCode Poisson3::setup(DM pism_da) {
   PetscErrorCode ierr;
   // DM
@@ -601,6 +633,8 @@ PetscErrorCode Poisson3::setup(DM pism_da) {
                                     &m_callback_data); CHKERRQ(ierr);
 
     ierr = SNESSetFromOptions(m_snes); CHKERRQ(ierr);
+
+    ierr = SNESMonitorSet(m_snes, p3_monitor, NULL, NULL); CHKERRQ(ierr);
   }
 
   // set the initial guess