Add hyper volume continuators & metrics handling feasibility constraint on objectives

moeo/src/continue/moeoHypContinue.h

@@ -46,11 +46,11 @@
 #include <archive/moeoUnboundedArchive.h>
 
 /**
-  Continues until the optimum ParetoSet level is reached.
+  Continues until the given ParetoSet level is reached.
 
   @ingroup Continuators
   */
-template< class MOEOT>
+template< class MOEOT, class MetricT = moeoHyperVolumeDifferenceMetric<typename MOEOT::ObjectiveVector> >
 class moeoHypContinue: public eoContinue<MOEOT>
 {
 public:
@@ -71,8 +71,20 @@ public:
         vectorToParetoSet(_OptimVec);
     }
 
+
     /** Returns false when a ParetoSet is reached. */
     virtual bool operator() ( const eoPop<MOEOT>& _pop )
+    {
+        std::vector<ObjectiveVector> bestCurrentParetoSet = pareto( arch );
+
+        return is_null_hypervolume( bestCurrentParetoSet );
+    }
+
+    virtual std::string className(void) const { return "moeoHypContinue"; }
+
+protected:
+
+    std::vector<ObjectiveVector> pareto( moeoArchive<MOEOT> & _archive )
     {
         std::vector < ObjectiveVector > bestCurrentParetoSet;
 
@@ -80,6 +92,11 @@ public:
             bestCurrentParetoSet.push_back(arch[i].objectiveVector());
         }
 
+        return bestCurrentParetoSet;
+    }
+
+    bool is_null_hypervolume( std::vector<ObjectiveVector>& bestCurrentParetoSet )
+    {
         double hypervolume= metric( bestCurrentParetoSet, OptimSet );
 
         if (hypervolume==0) {
@@ -91,7 +108,7 @@ public:
     }
 
     /** Translate a vector given as param to the ParetoSet that should be reached. */
-    void vectorToParetoSet(const std::vector<AtomType> & _OptimVec)
+    virtual void vectorToParetoSet(const std::vector<AtomType> & _OptimVec)
     {
         unsigned dim = (unsigned)(_OptimVec.size()/ObjectiveVector::Traits::nObjectives());
         OptimSet.resize(dim);
@@ -104,12 +121,98 @@ public:
         }
     }
 
-    virtual std::string className(void) const { return "moeoHypContinue"; }
+protected:
-
-private:
     moeoArchive <MOEOT> & arch;
-    moeoHyperVolumeDifferenceMetric <ObjectiveVector> metric;
+    MetricT metric;
     std::vector <ObjectiveVector> OptimSet;
 };
 
+
+/**
+  Continues until the (feasible or unfeasible) given Pareto set is reached.
+
+
+  @ingroup Continuators
+  */
+template< class MOEOT, class MetricT = moeoDualHyperVolumeDifferenceMetric<typename MOEOT::ObjectiveVector> >
+class moeoDualHypContinue: public moeoHypContinue<MOEOT, MetricT >
+{
+protected:
+    bool is_feasible;
+
+    using moeoHypContinue<MOEOT, MetricT>::arch;
+    using moeoHypContinue<MOEOT, MetricT>::OptimSet;
+
+public:
+    typedef typename MOEOT::ObjectiveVector ObjectiveVector;
+    typedef typename ObjectiveVector::Type AtomType;
+
+    /** A continuator that stops once a given Pareto front has been reached
+     *
+     * You should specify the feasibility of the targeted front.
+     * NOTE: the MOEOT::ObjectiveVector is supposed to implement the moeoDualRealObjectiveVector interface.
+     *
+     */
+    moeoDualHypContinue<MOEOT, MetricT>( const std::vector<AtomType> & _OptimVec, bool _is_feasible, moeoArchive < MOEOT > & _archive,  bool _normalize=true, double _rho=1.1 )
+        : moeoHypContinue<MOEOT, MetricT>( _OptimVec, _archive, _normalize, _rho ), is_feasible(_is_feasible)
+    {
+        assert( _OptimVec.size() > 0);
+        vectorToParetoSet(_OptimVec);
+    }
+
+    /** A continuator that stops once a given Pareto front has been reached
+     *
+     * You should specify the feasibility of the targeted front.
+     * NOTE: the MOEOT::ObjectiveVector is supposed to implement the moeoDualRealObjectiveVector interface.
+     *
+     */
+    moeoDualHypContinue<MOEOT, MetricT>( const std::vector<AtomType> & _OptimVec, bool _is_feasible, moeoArchive < MOEOT > & _archive,  bool _normalize=true, ObjectiveVector& _ref_point=NULL )
+        : moeoHypContinue<MOEOT, MetricT>( _OptimVec, _archive, _normalize, _ref_point ), is_feasible(_is_feasible)
+    {
+        assert( _OptimVec.size() > 0);
+        vectorToParetoSet(_OptimVec);
+    }
+
+    /** Returns false when a ParetoSet is reached. */
+    virtual bool operator() ( const eoPop<MOEOT>& _pop )
+    {
+        std::vector<ObjectiveVector> bestCurrentParetoSet = pareto( arch );
+
+#ifndef NDEBUG
+        assert( bestCurrentParetoSet.size() > 0 );
+        for( unsigned int i=1; i<bestCurrentParetoSet.size(); ++i ) {
+            assert( bestCurrentParetoSet[i].is_feasible() == bestCurrentParetoSet[0].is_feasible() );
+        }
+#endif
+
+        // The current Pareto front is either feasible or unfeasible.
+        // It could not contains both kind of objective vectors, because a feasible solution always dominates an unfeasible front.
+        if( bestCurrentParetoSet[0].is_feasible() != OptimSet[0].is_feasible() ) {
+            return false;
+        }
+
+        return is_null_hypervolume( bestCurrentParetoSet );
+    }
+
+protected:
+
+    using moeoHypContinue<MOEOT, MetricT>::pareto;
+    using moeoHypContinue<MOEOT, MetricT>::is_null_hypervolume;
+
+    /** Translate a vector given as param to the ParetoSet that should be reached. */
+    virtual void vectorToParetoSet(const std::vector<AtomType> & _OptimVec)
+    {
+        unsigned dim = (unsigned)(_OptimVec.size()/ObjectiveVector::Traits::nObjectives());
+        OptimSet.resize(dim);
+
+        unsigned k=0;
+        for(size_t i=0; i < dim; i++) {
+            for (size_t j=0; j < ObjectiveVector::Traits::nObjectives(); j++) {
+                // Use the feasibility declaration of an eoDualFitness
+                OptimSet[i][j] = AtomType(_OptimVec[k++], is_feasible);
+            }
+        }
+    }
+};
+
 #endif

