Add a copy of ConvexProblem as MyConvexProblem

src/myblockproblem.hh

@@ -4,7 +4,7 @@
 #include <dune/common/bitsetvector.hh>
 
 #include <dune/tnnmg/problem-classes/bisection.hh>
-#include <dune/tnnmg/problem-classes/convexproblem.hh>
+//#include <dune/tnnmg/problem-classes/convexproblem.hh>
 #include <dune/tnnmg/problem-classes/nonlinearity.hh>
 #include <dune/tnnmg/problem-classes/onedconvexfunction.hh>
 
@@ -14,21 +14,21 @@
 
 /** \brief Base class for problems where each block can be solved with a scalar
  * Gauss-Seidel method */
-template <class ConvexProblemTypeTEMPLATE> class MyBlockProblem {
+template <class MyConvexProblemTypeTEMPLATE> class MyBlockProblem {
 public:
-  typedef ConvexProblemTypeTEMPLATE ConvexProblemType;
-  typedef typename ConvexProblemType::NonlinearityType NonlinearityType;
-  typedef typename ConvexProblemType::VectorType VectorType;
-  typedef typename ConvexProblemType::MatrixType MatrixType;
-  typedef typename ConvexProblemType::LocalVectorType LocalVectorType;
-  typedef typename ConvexProblemType::LocalMatrixType LocalMatrixType;
+  typedef MyConvexProblemTypeTEMPLATE MyConvexProblemType;
+  typedef typename MyConvexProblemType::NonlinearityType NonlinearityType;
+  typedef typename MyConvexProblemType::VectorType VectorType;
+  typedef typename MyConvexProblemType::MatrixType MatrixType;
+  typedef typename MyConvexProblemType::LocalVectorType LocalVectorType;
+  typedef typename MyConvexProblemType::LocalMatrixType LocalMatrixType;
 
-  static const int block_size = ConvexProblemType::block_size;
+  static const int block_size = MyConvexProblemType::block_size;
 
   /** \brief Solves one local system using a scalar Gauss-Seidel method */
   class IterateObject;
 
-  MyBlockProblem(ConvexProblemType& problem) : problem(problem) {
+  MyBlockProblem(MyConvexProblemType& problem) : problem(problem) {
     bisection = Bisection(0.0, 1.0, 1e-15, true, 1e-14);
   };
 
@@ -37,15 +37,15 @@ template <class ConvexProblemTypeTEMPLATE> class MyBlockProblem {
 
 private:
   // problem data
-  ConvexProblemType& problem;
+  MyConvexProblemType& problem;
 
   // commonly used minimization stuff
   Bisection bisection;
 };
 
 /** \brief Solves one local system using a scalar Gauss-Seidel method */
-template <class ConvexProblemTypeTEMPLATE>
-class MyBlockProblem<ConvexProblemTypeTEMPLATE>::IterateObject {
+template <class MyConvexProblemTypeTEMPLATE>
+class MyBlockProblem<MyConvexProblemTypeTEMPLATE>::IterateObject {
   friend class MyBlockProblem;
 
 protected:
@@ -54,7 +54,7 @@ class MyBlockProblem<ConvexProblemTypeTEMPLATE>::IterateObject {
    * \param problem The problem including quadratic part and nonlinear/nonsmooth
    * part
    */
-  IterateObject(const Bisection& bisection, ConvexProblemType& problem)
+  IterateObject(const Bisection& bisection, MyConvexProblemType& problem)
       : problem(problem),
         bisection(bisection),
         local_J(1.0, 0.0, problem.phi, 0, 0) {};
@@ -126,7 +126,7 @@ class MyBlockProblem<ConvexProblemTypeTEMPLATE>::IterateObject {
 
 private:
   // problem data
-  ConvexProblemType& problem;
+  MyConvexProblemType& problem;
 
   // commonly used minimization stuff
   Bisection bisection;

src/myconvexproblem.hh

+// Based on dune/tnnmg/problem-classes/convexproblem.hh
+
+#include <dune/istl/bcrsmatrix.hh>
+
+#include <dune/tnnmg/problem-classes/nonlinearity.hh>
+
+/** \brief General convex problem for a Truncated Nonsmooth Newton Multigrid
+   (TNNMG) solver
+
+    \tparam NonlinearityTypeTEMPLATE The type used to implement the nonlinearity
+    \tparam MatrixTypeTEMPLATE The type used for the matrix of the quadratic
+   part
+*/
+template <class NonlinearityTypeTEMPLATE = Nonlinearity<
+              Dune::FieldVector<double, 1>, Dune::FieldMatrix<double, 1, 1>>,
+          class MatrixTypeTEMPLATE =
+              Dune::BCRSMatrix<Dune::FieldMatrix<double, 1, 1>>>
+class MyConvexProblem {
+public:
+  typedef NonlinearityTypeTEMPLATE NonlinearityType;
+
+  typedef typename NonlinearityType::VectorType VectorType;
+  typedef MatrixTypeTEMPLATE MatrixType;
+  typedef typename NonlinearityType::LocalVectorType LocalVectorType;
+  typedef typename MatrixType::block_type LocalMatrixType;
+
+  static const int block_size = NonlinearityType::block_size;
+
+  /** \brief Constructor with the problem components
+
+      \param a A scalar factor in front of the quadratic part (the quadratic
+     part includes a factor of 1/2 already)
+      \param A The matrix of the quadratic part
+      \param am A scalar factor in front of the optional rank-one matrix
+      \param Am A rank-one matrix given by a single vector.  The matrix is
+     AmAm^T
+      \param phi The nonlinearity
+      \param f The linear functional
+      \param u The solution vector
+  */
+  MyConvexProblem(double a, const MatrixType& A, double am,
+                  const VectorType& Am, NonlinearityType& phi,
+                  const VectorType& f, VectorType& u)
+      : a(a), A(A), am(am), Am(Am), phi(phi), f(f), u(u) {};
+
+  double a;
+  const MatrixType& A;
+  double am;
+  const VectorType& Am;
+  NonlinearityType& phi;
+
+  const VectorType& f;
+
+  VectorType& u;
+};

src/one-body-sample.cc

@@ -28,13 +28,14 @@
 
 #include <dune/grid/common/mcmgmapper.hh>
 
-#include <dune/tnnmg/problem-classes/convexproblem.hh>
+//#include <dune/tnnmg/problem-classes/convexproblem.hh>
 #include <dune/tnnmg/problem-classes/nonlinearity.hh>
 #include <dune/tnnmg/nonlinearities/zerononlinearity.hh>
 
 #include <exception>
 #include <iostream>
 
+#include "myconvexproblem.hh"
 #include "myblockproblem.hh"
 
 int const dim = 2;
@@ -141,13 +142,14 @@ int main() {
       typedef ZeroNonlinearity<SmallVector, SmallMatrix> NonlinearityType;
       NonlinearityType phi;
 
-      typedef ConvexProblem<NonlinearityType, OperatorType> ConvexProblemType;
+      typedef MyConvexProblem<NonlinearityType, OperatorType>
+      MyConvexProblemType;
 
       VectorType unused_vector(grid.size(grid.maxLevel(), dim), 0);
-      ConvexProblemType myConvexProblem(1, stiffnessMatrix, 0, unused_vector,
-                                        phi, f, u1);
+      MyConvexProblemType myConvexProblem(1, stiffnessMatrix, 0, unused_vector,
+                                          phi, f, u1);
 
-      typedef MyBlockProblem<ConvexProblemType> MyBlockProblemType;
+      typedef MyBlockProblem<MyConvexProblemType> MyBlockProblemType;
       MyBlockProblemType myBlockProblem(myConvexProblem);
 
       GenericNonlinearGS<MyBlockProblemType> nonlinearGSStep;