Modified the multiplication of a sparse matrix with a dense vector to only use the parallel version if available and the number of nonzero entries exceeds a certain threshold.

Former-commit-id: 07957b2d29

src/storage/SparseMatrix.cpp

@@ -770,6 +770,36 @@ namespace storm {
         template<typename T>
         void SparseMatrix<T>::multiplyWithVector(std::vector<T> const& vector, std::vector<T>& result) const {
 #ifdef STORM_HAVE_INTELTBB
+            if (this->getNonzeroEntryCount() > 10000) {
+                return this->multiplyWithVectorParallel(vector, result);
+            } else {
+                return this->multiplyWithVectorSequential(vector, result);
+            }
+#else
+            return multiplyWithVectorSequential(vector, result);
+#endif
+        }
+        
+        template<typename T>
+        void SparseMatrix<T>::multiplyWithVectorSequential(std::vector<T> const& vector, std::vector<T>& result) const {
+            const_iterator it = this->begin();
+            const_iterator ite;
+            typename std::vector<uint_fast64_t>::const_iterator rowIterator = rowIndications.begin();
+            typename std::vector<T>::iterator resultIterator = result.begin();
+            typename std::vector<T>::iterator resultIteratorEnd = result.end();
+            
+            for (; resultIterator != resultIteratorEnd; ++rowIterator, ++resultIterator) {
+                *resultIterator = storm::utility::constantZero<T>();
+                
+                for (ite = this->begin() + *(rowIterator + 1); it != ite; ++it) {
+                    *resultIterator += it->getValue() * vector[it->getColumn()];
+                }
+            }
+        }
+
+#ifdef STORM_HAVE_INTELTBB
+        template<typename T>
+        void SparseMatrix<T>::multiplyWithVectorParallel(std::vector<T> const& vector, std::vector<T>& result) const {
             tbb::parallel_for(tbb::blocked_range<uint_fast64_t>(0, result.size(), 10),
                               [&] (tbb::blocked_range<uint_fast64_t> const& range) {
                                   uint_fast64_t startRow = range.begin();
@@ -778,8 +808,8 @@ namespace storm {
                                   const_iterator ite;
                                   std::vector<uint_fast64_t>::const_iterator rowIterator = this->rowIndications.begin() + startRow;
                                   std::vector<uint_fast64_t>::const_iterator rowIteratorEnd = this->rowIndications.begin() + endRow;
-                                  std::vector<T>::iterator resultIterator = result.begin() + startRow;
-                                  std::vector<T>::iterator resultIteratorEnd = result.begin() + endRow;
+                                  typename std::vector<T>::iterator resultIterator = result.begin() + startRow;
+                                  typename std::vector<T>::iterator resultIteratorEnd = result.begin() + endRow;
                                   
                                   for (; resultIterator != resultIteratorEnd; ++rowIterator, ++resultIterator) {
                                       *resultIterator = storm::utility::constantZero<T>();
@@ -789,23 +819,9 @@ namespace storm {
                                       }
                                   }
                               });
-#else
-            const_iterator it = this->begin();
-            const_iterator ite;
-            typename std::vector<uint_fast64_t>::const_iterator rowIterator = rowIndications.begin();
-            typename std::vector<T>::iterator resultIterator = result.begin();
-            typename std::vector<T>::iterator resultIteratorEnd = result.end();
-            
-            for (; resultIterator != resultIteratorEnd; ++rowIterator, ++resultIterator) {
-                *resultIterator = storm::utility::constantZero<T>();
-                
-                for (ite = this->begin() + *(rowIterator + 1); it != ite; ++it) {
-                    *resultIterator += it->getValue() * vector[it->getColumn()];
-                }
-            }
-#endif
         }
-        
+#endif
+
         template<typename T>
         uint_fast64_t SparseMatrix<T>::getSizeInMemory() const {
             uint_fast64_t size = sizeof(*this);

src/storage/SparseMatrix.h

@@ -576,7 +576,9 @@ namespace storm {
             std::vector<T> getPointwiseProductRowSumVector(storm::storage::SparseMatrix<T> const& otherMatrix) const;
             
             /*!
-             * Multiplies the matrix with the given vector and writes the result to given result vector.
+             * Multiplies the matrix with the given vector and writes the result to the given result vector. If a
+             * parallel implementation is available and it is considered worthwhile (heuristically, based on the metrics
+             * of the matrix), the multiplication is carried out in parallel.
              *
              * @param vector The vector with which to multiply the matrix.
              * @param result The vector that is supposed to hold the result of the multiplication after the operation.
@@ -584,6 +586,28 @@ namespace storm {
              */
             void multiplyWithVector(std::vector<T> const& vector, std::vector<T>& result) const;
             
+            /*!
+             * Multiplies the matrix with the given vector in a sequential way and writes the result to the given result
+             * vector.
+             *
+             * @param vector The vector with which to multiply the matrix.
+             * @param result The vector that is supposed to hold the result of the multiplication after the operation.
+             * @return The product of the matrix and the given vector as the content of the given result vector.
+             */
+            void multiplyWithVectorSequential(std::vector<T> const& vector, std::vector<T>& result) const;
+
+#ifdef STORM_HAVE_INTELTBB
+            /*!
+             * Multiplies the matrix with the given vector in a parallel fashion using Intel's TBB and writes the result
+             * to the given result vector.
+             *
+             * @param vector The vector with which to multiply the matrix.
+             * @param result The vector that is supposed to hold the result of the multiplication after the operation.
+             * @return The product of the matrix and the given vector as the content of the given result vector.
+             */
+            void multiplyWithVectorParallel(std::vector<T> const& vector, std::vector<T>& result) const;
+#endif
+            
             /*!
              * Computes the sum of the entries in a given row.
              *