From 8193c8883e194451561d3558fab206c0619eb8e2 Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Fri, 21 Apr 2023 18:08:37 -0600 Subject: [PATCH] Also update some commentary. --- examples/step-57/step-57.cc | 28 +++++++++++++++++++++------- 1 file changed, 21 insertions(+), 7 deletions(-) diff --git a/examples/step-57/step-57.cc b/examples/step-57/step-57.cc index a282ea66b5..d8f139e0ee 100644 --- a/examples/step-57/step-57.cc +++ b/examples/step-57/step-57.cc @@ -487,18 +487,32 @@ namespace Step57 } } + // If we were asked to assemble the Newton matrix, then we also built + // a pressure mass matrix in the bottom right block of the matrix. + // We only need this for the preconditioner, so we need to copy it + // in into a separate matrix object, followed by zeroing out this + // block in the Newton matrix. + // + // Note that settings this bottom right block to zero is not identical to + // not assembling anything in this block, because applying boundary values + // and hanging node constraints (in the + // `constraints_used.distribute_local_to_global()` call above) puts entries + // into this block. As a consequence, setting the $(1,1)$ block to zero + // below does not result in what would have happened if we had just not + // assembled a pressure mass matrix in that block to begin with. + // + // The difference is that if we had not assembled anything in this block, + // dealing with constraint degrees of freedom would have put entries on + // the diagonal of the $(1,1)$ block whereas the last operation below, + // zeroing out the entire block, results in a system matrix with + // rows and columns that are completely empty. In other words, the + // linear problem is singular. Luckily, however, the FGMRES solver we + // use appears to handle these rows and columns without any problem. if (assemble_matrix) { - // Finally we move pressure mass matrix into a separate matrix: pressure_mass_matrix.reinit(sparsity_pattern.block(1, 1)); pressure_mass_matrix.copy_from(system_matrix.block(1, 1)); - // Note that settings this pressure block to zero is not identical to - // not assembling anything in this block, because this operation here - // will (incorrectly) delete diagonal entries that come in from - // hanging node constraints for pressure DoFs. This means that our - // whole system matrix will have rows that are completely - // zero. Luckily, FGMRES handles these rows without any problem. system_matrix.block(1, 1) = 0; } } -- 2.39.5