moeo/src/fitness/moeoExpBinaryIndicatorBasedDualFitnessAssignment.h

@@ -10,6 +10,9 @@ protected:
 
 public:
     typedef typename MOEOT::ObjectiveVector ObjectiveVector;
+    typedef typename ObjectiveVector::Type Type;
+
+    using moeoExpBinaryIndicatorBasedFitnessAssignment<MOEOT>::values;
 
     moeoExpBinaryIndicatorBasedDualFitnessAssignment(
             moeoNormalizedSolutionVsSolutionBinaryMetric<ObjectiveVector,double> & metric,
@@ -54,11 +57,33 @@ public:
         this->setFitnesses(*ppop);
     }
 
+    /**
+     * Compute every indicator value in values (values[i] = I(_v[i], _o))
+     * @param _pop the population
+     */
+    void computeValues(const eoPop < MOEOT > & _pop)
+    {
+      values.clear();
+      values.resize(_pop.size());
+      for (unsigned int i=0; i<_pop.size(); i++)
+        {
+          values[i].resize(_pop.size());
+          // the metric may not be symetric, thus neither is the matrix
+          for (unsigned int j=0; j<_pop.size(); j++)
+            {
+              if (i != j)
+                {
+                  values[i][j] = Type( metric(_pop[i].objectiveVector(), _pop[j].objectiveVector()), _pop[i].objectiveVector().is_feasible() );
+                }
+            }
+        }
+    }
+
     virtual void setFitnesses(eoPop < MOEOT > & pop)
     {
         for (unsigned int i=0; i<pop.size(); i++) {
             // We should maintain the feasibility of the fitness across computations
-            pop[i].fitness( this->computeFitness(i), pop[i].is_feasible() );
+            pop[i].fitness( this->computeFitness(i), pop[i].fitness().is_feasible() );
         }
     }
 

moeo/src/fitness/moeoExpBinaryIndicatorBasedFitnessAssignment.h

@@ -187,7 +187,7 @@ class moeoExpBinaryIndicatorBasedFitnessAssignment : public moeoBinaryIndicatorB
             {
               if (i != j)
                 {
-                  values[i][j] = metric(_pop[i].objectiveVector(), _pop[j].objectiveVector());
+                  values[i][j] = Type( metric(_pop[i].objectiveVector(), _pop[j].objectiveVector()) );
                 }
             }
         }

moeo/src/metric/moeoHyperVolumeDifferenceMetric.h

