From 82a3be3d74feeb47c2c11890be870819eaf6695c Mon Sep 17 00:00:00 2001
From: TimQu <tim.quatmann@rwth-aachen.de>
Date: Tue, 8 Nov 2016 22:20:41 +0100
Subject: [PATCH] .. missing files

Former-commit-id: f05bc337a5a514fa94f8481894f638770d76c48f
---
 .../pcaa/SparsePcaaAchievabilityQuery.cpp     | 101 ++++++++++++++++++
 .../pcaa/SparsePcaaAchievabilityQuery.h       |  58 ++++++++++
 .../pcaa/SparsePcaaParetoQuery.cpp            | 101 ++++++++++++++++++
 .../pcaa/SparsePcaaParetoQuery.h              |  47 ++++++++
 .../EndComponentEliminatorTest.cpp            |  14 +--
 5 files changed, 315 insertions(+), 6 deletions(-)
 create mode 100644 src/modelchecker/multiobjective/pcaa/SparsePcaaAchievabilityQuery.cpp
 create mode 100644 src/modelchecker/multiobjective/pcaa/SparsePcaaAchievabilityQuery.h
 create mode 100644 src/modelchecker/multiobjective/pcaa/SparsePcaaParetoQuery.cpp
 create mode 100644 src/modelchecker/multiobjective/pcaa/SparsePcaaParetoQuery.h

