Expose box constraints as intervals

dune/solvers/common/boxconstraint.hh

@@ -8,63 +8,67 @@
 
 #include <dune/common/fvector.hh>
 
+#include "interval.hh"
+
 template <class T, int dim>
 class BoxConstraint {
 public:
+    using Interval = typename Dune::Solvers::Interval<T>;
 
     /** \brief Default constructor, makes the obstacle as large as possible */
     BoxConstraint() {
-        for (int i=0; i<dim; i++) {
-            val[2*i]   = -std::numeric_limits<T>::max();
-            val[2*i+1] =  std::numeric_limits<T>::max();
-        }
+        clear();
     }
 
     /** \brief Set obstacle values to +/- std::numeric_limits<T>::max() */
     void clear() {
-        for (int i=0; i<dim; i++) {
-            val[2*i]   = -std::numeric_limits<T>::max();
-            val[2*i+1] =  std::numeric_limits<T>::max();
-        }
+        for (int i=0; i<dim; i++)
+            val[i] = { -std::numeric_limits<T>::max(),
+                       std::numeric_limits<T>::max() };
     }
 
     //! Subtract vector from box
     BoxConstraint<T,dim>& operator-= (const Dune::FieldVector<T, dim>& v)
     {
         for (int i=0; i<dim; i++) {
-            val[2*i]   -= v[i];
-            val[2*i+1] -= v[i];
+            val[i][0] -= v[i];
+            val[i][1] -= v[i];
         }
         return *this;
     }
 
+    Interval const &operator[](size_t i) const {
+        return val[i];
+    }
+    Interval &operator[](size_t i) {
+        return val[i];
+    }
+
     //! Access to the lower obstacles
-    T& lower(size_t i) {return val[2*i];}
+    T& lower(size_t i) {return val[i][0];}
 
     //! Const access to the lower obstacles
-    const T& lower(size_t i) const {return val[2*i];}
+    const T& lower(size_t i) const {return val[i][0];}
 
     //! Access to the upper obstacles
-    T& upper(size_t i) {return val[2*i+1];}
+    T& upper(size_t i) {return val[i][1];}
 
     //! Const access to the upper obstacles
-    const T& upper(size_t i) const {return val[2*i+1];}
+    const T& upper(size_t i) const {return val[i][1];}
 
     //! Send box constraint to an output stream
     friend
     std::ostream& operator<< (std::ostream& s, const BoxConstraint<T,dim>& v)
     {
         for (int i=0; i<dim; i++)
-            s << "Direction: " << i <<",  val[0] = " << v.val[2*i]
-              << ", val[1] = " << v.val[2*i+1] << std::endl;
+            s << "Direction: " << i <<",  val[0] = " << v[i][0]
+              << ", val[1] = " << v[i][1] << std::endl;
 
         return s;
     }
 
 protected:
-
-    std::array<T, 2*dim> val;
-
+    std::array<Interval, dim> val;
 };
 
 #endif

dune/solvers/common/interval.hh

@@ -24,8 +24,9 @@ public:
 
     /** \brief Construct from an initializer list */
     Interval(std::initializer_list<field_type> const &input)
-    : data_(input)
-    {}
+    {
+        std::copy(input.begin(), input.end(), data_.begin());
+    }
 
     /** \brief Array-like access
      */

dune/solvers/iterationsteps/obstacletnnmgstep.hh

@@ -265,12 +265,8 @@ class ObstacleTNNMGStep
 
                 // project correction into local defect constraints
                 for (int j = 0; j < blockSize; ++j)
-                {
-                    if (coarseCorrection_[i][j] < defectConstraint.lower(j))
-                        coarseCorrection_[i][j] = defectConstraint.lower(j);
-                    if (coarseCorrection_[i][j] > defectConstraint.upper(j))
-                        coarseCorrection_[i][j] = defectConstraint.upper(j);
-                }
+                    coarseCorrection_[i][j]
+                        = defectConstraint[j].projectIn(coarseCorrection_[i][j]);
 
                 for (int j = 0; j < blockSize; ++j)
                 {
@@ -299,8 +295,7 @@ class ObstacleTNNMGStep
                 damping = (residual_*coarseCorrection_) / coarseCorrectionNorm2;
             else
                 damping = 0;
-            damping = std::max(damping, directionalConstraints.lower(0));
-            damping = std::min(damping, directionalConstraints.upper(0));
+            damping = directionalConstraints[0].projectIn(damping);
 
             // apply correction
             x.axpy(damping, coarseCorrection_);