@@ -86,6 +86,7 @@ class moeoHyperVolumeDifferenceMetric : public moeoVectorVsVectorBinaryMetric <
      */
     double operator()(const std::vector < ObjectiveVector > & _set1, const std::vector < ObjectiveVector > & _set2)
     {
+
         double hypervolume_set1;
         double hypervolume_set2;
 
@@ -132,7 +133,7 @@ class moeoHyperVolumeDifferenceMetric : public moeoVectorVsVectorBinaryMetric <
     }
 
     /**
-     * method caclulate bounds for the normalization
+     * method calculate bounds for the normalization
      * @param _set1 the vector contains all objective Vector of the first pareto front
      * @param _set2 the vector contains all objective Vector of the second pareto front
      */
@@ -182,7 +183,7 @@ class moeoHyperVolumeDifferenceMetric : public moeoVectorVsVectorBinaryMetric <
       return 1e-6;
     }
 
-    private:
+    protected:
 
     /*boolean indicates if data must be normalized or not*/
     bool normalize;
@@ -196,4 +197,90 @@ class moeoHyperVolumeDifferenceMetric : public moeoVectorVsVectorBinaryMetric <
 
   };
 
+
+template<class ObjectiveVector>
+class moeoDualHyperVolumeDifferenceMetric : public moeoHyperVolumeDifferenceMetric<ObjectiveVector>
+{
+protected:
+    using moeoHyperVolumeDifferenceMetric<ObjectiveVector>::rho;
+    using moeoHyperVolumeDifferenceMetric<ObjectiveVector>::normalize;
+    using moeoHyperVolumeDifferenceMetric<ObjectiveVector>::ref_point;
+    using moeoHyperVolumeDifferenceMetric<ObjectiveVector>::bounds;
+
+public:
+
+    typedef typename ObjectiveVector::Type Type;
+
+    moeoDualHyperVolumeDifferenceMetric( bool _normalize=true, double _rho=1.1)
+        : moeoHyperVolumeDifferenceMetric<ObjectiveVector>(_normalize, _rho)
+    {
+
+    }
+
+    moeoDualHyperVolumeDifferenceMetric( bool _normalize/*=true*/, ObjectiveVector& _ref_point/*=NULL*/ )
+        : moeoHyperVolumeDifferenceMetric<ObjectiveVector>( _normalize, _ref_point )
+    {
+
+    }
+
+    /**
+     * calculates and returns the HyperVolume value of a pareto front
+     * @param _set1 the vector contains all objective Vector of the first pareto front
+     * @param _set2 the vector contains all objective Vector of the second pareto front
+     */
+    double operator()(const std::vector < ObjectiveVector > & _set1, const std::vector < ObjectiveVector > & _set2)
+    {
+#ifndef NDEBUG
+        // the two sets must be homogeneous in feasibility
+        assert( _set1.size() > 0 );
+        for( unsigned int i=1; i<_set1.size(); ++i ) {
+            assert( _set1[i].is_feasible() == _set1[0].is_feasible() );
+        }
+        assert( _set2.size() > 0 );
+        for( unsigned int i=1; i<_set2.size(); ++i ) {
+            assert( _set2[i].is_feasible() == _set2[0].is_feasible() );
+        }
+        // and they must have the same feasibility
+        assert( _set1[0].is_feasible() == _set2[0].is_feasible() );
+#endif
+        bool feasible = _set1[0].is_feasible();
+
+        double hypervolume_set1;
+        double hypervolume_set2;
+
+        if(rho >= 1.0){
+            //determine bounds
+            setup(_set1, _set2);
+            //determine reference point
+            for (unsigned int i=0; i<ObjectiveVector::Traits::nObjectives(); i++){
+                if(normalize){
+                    if (ObjectiveVector::Traits::minimizing(i))
+                        ref_point[i]= Type(rho, feasible);
+                    else
+                        ref_point[i]= Type(1-rho, feasible);
+                }
+                else{
+                    if (ObjectiveVector::Traits::minimizing(i))
+                        ref_point[i]= Type(bounds[i].maximum() * rho, feasible);
+                    else
+                        ref_point[i]= Type(bounds[i].maximum() * (1-rho), feasible);
+                }
+            }
+            //if no normalization, reinit bounds to O..1 for
+            if(!normalize)
+                for (unsigned int i=0; i<ObjectiveVector::Traits::nObjectives(); i++)
+                    bounds[i] = eoRealInterval(0,1);
+
+        }
+        else if(normalize)
+            setup(_set1, _set2);
+
+        moeoHyperVolumeMetric <ObjectiveVector> unaryMetric(ref_point, bounds);
+        hypervolume_set1 = unaryMetric(_set1);
+        hypervolume_set2 = unaryMetric(_set2);
+
+        return hypervolume_set1 - hypervolume_set2;
+    }
+};
+
 #endif /*MOEOHYPERVOLUMEMETRIC_H_*/