diff --git a/src/modelchecker/multiobjective/pcaa/SparsePcaaAchievabilityQuery.cpp b/src/modelchecker/multiobjective/pcaa/SparsePcaaAchievabilityQuery.cpp
new file mode 100644
index 000000000..a330208db
--- /dev/null
+++ b/src/modelchecker/multiobjective/pcaa/SparsePcaaAchievabilityQuery.cpp
@@ -0,0 +1,101 @@
+#include "src/modelchecker/multiobjective/pcaa/SparsePcaaAchievabilityQuery.h"
+
+#include "src/adapters/CarlAdapter.h"
+#include "src/models/sparse/Mdp.h"
+#include "src/models/sparse/MarkovAutomaton.h"
+#include "src/models/sparse/StandardRewardModel.h"
+#include "src/modelchecker/results/ExplicitQualitativeCheckResult.h"
+#include "src/utility/constants.h"
+#include "src/utility/vector.h"
+#include "src/settings//SettingsManager.h"
+#include "src/settings/modules/GeneralSettings.h"
+#include "src/settings/modules/MultiObjectiveSettings.h"
+
+namespace storm {
+    namespace modelchecker {
+        namespace multiobjective {
+            
+            
+            template <class SparseModelType, typename GeometryValueType>
+            SparsePcaaAchievabilityQuery<SparseModelType, GeometryValueType>::SparsePcaaAchievabilityQuery(SparsePcaaPreprocessorReturnType<SparseModelType>& preprocessorResult) : SparsePcaaQuery<SparseModelType, GeometryValueType>(preprocessorResult) {
+                STORM_LOG_ASSERT(preprocessorResult.queryType==SparsePcaaPreprocessorReturnType<SparseModelType>::QueryType::Achievability, "Invalid query Type");
+                initializeThresholdData();
+                
+                // Set the maximum gap between lower and upper bound of the weightVectorChecker result.
+                // This is the maximal edge length of the box we have to consider around each computed point
+                // We pick the gap such that the maximal distance between two points within this box is less than the given precision divided by two.
+                typename SparseModelType::ValueType gap = storm::utility::convertNumber<typename SparseModelType::ValueType>(storm::settings::getModule<storm::settings::modules::MultiObjectiveSettings>().getPrecision());
+                gap /= (storm::utility::one<typename SparseModelType::ValueType>() + storm::utility::one<typename SparseModelType::ValueType>());
+                gap /= storm::utility::sqrt(static_cast<typename SparseModelType::ValueType>(this->objectives.size()));
+                this->weightVectorChecker->setMaximumLowerUpperBoundGap(gap);
+                
+            }
+            
+            template <class SparseModelType, typename GeometryValueType>
+            void SparsePcaaAchievabilityQuery<SparseModelType, GeometryValueType>::initializeThresholdData() {
+                thresholds.reserve(this->objectives.size());
+                strictThresholds = storm::storage::BitVector(this->objectives.size(), false);
+                for(uint_fast64_t objIndex = 0; objIndex < this->objectives.size(); ++objIndex) {
+                    thresholds.push_back(storm::utility::convertNumber<GeometryValueType>(*this->objectives[objIndex].threshold));
+                    strictThresholds.set(objIndex, this->objectives[objIndex].thresholdIsStrict);
+                }
+            }
+            
+            template <class SparseModelType, typename GeometryValueType>
+            std::unique_ptr<CheckResult> SparsePcaaAchievabilityQuery<SparseModelType, GeometryValueType>::check() {
+               
+                bool result = this->checkAchievability();
+                
+                return std::unique_ptr<CheckResult>(new ExplicitQualitativeCheckResult(this->originalModel.getInitialStates().getNextSetIndex(0), result));
+                
+            }
+            
+            template <class SparseModelType, typename GeometryValueType>
+            bool SparsePcaaAchievabilityQuery<SparseModelType, GeometryValueType>::checkAchievability() {
+                // repeatedly refine the over/ under approximation until the threshold point is either in the under approx. or not in the over approx.
+                while(!this->maxStepsPerformed()){
+                    WeightVector separatingVector = this->findSeparatingVector(thresholds);
+                    this->performRefinementStep(std::move(separatingVector));
+                    if(!checkIfThresholdsAreSatisfied(this->overApproximation)){
+                        return false;
+                    }
+                    if(checkIfThresholdsAreSatisfied(this->underApproximation)){
+                        return true;
+                    }
+                }
+                STORM_LOG_ERROR("Could not check whether thresholds are achievable: Exceeded maximum number of refinement steps");
+                return false;
+            }
+
+            
+            template <class SparseModelType, typename GeometryValueType>
+            bool SparsePcaaAchievabilityQuery<SparseModelType, GeometryValueType>::checkIfThresholdsAreSatisfied(std::shared_ptr<storm::storage::geometry::Polytope<GeometryValueType>> const& polytope) {
+                std::vector<storm::storage::geometry::Halfspace<GeometryValueType>> halfspaces = polytope->getHalfspaces();
+                for(auto const& h : halfspaces) {
+                    GeometryValueType distance = h.distance(thresholds);
+                    if(distance < storm::utility::zero<GeometryValueType>()) {
+                        return false;
+                    }
+                    if(distance == storm::utility::zero<GeometryValueType>()) {
+                        // In this case, the thresholds point is on the boundary of the polytope.
+                        // Check if this is problematic for the strict thresholds
+                        for(auto strictThreshold : strictThresholds) {
+                            if(h.normalVector()[strictThreshold] > storm::utility::zero<GeometryValueType>()) {
+                                return false;
+                            }
+                        }
+                    }
+                }
+                return true;
+            }
+
+#ifdef STORM_HAVE_CARL
+            template class SparsePcaaAchievabilityQuery<storm::models::sparse::Mdp<double>, storm::RationalNumber>;
+            template class SparsePcaaAchievabilityQuery<storm::models::sparse::MarkovAutomaton<double>, storm::RationalNumber>;
+            
+            template class SparsePcaaAchievabilityQuery<storm::models::sparse::Mdp<storm::RationalNumber>, storm::RationalNumber>;
+         //   template class SparsePcaaAchievabilityQuery<storm::models::sparse::MarkovAutomaton<storm::RationalNumber>, storm::RationalNumber>;
+#endif
+        }
+    }
+}
diff --git a/src/modelchecker/multiobjective/pcaa/SparsePcaaAchievabilityQuery.h b/src/modelchecker/multiobjective/pcaa/SparsePcaaAchievabilityQuery.h
new file mode 100644
index 000000000..0b4ca9417
--- /dev/null
+++ b/src/modelchecker/multiobjective/pcaa/SparsePcaaAchievabilityQuery.h
@@ -0,0 +1,58 @@
+#ifndef STORM_MODELCHECKER_MULTIOBJECTIVE_PCAA_SPARSEPCAAACHIEVABILITYQUERY_H_
+#define STORM_MODELCHECKER_MULTIOBJECTIVE_PCAA_SPARSEPCAAACHIEVABILITYQUERY_H_
+
+#include "src/modelchecker/multiobjective/pcaa/SparsePcaaQuery.h"
+
+namespace storm {
+    namespace modelchecker {
+        namespace multiobjective {
+            
+            /*
+             * This class represents a query for the Pareto curve approximation algorithm (Pcaa).
+             * It implements the necessary computations for the different query types.
+             */
+            template <class SparseModelType, typename GeometryValueType>
+            class SparsePcaaAchievabilityQuery : public SparsePcaaQuery<SparseModelType, GeometryValueType> {
+            public:
+                
+                // Typedefs for simple geometric objects
+                typedef std::vector<GeometryValueType> Point;
+                typedef std::vector<GeometryValueType> WeightVector;
+                
+                /*
+                 * Creates a new query for the Pareto curve approximation algorithm (Pcaa)
+                 * @param preprocessorResult the result from preprocessing
+                 */
+                SparsePcaaAchievabilityQuery(SparsePcaaPreprocessorReturnType<SparseModelType>& preprocessorResult);
+                
+                
+                /*
+                 * Invokes the computation and retrieves the result
+                 */
+                virtual std::unique_ptr<CheckResult> check() override;
+                
+            private:
+                
+                void initializeThresholdData();
+                
+                /* 
+                 * Returns whether the given thresholds are achievable.
+                 */
+                bool checkAchievability();
+                
+                /*
+                 * Returns true iff there is one point in the given polytope that satisfies the given thresholds.
+                 * It is assumed that the given polytope contains the downward closure of its vertices.
+                 */
+                bool checkIfThresholdsAreSatisfied(std::shared_ptr<storm::storage::geometry::Polytope<GeometryValueType>> const& polytope);
+
+                
+                Point thresholds;
+                storm::storage::BitVector strictThresholds;
+            };
+            
+        }
+    }
+}
+
+#endif /* STORM_MODELCHECKER_MULTIOBJECTIVE_PCAA_SPARSEPCAAACHIEVABILITYQUERY_H_ */
diff --git a/src/modelchecker/multiobjective/pcaa/SparsePcaaParetoQuery.cpp b/src/modelchecker/multiobjective/pcaa/SparsePcaaParetoQuery.cpp
new file mode 100644
index 000000000..14605ad28
--- /dev/null
+++ b/src/modelchecker/multiobjective/pcaa/SparsePcaaParetoQuery.cpp
@@ -0,0 +1,101 @@
+#include "src/modelchecker/multiobjective/pcaa/SparsePcaaParetoQuery.h"
+
+#include "src/adapters/CarlAdapter.h"
+#include "src/models/sparse/Mdp.h"
+#include "src/models/sparse/MarkovAutomaton.h"
+#include "src/models/sparse/StandardRewardModel.h"
+#include "src/modelchecker/results/ParetoCurveCheckResult.h"
+#include "src/utility/constants.h"
+#include "src/utility/vector.h"
+#include "src/settings//SettingsManager.h"
+#include "src/settings/modules/MultiObjectiveSettings.h"
+#include "src/settings/modules/GeneralSettings.h"
+
+
+namespace storm {
+    namespace modelchecker {
+        namespace multiobjective {
+            
+            template <class SparseModelType, typename GeometryValueType>
+            SparsePcaaParetoQuery<SparseModelType, GeometryValueType>::SparsePcaaParetoQuery(SparsePcaaPreprocessorReturnType<SparseModelType>& preprocessorResult) : SparsePcaaQuery<SparseModelType, GeometryValueType>(preprocessorResult) {
+                STORM_LOG_ASSERT(preprocessorResult.queryType==SparsePcaaPreprocessorReturnType<SparseModelType>::QueryType::Pareto, "Invalid query Type");
+                
+                // Set the maximum gap between lower and upper bound of the weightVectorChecker result.
+                // This is the maximal edge length of the box we have to consider around each computed point
+                // We pick the gap such that the maximal distance between two points within this box is less than the given precision divided by two.
+                typename SparseModelType::ValueType gap = storm::utility::convertNumber<typename SparseModelType::ValueType>(storm::settings::getModule<storm::settings::modules::MultiObjectiveSettings>().getPrecision());
+                gap /= (storm::utility::one<typename SparseModelType::ValueType>() + storm::utility::one<typename SparseModelType::ValueType>());
+                gap /= storm::utility::sqrt(static_cast<typename SparseModelType::ValueType>(this->objectives.size()));
+                this->weightVectorChecker->setMaximumLowerUpperBoundGap(gap);
+                
+            }
+            
+            template <class SparseModelType, typename GeometryValueType>
+            std::unique_ptr<CheckResult> SparsePcaaParetoQuery<SparseModelType, GeometryValueType>::check() {
+                
+                // refine the approximation
+                exploreSetOfAchievablePoints();
+                
+                // obtain the data for the checkresult
+                std::vector<std::vector<typename SparseModelType::ValueType>> paretoOptimalPoints;
+                paretoOptimalPoints.reserve(this->refinementSteps.size());
+                for(auto const& step : this->refinementSteps) {
+                    paretoOptimalPoints.push_back(storm::utility::vector::convertNumericVector<typename SparseModelType::ValueType>(this->transformPointToOriginalModel(step.lowerBoundPoint)));
+                }
+                return std::unique_ptr<CheckResult>(new ParetoCurveCheckResult<typename SparseModelType::ValueType>(this->originalModel.getInitialStates().getNextSetIndex(0),
+                                                                        std::move(paretoOptimalPoints),
+                                                                        this->transformPolytopeToOriginalModel(this->underApproximation)->template convertNumberRepresentation<typename SparseModelType::ValueType>(),
+                                                                        this->transformPolytopeToOriginalModel(this->overApproximation)->template convertNumberRepresentation<typename SparseModelType::ValueType>()));
+                
+                
+            }
+            
+            template <class SparseModelType, typename GeometryValueType>
+            void SparsePcaaParetoQuery<SparseModelType, GeometryValueType>::exploreSetOfAchievablePoints() {
+            
+                //First consider the objectives individually
+                for(uint_fast64_t objIndex = 0; objIndex<this->objectives.size() && !this->maxStepsPerformed(); ++objIndex) {
+                    WeightVector direction(this->objectives.size(), storm::utility::zero<GeometryValueType>());
+                    direction[objIndex] = storm::utility::one<GeometryValueType>();
+                    this->performRefinementStep(std::move(direction));
+                }
+                
+                while(!this->maxStepsPerformed()) {
+                    // Get the halfspace of the underApproximation with maximal distance to a vertex of the overApproximation
+                    std::vector<storm::storage::geometry::Halfspace<GeometryValueType>> underApproxHalfspaces = this->underApproximation->getHalfspaces();
+                    std::vector<Point> overApproxVertices = this->overApproximation->getVertices();
+                    uint_fast64_t farestHalfspaceIndex = underApproxHalfspaces.size();
+                    GeometryValueType farestDistance = storm::utility::zero<GeometryValueType>();
+                    for(uint_fast64_t halfspaceIndex = 0; halfspaceIndex < underApproxHalfspaces.size(); ++halfspaceIndex) {
+                        for(auto const& vertex : overApproxVertices) {
+                            GeometryValueType distance = -underApproxHalfspaces[halfspaceIndex].euclideanDistance(vertex);
+                            if(distance > farestDistance) {
+                                farestHalfspaceIndex = halfspaceIndex;
+                                farestDistance = distance;
+                            }
+                        }
+                    }
+                    if(farestDistance < storm::utility::convertNumber<GeometryValueType>(storm::settings::getModule<storm::settings::modules::MultiObjectiveSettings>().getPrecision())) {
+                        // Goal precision reached!
+                        return;
+                    }
+                    STORM_LOG_DEBUG("Current precision of the approximation of the pareto curve is ~" << storm::utility::convertNumber<double>(farestDistance));
+                    WeightVector direction = underApproxHalfspaces[farestHalfspaceIndex].normalVector();
+                    this->performRefinementStep(std::move(direction));
+                }
+                STORM_LOG_ERROR("Could not reach the desired precision: Exceeded maximum number of refinement steps");
+            }
+            
+            
+
+            
+#ifdef STORM_HAVE_CARL
+            template class SparsePcaaParetoQuery<storm::models::sparse::Mdp<double>, storm::RationalNumber>;
+            template class SparsePcaaParetoQuery<storm::models::sparse::MarkovAutomaton<double>, storm::RationalNumber>;
+            
+            template class SparsePcaaParetoQuery<storm::models::sparse::Mdp<storm::RationalNumber>, storm::RationalNumber>;
+          //  template class SparsePcaaParetoQuery<storm::models::sparse::MarkovAutomaton<storm::RationalNumber>, storm::RationalNumber>;
+#endif
+        }
+    }
+}
diff --git a/src/modelchecker/multiobjective/pcaa/SparsePcaaParetoQuery.h b/src/modelchecker/multiobjective/pcaa/SparsePcaaParetoQuery.h
new file mode 100644
index 000000000..bb47d5158
--- /dev/null
+++ b/src/modelchecker/multiobjective/pcaa/SparsePcaaParetoQuery.h
@@ -0,0 +1,47 @@
+#ifndef STORM_MODELCHECKER_MULTIOBJECTIVE_PCAA_SPARSEPCAAPARETOQUERY_H_
+#define STORM_MODELCHECKER_MULTIOBJECTIVE_PCAA_SPARSEPCAAPARETOQUERY_H_
+
+#include "src/modelchecker/multiobjective/pcaa/SparsePcaaQuery.h"
+
+namespace storm {
+    namespace modelchecker {
+        namespace multiobjective {
+            
+            /*
+             * This class represents a query for the Pareto curve approximation algorithm (Pcaa).
+             * It implements the necessary computations for the different query types.
+             */
+            template <class SparseModelType, typename GeometryValueType>
+            class SparsePcaaParetoQuery : public SparsePcaaQuery<SparseModelType, GeometryValueType> {
+            public:
+                
+                // Typedefs for simple geometric objects
+                typedef std::vector<GeometryValueType> Point;
+                typedef std::vector<GeometryValueType> WeightVector;
+                
+                /*
+                 * Creates a new query for the Pareto curve approximation algorithm (Pcaa)
+                 * @param preprocessorResult the result from preprocessing
+                 */
+                SparsePcaaParetoQuery(SparsePcaaPreprocessorReturnType<SparseModelType>& preprocessorResult);
+                
+                
+                /*
+                 * Invokes the computation and retrieves the result
+                 */
+                virtual std::unique_ptr<CheckResult> check() override;
+                
+            private:
+                
+                
+                /* 
+                 * Performs refinement steps until the approximation is sufficiently precise
+                 */
+                void exploreSetOfAchievablePoints();
+            };
+            
+        }
+    }
+}
+
+#endif /* STORM_MODELCHECKER_MULTIOBJECTIVE_PCAA_SPARSEPCAAPARETOQUERY_H_ */
diff --git a/test/functional/transformer/EndComponentEliminatorTest.cpp b/test/functional/transformer/EndComponentEliminatorTest.cpp
index 8785df336..2210eafc8 100644
--- a/test/functional/transformer/EndComponentEliminatorTest.cpp
+++ b/test/functional/transformer/EndComponentEliminatorTest.cpp
@@ -35,17 +35,19 @@ TEST(NeutralECRemover, SimpleModelTest) {
     ASSERT_NO_THROW(matrix = builder.build());
     
     storm::storage::BitVector possibleEcRows(12, true);
-    consideredRows.set(3, false);
-    consideredRows.set(6, false);
-    consideredRows.set(7, false);
-    consideredRows.set(8, false);
-    consideredRows.set(11, false);
+    possibleEcRows.set(3, false);
+    possibleEcRows.set(6, false);
+    possibleEcRows.set(7, false);
+    possibleEcRows.set(8, false);
+    possibleEcRows.set(11, false);
     
     storm::storage::BitVector allowEmptyRows(5, true);
     allowEmptyRows.set(1, false);
     allowEmptyRows.set(4, false);
+    storm::storage::BitVector subsystem(5, true);
+    allowEmptyRows.set(2, false);
     
-    auto res = storm::transformer::EndComponentEliminator<double>::transform(matrix, possibleEcRows, allowEmptyRows);
+    auto res = storm::transformer::EndComponentEliminator<double>::transform(matrix, subsystem, possibleEcRows, allowEmptyRows);
     
     // Expected data
     storm::storage::SparseMatrixBuilder<double> expectedBuilder(8, 3, 8, true, true, 3);