Use BCGS_AMG_SSOR linear solver for FV computations

This selects the BCGS_AMG_SSOR linear solver backend for FV computations
in both the Richards and the Transport model.

CHANGELOG.md

@@ -2,7 +2,10 @@
 
 ## Unreleased
 
-<!-- Add new entries here! -->
+### Changed
+
+* Linerar solver for finite volumes changed from `AMG_4_DG` to
+    `BCGS_AMG_SSOR` !204
 
 ## 2.0.0 (2020-05-14)
 

dune/dorie/model/richards/richards.cc

@@ -154,9 +154,20 @@ void ModelRichards<Traits>::operator_setup()
 	igo = std::make_unique<IGO>(*go0,*go1);
 
 	// --- Solvers ---
+
+	// FV linear solver
+	if constexpr (order == 0) {
+		// Solver args: GFS, maxiter, verbose, reuse, superlu
+		ls = std::make_unique<LS>(*gfs, 1000, 0, true, true);
+	}
+	// DG linear solver
+	else {
 		lsgfs = std::make_unique<LSGFS>(LSGFSHelper::create(gv));
 		lscc = std::make_unique<LSCC>();
-	ls = std::make_unique<LS>(*igo,*cc,*lsgfs,*lscc,1000,0,true,true);
+		// Solver args: GO, constraints, LSGFS, LSGFS constraints, maxiter, verbose,
+		//              reuse, superlu
+		ls = std::make_unique<LS>(*igo, *cc, *lsgfs, *lscc, 1000, 0, true, true);
+	}
 
 	pdesolver = std::make_unique<PDESOLVER>(*igo,*ls);
 	pdesolver->setParameters(inifile.sub("NewtonParameters"));

dune/dorie/model/richards/richards.hh

@@ -130,8 +130,10 @@ struct ModelRichardsTraits : public BaseTraits
 	/// Solution vector type
 	using U = typename IGO::Traits::Domain;
 	/// Linear solver types
-	using LS = Dune::PDELab::ISTLBackend_OVLP_AMG_4_DG<IGO,CC,LSGFS,LSCC,
-		Dune::PDELab::CG2DGProlongation,Dune::SeqSSOR,Dune::BiCGSTABSolver>;
+	using LS =
+		std::conditional_t<order == 0,
+			Dune::PDELab::ISTLBackend_BCGS_AMG_SSOR<IGO>,
+			Dune::PDELab::ISTLBackend_OVLP_AMG_4_DG<IGO, CC, LSGFS, LSCC, Dune::PDELab::CG2DGProlongation, Dune::SeqSSOR, Dune::BiCGSTABSolver>>;
 	/// Non-linear solver types
 	using PDESOLVER = Dune::PDELab::Newton<IGO,LS,U>;
 	/// Time stepping scheme

dune/dorie/model/transport/transport.cc

@@ -182,16 +182,33 @@ void ModelTransport<Traits>::operator_setup()
     explicit_igo = std::make_unique<ExplicitIGO>(*go0,*go1);
 
   // --- Solvers ---
+  // Initialize helper spaces for DG linear solver only
+  if constexpr (order > 0)
+  {
     lsgfs = std::make_unique<LSGFS>(LSGFSHelper::create(gv));
     lscc = std::make_unique<LSCC>();
+  }
+
   if (ts_param->implicit())
+  {
+    if constexpr (order == 0)
+      implicit_ls = std::make_unique<ImplicitLS>(*gfs, 1000, 0, true, true);
+    else
+    {
       implicit_ls = std::make_unique<ImplicitLS>(
-      *implicit_igo, *cc, *lsgfs, *lscc, 1000, 0, true, true
-    );
+        *implicit_igo, *cc, *lsgfs, *lscc, 1000, 0, true, true);
+    }
+  }
+  else
+  {
+    if constexpr (order == 0)
+      explicit_ls = std::make_unique<ExplicitLS>(*gfs, 1000, 0, true, true);
     else
+    {
       explicit_ls = std::make_unique<ExplicitLS>(
-      *explicit_igo, *cc, *lsgfs, *lscc, 1000, 0, true, true
-    );
+        *explicit_igo, *cc, *lsgfs, *lscc, 1000, 0, true, true);
+    }
+  }
 
   // --- Time Step Operators ---
   if (ts_param->implicit()){

dune/dorie/model/transport/transport.hh

@@ -158,14 +158,26 @@ struct ModelTransportTraits : public BaseTraits
   /// Solution vector type
   using U = typename ImplicitIGO::Traits::Domain;
   /// Linear solver type
-  using ImplicitLS = Dune::PDELab::ISTLBackend_OVLP_AMG_4_DG<
-    ImplicitIGO, CC, LSGFS, LSCC, Dune::PDELab::CG2DGProlongation,
-    Dune::SeqSSOR, Dune::BiCGSTABSolver
-  >;
-  using ExplicitLS = Dune::PDELab::ISTLBackend_OVLP_AMG_4_DG<
-    ExplicitIGO, CC, LSGFS, LSCC, Dune::PDELab::CG2DGProlongation,
-    Dune::SeqSSOR, Dune::BiCGSTABSolver
-  >;
+  using ImplicitLS = std::conditional_t<
+    order == 0,
+    Dune::PDELab::ISTLBackend_BCGS_AMG_SSOR<ImplicitIGO>,
+    Dune::PDELab::ISTLBackend_OVLP_AMG_4_DG<ImplicitIGO,
+                                            CC,
+                                            LSGFS,
+                                            LSCC,
+                                            Dune::PDELab::CG2DGProlongation,
+                                            Dune::SeqSSOR,
+                                            Dune::BiCGSTABSolver>>;
+  using ExplicitLS = std::conditional_t<
+    order == 0,
+    Dune::PDELab::ISTLBackend_BCGS_AMG_SSOR<ExplicitIGO>,
+    Dune::PDELab::ISTLBackend_OVLP_AMG_4_DG<ExplicitIGO,
+                                            CC,
+                                            LSGFS,
+                                            LSCC,
+                                            Dune::PDELab::CG2DGProlongation,
+                                            Dune::SeqSSOR,
+                                            Dune::BiCGSTABSolver>>;
   /// Methods computing the time step
   using ImplicitSLPS
     = Dune::PDELab::StationaryLinearProblemSolver<ImplicitIGO, ImplicitLS, U>;