[NL] Copy current Newton applyBCs for PETSc-SNES

The system solved by PETSc-SNES is different one and so far
the methods could be reused. With upcoming changes to the
Newton's version of applying Dirichlet-BCs this is not
going to work anymore.

NumLib/ODESolver/NonlinearSystem.h

@@ -78,6 +78,13 @@ public:
         GlobalMatrix& Jac, GlobalVector& res,
         GlobalVector& minus_delta_x) const = 0;
 
+    //! Apply known solutions to the linearized equation system
+    //! \f$ \mathit{Jac} \cdot (-\Delta x) = \mathit{res} \f$.
+    //! \pre computeKnownSolutions() must have been called before.
+    virtual void applyKnownSolutionsPETScSNES(GlobalMatrix& Jac,
+                                              GlobalVector& res,
+                                              GlobalVector& x) const = 0;
+
     virtual void updateConstraints(GlobalVector& /*lower*/,
                                    GlobalVector& /*upper*/,
                                    int /*process_id*/) = 0;

NumLib/ODESolver/PETScNonlinearSolver.cpp

@@ -67,7 +67,7 @@ PetscErrorCode updateResidual(SNES /*snes*/, Vec x, Vec petsc_r,
     context->J->finalizeAssembly();
 
     context->system->getJacobian(*context->J);
-    context->system->applyKnownSolutionsNewton(
+    context->system->applyKnownSolutionsPETScSNES(
         *context->J, *context->r, *context->x[context->process_id]);
 
     VecCopy(context->r->getRawVector(), petsc_r);

NumLib/ODESolver/TimeDiscretizedODESystem.cpp

@@ -152,6 +152,28 @@ void TimeDiscretizedODESystem<ODESystemTag::FirstOrderImplicitQuasilinear,
     MathLib::applyKnownSolution(Jac, res, minus_delta_x, ids, values);
 }
 
+void TimeDiscretizedODESystem<ODESystemTag::FirstOrderImplicitQuasilinear,
+                              NonlinearSolverTag::Newton>::
+    applyKnownSolutionsPETScSNES(GlobalMatrix& Jac, GlobalVector& res,
+                                 GlobalVector& x) const
+{
+    if (!_known_solutions)
+    {
+        return;
+    }
+
+    using IndexType = MathLib::MatrixVectorTraits<GlobalMatrix>::Index;
+    std::vector<IndexType> ids;
+    for (auto const& bc : *_known_solutions)
+    {
+        ids.insert(end(ids), begin(bc.ids), end(bc.ids));
+    }
+
+    // For the Newton method the values must be zero
+    std::vector<double> values(ids.size(), 0);
+    MathLib::applyKnownSolution(Jac, res, x, ids, values);
+}
+
 TimeDiscretizedODESystem<ODESystemTag::FirstOrderImplicitQuasilinear,
                          NonlinearSolverTag::Picard>::
     TimeDiscretizedODESystem(const int process_id, ODE& ode,

NumLib/ODESolver/TimeDiscretizedODESystem.h

@@ -108,6 +108,9 @@ public:
     void applyKnownSolutionsNewton(GlobalMatrix& Jac, GlobalVector& res,
                                    GlobalVector& minus_delta_x) const override;
 
+    void applyKnownSolutionsPETScSNES(GlobalMatrix& Jac, GlobalVector& res,
+                                      GlobalVector& x) const override;
+
     void updateConstraints(GlobalVector& lower, GlobalVector& upper,
                            int const process_id) override
     {