diff --git a/src/adapters/GmmxxAdapter.h b/src/adapters/GmmxxAdapter.h
index 43a3f8302..748056e4b 100644
--- a/src/adapters/GmmxxAdapter.h
+++ b/src/adapters/GmmxxAdapter.h
@@ -27,7 +27,7 @@ public:
*/
template<class T>
static gmm::csr_matrix<T>* toGmmxxSparseMatrix(storm::storage::SquareSparseMatrix<T> const& matrix) {
- uint_fast64_t realNonZeros = matrix.getNonZeroEntryCount() + matrix.getDiagonalNonZeroEntryCount();
+ uint_fast64_t realNonZeros = matrix.getNonZeroEntryCount();
LOG4CPLUS_DEBUG(logger, "Converting matrix with " << realNonZeros << " non-zeros to gmm++ format.");
// Prepare the resulting matrix.
@@ -44,48 +44,23 @@ public:
T* tmpValueArray = new T[realNonZeros];
T zero(0);
uint_fast64_t currentPosition = 0;
- uint_fast64_t insertedDiagonalElements = 0;
for (uint_fast64_t i = 0; i < matrix.rowCount; ++i) {
// Compute correct start index of row.
- result->jc[i] = matrix.rowIndications[i] + insertedDiagonalElements;
-
- // If the current row has no non-zero which is not on the diagonal, we have to check the
- // diagonal element explicitly.
- if (matrix.rowIndications[i + 1] - matrix.rowIndications[i] == 0) {
- if (matrix.diagonalStorage[i] != zero) {
- tmpColumnIndicationsArray[currentPosition] = i;
- tmpValueArray[currentPosition] = matrix.diagonalStorage[i];
- ++currentPosition; ++insertedDiagonalElements;
- }
- } else {
- // Otherwise, we can just enumerate the non-zeros which are not on the diagonal
- // and fit in the diagonal element where appropriate.
- bool includedDiagonal = false;
- for (uint_fast64_t j = matrix.rowIndications[i]; j < matrix.rowIndications[i + 1]; ++j) {
- if (matrix.diagonalStorage[i] != zero && !includedDiagonal && matrix.columnIndications[j] > i) {
- includedDiagonal = true;
- tmpColumnIndicationsArray[currentPosition] = i;
- tmpValueArray[currentPosition] = matrix.diagonalStorage[i];
- ++currentPosition; ++insertedDiagonalElements;
- }
- tmpColumnIndicationsArray[currentPosition] = matrix.columnIndications[j];
- tmpValueArray[currentPosition] = matrix.valueStorage[j];
- ++currentPosition;
- }
-
- // If the diagonal element is non-zero and was not inserted until now (i.e. all
- // off-diagonal elements in the row are before the diagonal element.
- if (!includedDiagonal && matrix.diagonalStorage[i] != zero) {
- tmpColumnIndicationsArray[currentPosition] = i;
- tmpValueArray[currentPosition] = matrix.diagonalStorage[i];
- ++currentPosition; ++insertedDiagonalElements;
- }
+ result->jc[i] = matrix.rowIndications[i];
+
+ // Otherwise, we can just enumerate the non-zeros which are not on the diagonal
+ // and fit in the diagonal element where appropriate.
+ for (uint_fast64_t j = matrix.rowIndications[i]; j < matrix.rowIndications[i + 1]; ++j) {
+ tmpColumnIndicationsArray[currentPosition] = matrix.columnIndications[j];
+ tmpValueArray[currentPosition] = matrix.valueStorage[j];
+ ++currentPosition;
}
}
// Fill in sentinel element at the end.
result->jc[matrix.rowCount] = static_cast<unsigned int>(realNonZeros);
// Now, we can copy the temporary array to the GMMXX format.
+ // FIXME Now everything can just be copied. No TMP needed anymore
result->ir.resize(realNonZeros);
std::copy(tmpColumnIndicationsArray, tmpColumnIndicationsArray + realNonZeros, result->ir.begin());
delete[] tmpColumnIndicationsArray;
diff --git a/src/formula/BoundedEventually.h b/src/formula/BoundedEventually.h
index 58973456d..0a66f6b83 100644
--- a/src/formula/BoundedEventually.h
+++ b/src/formula/BoundedEventually.h
@@ -120,7 +120,7 @@ public:
BoundedEventually<T>* result = new BoundedEventually<T>();
result->setBound(bound);
if (child != nullptr) {
- result->setRight(child->clone());
+ result->setChild(child->clone());
}
return result;
}
diff --git a/src/models/GraphTransitions.h b/src/models/GraphTransitions.h
index 3ccc6bed5..122f20dfd 100644
--- a/src/models/GraphTransitions.h
+++ b/src/models/GraphTransitions.h
@@ -93,12 +93,12 @@ private:
// First, we copy the index list from the sparse matrix as this will
// stay the same.
- std::copy(transitionMatrix->getRowIndicationsPointer(), transitionMatrix->getRowIndicationsPointer() + numberOfStates + 1, this->stateIndications);
+ std::copy(transitionMatrix->getRowIndicationsPointer().begin(), transitionMatrix->getRowIndicationsPointer().end(), this->stateIndications);
// Now we can iterate over all rows of the transition matrix and record
// the target state.
for (uint_fast64_t i = 0, currentNonZeroElement = 0; i < numberOfStates; i++) {
- for (auto rowIt = transitionMatrix->beginConstColumnNoDiagIterator(i); rowIt != transitionMatrix->endConstColumnNoDiagIterator(i); ++rowIt) {
+ for (auto rowIt = transitionMatrix->beginConstColumnIterator(i); rowIt != transitionMatrix->endConstColumnIterator(i); ++rowIt) {
this->stateIndications[currentNonZeroElement++] = *rowIt;
}
}
@@ -117,7 +117,7 @@ private:
// NOTE: We disregard the diagonal here, as we only consider "true"
// predecessors.
for (uint_fast64_t i = 0; i < numberOfStates; i++) {
- for (auto rowIt = transitionMatrix->beginConstColumnNoDiagIterator(i); rowIt != transitionMatrix->endConstColumnNoDiagIterator(i); ++rowIt) {
+ for (auto rowIt = transitionMatrix->beginConstColumnIterator(i); rowIt != transitionMatrix->endConstColumnIterator(i); ++rowIt) {
this->stateIndications[*rowIt + 1]++;
}
}
@@ -140,7 +140,7 @@ private:
// Now we are ready to actually fill in the list of predecessors for
// every state. Again, we start by considering all but the last row.
for (uint_fast64_t i = 0; i < numberOfStates; i++) {
- for (auto rowIt = transitionMatrix->beginConstColumnNoDiagIterator(i); rowIt != transitionMatrix->endConstColumnNoDiagIterator(i); ++rowIt) {
+ for (auto rowIt = transitionMatrix->beginConstColumnIterator(i); rowIt != transitionMatrix->endConstColumnIterator(i); ++rowIt) {
this->successorList[nextIndicesList[*rowIt]++] = i;
}
}
diff --git a/src/parser/DeterministicSparseTransitionParser.cpp b/src/parser/DeterministicSparseTransitionParser.cpp
index f3d3b9d00..9c2e5c1b4 100644
--- a/src/parser/DeterministicSparseTransitionParser.cpp
+++ b/src/parser/DeterministicSparseTransitionParser.cpp
@@ -32,13 +32,10 @@ namespace parser{
* non-zero cells and maximum node id.
*
* This method does the first pass through the .tra file and computes
- * the number of non-zero elements that are not diagonal elements,
- * which correspondents to the number of transitions that are not
- * self-loops.
- * (Diagonal elements are treated in a special way).
+ * the number of non-zero elements.
* It also calculates the maximum node id and stores it in maxnode.
*
- * @return The number of non-zero elements that are not on the diagonal
+ * @return The number of non-zero elements
* @param buf Data to scan. Is expected to be some char array.
* @param maxnode Is set to highest id of all nodes.
*/
@@ -89,7 +86,8 @@ uint_fast64_t DeterministicSparseTransitionParser::firstPass(char* buf, uint_fas
LOG4CPLUS_ERROR(logger, "Expected a positive probability but got \"" << std::string(buf, 0, 16) << "\".");
return 0;
}
- if (row == col) non_zero--;
+ // not needed anymore
+ //if (row == col) non_zero--;
buf = trimWhitespaces(tmp);
}
diff --git a/src/storage/SquareSparseMatrix.h b/src/storage/SquareSparseMatrix.h
index ae48cf4ed..51d0f524e 100644
--- a/src/storage/SquareSparseMatrix.h
+++ b/src/storage/SquareSparseMatrix.h
@@ -5,6 +5,7 @@
#include <new>
#include <algorithm>
#include <iostream>
+#include <iterator>
#include "boost/integer/integer_mask.hpp"
#include "src/exceptions/InvalidStateException.h"
@@ -31,8 +32,7 @@ namespace storm {
namespace storage {
/*!
- * A sparse matrix class with a constant number of non-zero entries on the non-diagonal fields
- * and a separate dense storage for the diagonal elements.
+ * A sparse matrix class with a constant number of non-zero entries.
* NOTE: Addressing *is* zero-based, so the valid range for getValue and addNextValue is 0..(rows - 1)
* where rows is the first argument to the constructor.
*/
@@ -73,9 +73,25 @@ public:
* Constructs a sparse matrix object with the given number of rows.
* @param rows The number of rows of the matrix
*/
- SquareSparseMatrix(uint_fast64_t rows)
- : rowCount(rows), nonZeroEntryCount(0), valueStorage(nullptr),
- diagonalStorage(nullptr),columnIndications(nullptr), rowIndications(nullptr),
+ SquareSparseMatrix(uint_fast64_t rows, uint_fast64_t cols)
+ : rowCount(rows), colCount(cols), nonZeroEntryCount(0),
+ internalStatus(MatrixStatus::UnInitialized), currentSize(0), lastRow(0) { }
+
+ /* Sadly, Delegate Constructors are not yet available with MSVC2012
+ //! Constructor
+ /*!
+ * Constructs a square sparse matrix object with the given number rows
+ * @param size The number of rows and cols in the matrix
+ */ /*
+ SquareSparseMatrix(uint_fast64_t size) : SquareSparseMatrix(size, size) { }
+ */
+
+ //! Constructor
+ /*!
+ * Constructs a square sparse matrix object with the given number rows
+ * @param size The number of rows and cols in the matrix
+ */
+ SquareSparseMatrix(uint_fast64_t size) : rowCount(size), colCount(size), nonZeroEntryCount(0),
internalStatus(MatrixStatus::UnInitialized), currentSize(0), lastRow(0) { }
//! Copy Constructor
@@ -84,7 +100,7 @@ public:
* @param ssm A reference to the matrix to be copied.
*/
SquareSparseMatrix(const SquareSparseMatrix<T> &ssm)
- : rowCount(ssm.rowCount), nonZeroEntryCount(ssm.nonZeroEntryCount),
+ : rowCount(ssm.rowCount), colCount(ssm.colCount), nonZeroEntryCount(ssm.nonZeroEntryCount),
internalStatus(ssm.internalStatus), currentSize(ssm.currentSize), lastRow(ssm.lastRow) {
LOG4CPLUS_WARN(logger, "Invoking copy constructor.");
// Check whether copying the matrix is safe.
@@ -98,24 +114,12 @@ public:
LOG4CPLUS_ERROR(logger, "Unable to allocate internal storage.");
throw std::bad_alloc();
} else {
- // Now that all storages have been prepared, copy over all
- // elements. Start by copying the elements of type value and
- // copy them seperately in order to invoke copy their copy
- // constructor. This may not be necessary, but it is safer to
- // do so in any case.
- for (uint_fast64_t i = 0; i < nonZeroEntryCount; ++i) {
- // use T() to force use of the copy constructor for complex T types
- valueStorage[i] = T(ssm.valueStorage[i]);
- }
- for (uint_fast64_t i = 0; i < rowCount; ++i) {
- // use T() to force use of the copy constructor for complex T types
- diagonalStorage[i] = T(ssm.diagonalStorage[i]);
- }
+ std::copy(ssm.valueStorage.begin(), ssm.valueStorage.end(), std::back_inserter(valueStorage));
// The elements that are not of the value type but rather the
// index type may be copied directly.
- std::copy(ssm.columnIndications, ssm.columnIndications + nonZeroEntryCount, columnIndications);
- std::copy(ssm.rowIndications, ssm.rowIndications + rowCount + 1, rowIndications);
+ std::copy(ssm.columnIndications.begin(), ssm.columnIndications.end(), std::back_inserter(columnIndications));
+ std::copy(ssm.rowIndications.begin(), ssm.rowIndications.end(), std::back_inserter(rowIndications));
}
}
}
@@ -126,18 +130,9 @@ public:
*/
~SquareSparseMatrix() {
setState(MatrixStatus::UnInitialized);
- if (valueStorage != nullptr) {
- delete[] valueStorage;
- }
- if (columnIndications != nullptr) {
- delete[] columnIndications;
- }
- if (rowIndications != nullptr) {
- delete[] rowIndications;
- }
- if (diagonalStorage != nullptr) {
- delete[] diagonalStorage;
- }
+ valueStorage.resize(0);
+ columnIndications.resize(0);
+ rowIndications.resize(0);
}
/*!
@@ -155,11 +150,11 @@ public:
triggerErrorState();
LOG4CPLUS_ERROR(logger, "Trying to initialize matrix that is not uninitialized.");
throw storm::exceptions::InvalidStateException("Trying to initialize matrix that is not uninitialized.");
- } else if (rowCount == 0) {
+ } else if ((rowCount == 0) || (colCount == 0)) {
triggerErrorState();
- LOG4CPLUS_ERROR(logger, "Trying to create initialize a matrix with 0 rows.");
- throw storm::exceptions::InvalidArgumentException("Trying to create initialize a matrix with 0 rows.");
- } else if (((rowCount * rowCount) - rowCount) < nonZeroEntries) {
+ LOG4CPLUS_ERROR(logger, "Trying to create initialize a matrix with 0 rows or 0 columns.");
+ throw storm::exceptions::InvalidArgumentException("Trying to create initialize a matrix with 0 rows or 0 columns.");
+ } else if ((rowCount * colCount) < nonZeroEntries) {
triggerErrorState();
LOG4CPLUS_ERROR(logger, "Trying to initialize a matrix with more non-zero entries than there can be.");
throw storm::exceptions::InvalidArgumentException("Trying to initialize a matrix with more non-zero entries than there can be.");
@@ -196,37 +191,72 @@ public:
throw storm::exceptions::InvalidArgumentException("Trying to initialize from an Eigen matrix that is not in compressed form.");
}
- // Compute the actual (i.e. non-diagonal) number of non-zero entries.
- nonZeroEntryCount = getEigenSparseMatrixCorrectNonZeroEntryCount(eigenSparseMatrix);
+ if (eigenSparseMatrix.rows() > this->rowCount) {
+ triggerErrorState();
+ LOG4CPLUS_ERROR(logger, "Trying to initialize from an Eigen matrix that has more rows than the target matrix.");
+ throw storm::exceptions::InvalidArgumentException("Trying to initialize from an Eigen matrix that has more rows than the target matrix.");
+ }
+ if (eigenSparseMatrix.cols() > this->colCount) {
+ triggerErrorState();
+ LOG4CPLUS_ERROR(logger, "Trying to initialize from an Eigen matrix that has more columns than the target matrix.");
+ throw storm::exceptions::InvalidArgumentException("Trying to initialize from an Eigen matrix that has more columns than the target matrix.");
+ }
+
+ const _Index entryCount = eigenSparseMatrix.nonZeros();
+ nonZeroEntryCount = entryCount;
lastRow = 0;
// Try to prepare the internal storage and throw an error in case of
// failure.
- if (!prepareInternalStorage()) {
- triggerErrorState();
- LOG4CPLUS_ERROR(logger, "Unable to allocate internal storage.");
- throw std::bad_alloc();
- } else {
- // Get necessary pointers to the contents of the Eigen matrix.
- const T* valuePtr = eigenSparseMatrix.valuePtr();
- const _Index* indexPtr = eigenSparseMatrix.innerIndexPtr();
- const _Index* outerPtr = eigenSparseMatrix.outerIndexPtr();
-
- // If the given matrix is in RowMajor format, copying can simply
- // be done by adding all values in order.
- // Direct copying is, however, prevented because we have to
- // separate the diagonal entries from others.
- if (isEigenRowMajor(eigenSparseMatrix)) {
- // Because of the RowMajor format outerSize evaluates to the
- // number of rows.
- const _Index rowCount = eigenSparseMatrix.outerSize();
- for (_Index row = 0; row < rowCount; ++row) {
- for (_Index col = outerPtr[row]; col < outerPtr[row + 1]; ++col) {
- addNextValue(row, indexPtr[col], valuePtr[col]);
- }
+
+ // Get necessary pointers to the contents of the Eigen matrix.
+ const T* valuePtr = eigenSparseMatrix.valuePtr();
+ const _Index* indexPtr = eigenSparseMatrix.innerIndexPtr();
+ const _Index* outerPtr = eigenSparseMatrix.outerIndexPtr();
+
+ // If the given matrix is in RowMajor format, copying can simply
+ // be done by adding all values in order.
+ // Direct copying is, however, prevented because we have to
+ // separate the diagonal entries from others.
+ if (isEigenRowMajor(eigenSparseMatrix)) {
+ // Because of the RowMajor format outerSize evaluates to the
+ // number of rows.
+ if (!prepareInternalStorage(false)) {
+ triggerErrorState();
+ LOG4CPLUS_ERROR(logger, "Unable to allocate internal storage.");
+ throw std::bad_alloc();
+ } else {
+ if ((eigenSparseMatrix.innerSize() > nonZeroEntryCount) || (entryCount > nonZeroEntryCount)) {
+ triggerErrorState();
+ LOG4CPLUS_ERROR(logger, "Invalid internal composition of Eigen Sparse Matrix");
+ throw storm::exceptions::InvalidArgumentException("Invalid internal composition of Eigen Sparse Matrix");
+ }
+ std::vector<uint_fast64_t> eigenColumnTemp;
+ std::vector<uint_fast64_t> eigenRowTemp;
+ std::vector<T> eigenValueTemp;
+ uint_fast64_t outerSize = eigenSparseMatrix.outerSize() + 1;
+
+ for (uint_fast64_t i = 0; i < entryCount; ++i) {
+ eigenColumnTemp.push_back(indexPtr[i]);
+ eigenValueTemp.push_back(valuePtr[i]);
}
+ for (uint_fast64_t i = 0; i < outerSize; ++i) {
+ eigenRowTemp.push_back(outerPtr[i]);
+ }
+
+ std::copy(eigenRowTemp.begin(), eigenRowTemp.end(), std::back_inserter(this->rowIndications));
+ std::copy(eigenColumnTemp.begin(), eigenColumnTemp.end(), std::back_inserter(this->columnIndications));
+ std::copy(eigenValueTemp.begin(), eigenValueTemp.end(), std::back_inserter(this->valueStorage));
+
+ currentSize = entryCount;
+ lastRow = rowCount;
+ }
+ } else {
+ if (!prepareInternalStorage()) {
+ triggerErrorState();
+ LOG4CPLUS_ERROR(logger, "Unable to allocate internal storage.");
+ throw std::bad_alloc();
} else {
- const _Index entryCount = eigenSparseMatrix.nonZeros();
// Because of the ColMajor format outerSize evaluates to the
// number of columns.
const _Index colCount = eigenSparseMatrix.outerSize();
@@ -250,8 +280,7 @@ public:
// add it in case it is also in the current row.
if ((positions[currentColumn] < outerPtr[currentColumn + 1])
&& (indexPtr[positions[currentColumn]] == currentRow)) {
- addNextValue(currentRow, currentColumn,
- valuePtr[positions[currentColumn]]);
+ addNextValue(currentRow, currentColumn, valuePtr[positions[currentColumn]]);
// Remember that we found one more non-zero element.
++i;
// Mark this position as "used".
@@ -268,8 +297,8 @@ public:
}
delete[] positions;
}
- setState(MatrixStatus::Initialized);
}
+ setState(MatrixStatus::Initialized);
}
/*!
@@ -283,30 +312,27 @@ public:
void addNextValue(const uint_fast64_t row, const uint_fast64_t col, const T& value) {
// Check whether the given row and column positions are valid and throw
// error otherwise.
- if ((row > rowCount) || (col > rowCount)) {
+ if ((row > rowCount) || (col > colCount)) {
triggerErrorState();
LOG4CPLUS_ERROR(logger, "Trying to add a value at illegal position (" << row << ", " << col << ").");
throw storm::exceptions::OutOfRangeException("Trying to add a value at illegal position.");
}
- if (row == col) { // Set a diagonal element.
- diagonalStorage[row] = value;
- } else { // Set a non-diagonal element.
- // If we switched to another row, we have to adjust the missing
- // entries in the row_indications array.
- if (row != lastRow) {
- for (uint_fast64_t i = lastRow + 1; i <= row; ++i) {
- rowIndications[i] = currentSize;
- }
- lastRow = row;
+
+ // If we switched to another row, we have to adjust the missing
+ // entries in the row_indications array.
+ if (row != lastRow) {
+ for (uint_fast64_t i = lastRow + 1; i <= row; ++i) {
+ rowIndications[i] = currentSize;
}
+ lastRow = row;
+ }
- // Finally, set the element and increase the current size.
- valueStorage[currentSize] = value;
- columnIndications[currentSize] = col;
+ // Finally, set the element and increase the current size.
+ valueStorage[currentSize] = value;
+ columnIndications[currentSize] = col;
- ++currentSize;
- }
+ ++currentSize;
}
/*
@@ -355,18 +381,12 @@ public:
*/
inline bool getValue(uint_fast64_t row, uint_fast64_t col, T* const target) const {
// Check for illegal access indices.
- if ((row > rowCount) || (col > rowCount)) {
+ if ((row > rowCount) || (col > colCount)) {
LOG4CPLUS_ERROR(logger, "Trying to read a value from illegal position (" << row << ", " << col << ").");
throw storm::exceptions::OutOfRangeException("Trying to read a value from illegal position.");
return false;
}
- // Read elements on the diagonal directly.
- if (row == col) {
- *target = diagonalStorage[row];
- return true;
- }
-
// In case the element is not on the diagonal, we have to iterate
// over the accessed row to find the element.
uint_fast64_t rowStart = rowIndications[row];
@@ -405,17 +425,12 @@ public:
*/
inline T& getValue(uint_fast64_t row, uint_fast64_t col) {
// Check for illegal access indices.
- if ((row > rowCount) || (col > rowCount)) {
+ if ((row > rowCount) || (col > colCount)) {
LOG4CPLUS_ERROR(logger, "Trying to read a value from illegal position (" << row << ", " << col << ").");
throw storm::exceptions::OutOfRangeException("Trying to read a value from illegal position.");
}
- // Read elements on the diagonal directly.
- if (row == col) {
- return diagonalStorage[row];
- }
-
- // In case the element is not on the diagonal, we have to iterate
+ // we have to iterate
// over the accessed row to find the element.
uint_fast64_t rowStart = rowIndications[row];
uint_fast64_t rowEnd = rowIndications[row + 1];
@@ -445,20 +460,18 @@ public:
}
/*!
- * Returns a pointer to the value storage of the matrix. This storage does
- * *not* include elements on the diagonal.
- * @return A pointer to the value storage of the matrix.
+ * Returns the number of columns of the matrix.
*/
- T* getStoragePointer() const {
- return valueStorage;
+ uint_fast64_t getColumnCount() const {
+ return colCount;
}
/*!
- * Returns a pointer to the storage of elements on the diagonal.
- * @return A pointer to the storage of elements on the diagonal.
+ * Returns a pointer to the value storage of the matrix.
+ * @return A pointer to the value storage of the matrix.
*/
- T* getDiagonalStoragePointer() const {
- return diagonalStorage;
+ std::vector<T>& getStoragePointer() {
+ return valueStorage;
}
/*!
@@ -467,17 +480,17 @@ public:
* @return A pointer to the array that stores the start indices of non-zero
* entries in the value storage for each row.
*/
- uint_fast64_t* getRowIndicationsPointer() const {
+ std::vector<uint_fast64_t>& getRowIndicationsPointer() {
return rowIndications;
}
/*!
* Returns a pointer to an array that stores the column of each non-zero
- * element that is not on the diagonal.
+ * element.
* @return A pointer to an array that stores the column of each non-zero
- * element that is not on the diagonal.
+ * element.
*/
- uint_fast64_t* getColumnIndicationsPointer() const {
+ std::vector<uint_fast64_t>& getColumnIndicationsPointer() {
return columnIndications;
}
@@ -548,10 +561,6 @@ public:
#define STORM_USE_TRIPLETCONVERT
# ifdef STORM_USE_TRIPLETCONVERT
- // FIXME: Wouldn't it be more efficient to add the elements in
- // order including the diagonal elements? Otherwise, Eigen has to
- // perform some sorting.
-
// Prepare the triplet storage.
typedef Eigen::Triplet<T> IntTriplet;
std::vector<IntTriplet> tripletList;
@@ -572,12 +581,6 @@ public:
}
}
- // Then add the elements on the diagonal.
- for (uint_fast64_t i = 0; i < rowCount; ++i) {
- if (diagonalStorage[i] == 0) zeroCount++;
- tripletList.push_back(IntTriplet(static_cast<int_fast32_t>(i), static_cast<int_fast32_t>(i), diagonalStorage[i]));
- }
-
// Let Eigen create a matrix from the given list of triplets.
mat->setFromTriplets(tripletList.begin(), tripletList.end());
@@ -596,10 +599,6 @@ public:
rowStart = rowIndications[row];
rowEnd = rowIndications[row + 1];
- // Insert the element on the diagonal.
- mat->insert(row, row) = diagonalStorage[row];
- count++;
-
// Insert the elements that are not on the diagonal
while (rowStart < rowEnd) {
mat->insert(row, columnIndications[rowStart]) = valueStorage[rowStart];
@@ -628,19 +627,6 @@ public:
return nonZeroEntryCount;
}
- /*!
- * Returns the number of non-zero entries on the diagonal.
- * @return The number of non-zero entries on the diagonal.
- */
- uint_fast64_t getDiagonalNonZeroEntryCount() const {
- uint_fast64_t result = 0;
- T zero(0);
- for (uint_fast64_t i = 0; i < rowCount; ++i) {
- if (diagonalStorage[i] != zero) ++result;
- }
- return result;
- }
-
/*!
* This function makes the rows given by the bit vector absorbing.
* @param rows A bit vector indicating which rows to make absorbing.
@@ -658,7 +644,7 @@ public:
/*!
* This function makes the given row absorbing. This means that all
* entries in will be set to 0 and the value 1 will be written
- * to the element on the diagonal.
+ * to the element on the (pseudo-) diagonal.
* @param row The row to be made absorbing.
* @returns True iff the operation was successful.
*/
@@ -675,13 +661,31 @@ public:
uint_fast64_t rowStart = rowIndications[row];
uint_fast64_t rowEnd = rowIndications[row + 1];
+ if (rowStart >= rowEnd) {
+ LOG4CPLUS_ERROR(logger, "The row " << row << " can not be made absorbing, no state in row, would have to recreate matrix!");
+ throw storm::exceptions::InvalidStateException("A row can not be made absorbing, no state in row, would have to recreate matrix!");
+ }
+ uint_fast64_t pseudoDiagonal = row % colCount;
+
+ bool foundDiagonal = false;
while (rowStart < rowEnd) {
- valueStorage[rowStart] = storm::utility::constGetZero<T>();
+ if (!foundDiagonal && columnIndications[rowStart] >= pseudoDiagonal) {
+ foundDiagonal = true;
+ // insert/replace the diagonal entry
+ columnIndications[rowStart] = pseudoDiagonal;
+ valueStorage[rowStart] = storm::utility::constGetOne<T>();
+ } else {
+ valueStorage[rowStart] = storm::utility::constGetZero<T>();
+ }
++rowStart;
}
- // Set the element on the diagonal to one.
- diagonalStorage[row] = storm::utility::constGetOne<T>();
+ if (!foundDiagonal) {
+ --rowStart;
+ columnIndications[rowStart] = pseudoDiagonal;
+ valueStorage[rowStart] = storm::utility::constGetOne<T>();
+ }
+
return true;
}
@@ -763,8 +767,6 @@ public:
// Copy over selected entries.
uint_fast64_t rowCount = 0;
for (auto rowIndex : constraint) {
- result->addNextValue(rowCount, rowCount, diagonalStorage[rowIndex]);
-
for (uint_fast64_t i = rowIndications[rowIndex]; i < rowIndications[rowIndex + 1]; ++i) {
if (constraint.get(columnIndications[i])) {
result->addNextValue(rowCount, bitsSetBeforeIndex[columnIndications[i]], valueStorage[i]);
@@ -794,8 +796,17 @@ public:
*/
void invertDiagonal() {
T one(1);
- for (uint_fast64_t i = 0; i < rowCount; ++i) {
- diagonalStorage[i] = one - diagonalStorage[i];
+ for (uint_fast64_t row = 0; row < rowCount; ++row) {
+ uint_fast64_t rowStart = rowIndications[row];
+ uint_fast64_t rowEnd = rowIndications[row + 1];
+ uint_fast64_t pseudoDiagonal = row % colCount;
+ while (rowStart < rowEnd) {
+ if (columnIndications[rowStart] == pseudoDiagonal) {
+ valueStorage[rowStart] = one - valueStorage[rowStart];
+ break;
+ }
+ ++rowStart;
+ }
}
}
@@ -803,8 +814,16 @@ public:
* Negates all non-zero elements that are not on the diagonal.
*/
void negateAllNonDiagonalElements() {
- for (uint_fast64_t i = 0; i < nonZeroEntryCount; ++i) {
- valueStorage[i] = - valueStorage[i];
+ for (uint_fast64_t row = 0; row < rowCount; ++row) {
+ uint_fast64_t rowStart = rowIndications[row];
+ uint_fast64_t rowEnd = rowIndications[row + 1];
+ uint_fast64_t pseudoDiagonal = row % colCount;
+ while (rowStart < rowEnd) {
+ if (columnIndications[rowStart] != pseudoDiagonal) {
+ valueStorage[rowStart] = - valueStorage[rowStart];
+ }
+ ++rowStart;
+ }
}
}
@@ -844,7 +863,6 @@ public:
// in case the given matrix does not have a non-zero element at this column position, or
// multiply the two entries and add the result to the corresponding position in the vector.
for (uint_fast64_t row = 0; row < rowCount && row < otherMatrix.rowCount; ++row) {
- (*result)[row] += diagonalStorage[row] * otherMatrix.diagonalStorage[row];
for (uint_fast64_t element = rowIndications[row], nextOtherElement = otherMatrix.rowIndications[row]; element < rowIndications[row + 1] && nextOtherElement < otherMatrix.rowIndications[row + 1]; ++element) {
if (columnIndications[element] < otherMatrix.columnIndications[nextOtherElement]) {
continue;
@@ -868,25 +886,23 @@ public:
uint_fast64_t getSizeInMemory() const {
uint_fast64_t size = sizeof(*this);
// Add value_storage size.
- size += sizeof(T) * nonZeroEntryCount;
- // Add diagonal_storage size.
- size += sizeof(T) * (rowCount + 1);
+ size += sizeof(T) * valueStorage.capacity();
// Add column_indications size.
- size += sizeof(uint_fast64_t) * nonZeroEntryCount;
+ size += sizeof(uint_fast64_t) * columnIndications.capacity();
// Add row_indications size.
- size += sizeof(uint_fast64_t) * (rowCount + 1);
+ size += sizeof(uint_fast64_t) * rowIndications.capacity();
return size;
}
/*!
* Returns an iterator to the columns of the non-zero entries of the given
- * row that are not on the diagonal.
+ * row.
* @param row The row whose columns the iterator will return.
* @return An iterator to the columns of the non-zero entries of the given
- * row that are not on the diagonal.
+ * row.
*/
- constIndexIterator beginConstColumnNoDiagIterator(uint_fast64_t row) const {
- return this->columnIndications + this->rowIndications[row];
+ constIndexIterator beginConstColumnIterator(uint_fast64_t row) const {
+ return &(this->columnIndications[0]) + this->rowIndications[row];
}
/*!
@@ -894,18 +910,18 @@ public:
* @param row The row for which the iterator should point to the past-the-end
* element.
*/
- constIndexIterator endConstColumnNoDiagIterator(uint_fast64_t row) const {
- return this->columnIndications + this->rowIndications[row + 1];
+ constIndexIterator endConstColumnIterator(uint_fast64_t row) const {
+ return &(this->columnIndications[0]) + this->rowIndications[row + 1];
}
/*!
* Returns an iterator over the elements of the given row. The iterator
- * will include neither the diagonal element nor zero entries.
+ * will include no zero entries.
* @param row The row whose elements the iterator will return.
* @return An iterator over the elements of the given row.
*/
- constIterator beginConstNoDiagIterator(uint_fast64_t row) const {
- return this->valueStorage + this->rowIndications[row];
+ constIterator beginConstIterator(uint_fast64_t row) const {
+ return &(this->valueStorage[0]) + this->rowIndications[row];
}
/*!
* Returns an iterator pointing to the first element after the given
@@ -914,32 +930,28 @@ public:
* past-the-end element.
* @return An iterator to the element after the given row.
*/
- constIterator endConstNoDiagIterator(uint_fast64_t row) const {
- return this->valueStorage + this->rowIndications[row + 1];
+ constIterator endConstIterator(uint_fast64_t row) const {
+ return &(this->valueStorage[0]) + this->rowIndications[row + 1];
}
/*!
* @brief Calculate sum of all entries in given row.
*
- * Adds up all values in the given row (including the diagonal value)
+ * Adds up all values in the given row
* and returns the sum.
* @param row The row that should be added up.
* @return Sum of the row.
*/
T getRowSum(uint_fast64_t row) const {
- T sum = this->diagonalStorage[row];
- for (auto it = this->beginConstNoDiagIterator(row); it != this->endConstNoDiagIterator(row); it++) {
+ T sum = storm::utility::constGetZero<T>();
+ for (auto it = this->beginConstIterator(row); it != this->endConstIterator(row); it++) {
sum += *it;
}
return sum;
}
void print() const {
- std::cout << "diag: --------------------------------" << std::endl;
- for (uint_fast64_t i = 0; i < rowCount; ++i) {
- std::cout << "(" << i << "," << i << ") = " << diagonalStorage[i] << std::endl;
- }
- std::cout << "non diag: ----------------------------" << std::endl;
+ std::cout << "entries: ----------------------------" << std::endl;
for (uint_fast64_t i = 0; i < rowCount; ++i) {
for (uint_fast64_t j = rowIndications[i]; j < rowIndications[i + 1]; ++j) {
std::cout << "(" << i << "," << columnIndications[j] << ") = " << valueStorage[j] << std::endl;
@@ -955,31 +967,31 @@ private:
uint_fast64_t rowCount;
/*!
- * The number of non-zero elements that are not on the diagonal.
+ * The number of columns of the matrix.
*/
- uint_fast64_t nonZeroEntryCount;
+ uint_fast64_t colCount;
/*!
- * Stores all non-zero values that are not on the diagonal.
+ * The number of non-zero elements.
*/
- T* valueStorage;
+ uint_fast64_t nonZeroEntryCount;
/*!
- * Stores all elements on the diagonal, even the ones that are zero.
+ * Stores all non-zero values.
*/
- T* diagonalStorage;
+ std::vector<T> valueStorage;
/*!
- * Stores the column for each non-zero element that is not on the diagonal.
+ * Stores the column for each non-zero element.
*/
- uint_fast64_t* columnIndications;
+ std::vector<uint_fast64_t> columnIndications;
/*!
- * Array containing the boundaries (indices) in the value_storage array
+ * Vector containing the boundaries (indices) in the value_storage array
* for each row. All elements of value_storage with indices between the
* i-th and the (i+1)-st element of this array belong to row i.
*/
- uint_fast64_t* rowIndications;
+ std::vector<uint_fast64_t> rowIndications;
/*!
* The internal status of the matrix.
@@ -1017,24 +1029,37 @@ private:
/*!
* Prepares the internal CSR storage. For this, it requires
* non_zero_entry_count and row_count to be set correctly.
+ * @param alsoPerformAllocation If set to true, all entries are pre-allocated. This is the default.
* @return True on success, false otherwise (allocation failed).
*/
- bool prepareInternalStorage() {
- // Set up the arrays for the elements that are not on the diagonal.
- valueStorage = new (std::nothrow) T[nonZeroEntryCount]();
- columnIndications = new (std::nothrow) uint_fast64_t[nonZeroEntryCount]();
-
- // Set up the row_indications array and reserve one element more than
- // there are rows in order to put a sentinel element at the end,
- // which eases iteration process.
- rowIndications = new (std::nothrow) uint_fast64_t[rowCount + 1]();
-
- // Set up the array for the elements on the diagonal.
- diagonalStorage = new (std::nothrow) T[rowCount]();
+ bool prepareInternalStorage(const bool alsoPerformAllocation) {
+ if (alsoPerformAllocation) {
+ // Set up the arrays for the elements that are not on the diagonal.
+ valueStorage.resize(nonZeroEntryCount, storm::utility::constGetZero<T>());
+ columnIndications.resize(nonZeroEntryCount, 0);
+
+ // Set up the row_indications vector and reserve one element more than
+ // there are rows in order to put a sentinel element at the end,
+ // which eases iteration process.
+ rowIndications.resize(rowCount + 1, 0);
+
+ // Return whether all the allocations could be made without error.
+ return ((valueStorage.capacity() >= nonZeroEntryCount) && (columnIndications.capacity() >= nonZeroEntryCount)
+ && (rowIndications.capacity() >= (rowCount + 1)));
+ } else {
+ valueStorage.reserve(nonZeroEntryCount);
+ columnIndications.reserve(nonZeroEntryCount);
+ rowIndications.reserve(rowCount + 1);
+ return true;
+ }
+ }
- // Return whether all the allocations could be made without error.
- return ((valueStorage != NULL) && (columnIndications != NULL)
- && (rowIndications != NULL) && (diagonalStorage != NULL));
+ /*!
+ * Shorthand for prepareInternalStorage(true)
+ * @see prepareInternalStorage(const bool)
+ */
+ bool prepareInternalStorage() {
+ return this->prepareInternalStorage(true);
}
/*!
@@ -1060,53 +1085,6 @@ private:
return false;
}
- /*!
- * Helper function to determine the number of non-zero elements that are
- * not on the diagonal of the given Eigen matrix.
- * @param eigen_sparse_matrix The Eigen matrix to analyze.
- * @return The number of non-zero elements that are not on the diagonal of
- * the given Eigen matrix.
- */
- template<typename _Scalar, int _Options, typename _Index>
- _Index getEigenSparseMatrixCorrectNonZeroEntryCount(const Eigen::SparseMatrix<_Scalar, _Options, _Index>& eigen_sparse_matrix) const {
- const _Index* indexPtr = eigen_sparse_matrix.innerIndexPtr();
- const _Index* outerPtr = eigen_sparse_matrix.outerIndexPtr();
-
- const _Index entryCount = eigen_sparse_matrix.nonZeros();
- const _Index outerCount = eigen_sparse_matrix.outerSize();
-
- uint_fast64_t diagNonZeros = 0;
-
- // For RowMajor, row is the current row and col the column and for
- // ColMajor, row is the current column and col the row, but this is
- // not important as we are only looking for elements on the diagonal.
- _Index innerStart = 0;
- _Index innerEnd = 0;
- _Index innerMid = 0;
- for (_Index row = 0; row < outerCount; ++row) {
- innerStart = outerPtr[row];
- innerEnd = outerPtr[row + 1] - 1;
-
- // Now use binary search (but defer equality detection).
- while (innerStart < innerEnd) {
- innerMid = innerStart + ((innerEnd - innerStart) / 2);
-
- if (indexPtr[innerMid] < row) {
- innerStart = innerMid + 1;
- } else {
- innerEnd = innerMid;
- }
- }
-
- // Check whether we have found an element on the diagonal.
- if ((innerStart == innerEnd) && (indexPtr[innerStart] == row)) {
- ++diagNonZeros;
- }
- }
-
- return static_cast<_Index>(entryCount - diagNonZeros);
- }
-
};
} // namespace storage
diff --git a/src/utility/IoUtility.cpp b/src/utility/IoUtility.cpp
index e3fc7dbaf..4687e719e 100644
--- a/src/utility/IoUtility.cpp
+++ b/src/utility/IoUtility.cpp
@@ -1,9 +1,9 @@
/*
- * IoUtility.cpp
- *
- * Created on: 17.10.2012
- * Author: Thomas Heinemann
- */
+* IoUtility.cpp
+*
+* Created on: 17.10.2012
+* Author: Thomas Heinemann
+*/
#include "src/utility/IoUtility.h"
#include "src/parser/DeterministicSparseTransitionParser.h"
@@ -13,63 +13,58 @@
namespace storm {
-namespace utility {
-
-void dtmcToDot(storm::models::Dtmc<double> const &dtmc, std::string filename) {
- std::shared_ptr<storm::storage::SquareSparseMatrix<double>> matrix(dtmc.getTransitionProbabilityMatrix());
- double* diagonal_storage = matrix->getDiagonalStoragePointer();
-
- std::ofstream file;
- file.open(filename);
-
- file << "digraph dtmc {\n";
-
- //Specify the nodes and their labels
- for (uint_fast64_t i = 1; i < dtmc.getNumberOfStates(); i++) {
- file << "\t" << i << "[label=\"" << i << "\\n{";
- char komma=' ';
- std::set<std::string> propositions = dtmc.getPropositionsForState(i);
- for(auto it = propositions.begin();
- it != propositions.end();
- it++) {
- file << komma << *it;
- komma=',';
- }
-
- file << " }\"];\n";
-
- }
-
- for (uint_fast64_t row = 0; row < dtmc.getNumberOfStates(); row++ ) {
- //write diagonal entry/self loop first
- if (diagonal_storage[row] != 0) {
- file << "\t" << row << " -> " << row << " [label=" << diagonal_storage[row] <<"]\n";
- }
- //Then, iterate through the row and write each non-diagonal value into the file
- for ( auto it = matrix->beginConstColumnNoDiagIterator(row);
- it != matrix->endConstColumnNoDiagIterator(row);
- it++) {
- double value = 0;
- matrix->getValue(row,*it,&value);
- file << "\t" << row << " -> " << *it << " [label=" << value << "]\n";
- }
- }
-
- file << "}\n";
- file.close();
-}
+ namespace utility {
-//TODO: Should this stay here or be integrated in the new parser structure?
-/*storm::models::Dtmc<double>* parseDTMC(std::string const &tra_file, std::string const &lab_file) {
- storm::parser::DeterministicSparseTransitionParser tp(tra_file);
- uint_fast64_t node_count = tp.getMatrix()->getRowCount();
+ void dtmcToDot(storm::models::Dtmc<double> const &dtmc, std::string filename) {
+ std::shared_ptr<storm::storage::SquareSparseMatrix<double>> matrix(dtmc.getTransitionProbabilityMatrix());
+ std::ofstream file;
+ file.open(filename);
- storm::parser::AtomicPropositionLabelingParser lp(node_count, lab_file);
+ file << "digraph dtmc {\n";
- storm::models::Dtmc<double>* result = new storm::models::Dtmc<double>(tp.getMatrix(), lp.getLabeling());
- return result;
-}*/
+ //Specify the nodes and their labels
+ for (uint_fast64_t i = 1; i < dtmc.getNumberOfStates(); i++) {
+ file << "\t" << i << "[label=\"" << i << "\\n{";
+ char komma=' ';
+ std::set<std::string> propositions = dtmc.getPropositionsForState(i);
+ for(auto it = propositions.begin();
+ it != propositions.end();
+ it++) {
+ file << komma << *it;
+ komma=',';
+ }
-}
+ file << " }\"];\n";
+
+ }
+
+ for (uint_fast64_t row = 0; row < dtmc.getNumberOfStates(); row++ ) {
+
+ //Then, iterate through the row and write each non-diagonal value into the file
+ for ( auto it = matrix->beginConstColumnIterator(row);
+ it != matrix->endConstColumnIterator(row);
+ it++) {
+ double value = 0;
+ matrix->getValue(row,*it,&value);
+ file << "\t" << row << " -> " << *it << " [label=" << value << "]\n";
+ }
+ }
+
+ file << "}\n";
+ file.close();
+ }
+
+ //TODO: Should this stay here or be integrated in the new parser structure?
+ /*storm::models::Dtmc<double>* parseDTMC(std::string const &tra_file, std::string const &lab_file) {
+ storm::parser::DeterministicSparseTransitionParser tp(tra_file);
+ uint_fast64_t node_count = tp.getMatrix()->getRowCount();
+
+ storm::parser::AtomicPropositionLabelingParser lp(node_count, lab_file);
+
+ storm::models::Dtmc<double>* result = new storm::models::Dtmc<double>(tp.getMatrix(), lp.getLabeling());
+ return result;
+ }*/
+
+ }
}
diff --git a/test/parser/ReadTraFileTest.cpp b/test/parser/ReadTraFileTest.cpp
index 77377d857..752351e33 100644
--- a/test/parser/ReadTraFileTest.cpp
+++ b/test/parser/ReadTraFileTest.cpp
@@ -53,13 +53,13 @@ TEST(ReadTraFileTest, ParseFileTest1) {
ASSERT_TRUE(result->getValue(3,2,&val));
ASSERT_EQ(val,0.0806451612903225806451612903225812);
- ASSERT_TRUE(result->getValue(3,3,&val));
+ ASSERT_FALSE(result->getValue(3,3,&val));
ASSERT_EQ(val,0);
ASSERT_TRUE(result->getValue(3,4,&val));
ASSERT_EQ(val,0.080645161290322580645161290322581);
- ASSERT_TRUE(result->getValue(4,4,&val));
+ ASSERT_FALSE(result->getValue(4,4,&val));
ASSERT_EQ(val,0);
delete parser;
diff --git a/test/storage/SquareSparseMatrixTest.cpp b/test/storage/SquareSparseMatrixTest.cpp
index 61d64c933..289a52fc6 100644
--- a/test/storage/SquareSparseMatrixTest.cpp
+++ b/test/storage/SquareSparseMatrixTest.cpp
@@ -116,11 +116,7 @@ TEST(SquareSparseMatrixTest, Test) {
for (int row = 15; row < 24; ++row) {
for (int col = 1; col <= 25; ++col) {
target = 1;
- if (row != col) {
- ASSERT_FALSE(ssm->getValue(row, col, &target));
- } else {
- ASSERT_TRUE(ssm->getValue(row, col, &target));
- }
+ ASSERT_FALSE(ssm->getValue(row, col, &target));
ASSERT_EQ(0, target);
}
@@ -245,7 +241,7 @@ TEST(SquareSparseMatrixTest, ConversionFromSparseEigen_ColMajor_SparseMatrixTest
TEST(SquareSparseMatrixTest, ConversionFromSparseEigen_RowMajor_SparseMatrixTest) {
// 10 rows, 15 non zero entries
- storm::storage::SquareSparseMatrix<int> *ssm = new storm::storage::SquareSparseMatrix<int>(10);
+ storm::storage::SquareSparseMatrix<int> *ssm = new storm::storage::SquareSparseMatrix<int>(10, 10);
ASSERT_EQ(ssm->getState(), storm::storage::SquareSparseMatrix<int>::MatrixStatus::UnInitialized);
Eigen::SparseMatrix<int, Eigen::RowMajor> esm(10, 10);
@@ -253,7 +249,7 @@ TEST(SquareSparseMatrixTest, ConversionFromSparseEigen_RowMajor_SparseMatrixTest
typedef Eigen::Triplet<int> IntTriplet;
std::vector<IntTriplet> tripletList;
tripletList.reserve(15);
- tripletList.push_back(IntTriplet(1, 0, 0));
+ tripletList.push_back(IntTriplet(1, 0, 15));
tripletList.push_back(IntTriplet(1, 1, 1));
tripletList.push_back(IntTriplet(1, 2, 2));
tripletList.push_back(IntTriplet(1, 3, 3));
@@ -281,11 +277,15 @@ TEST(SquareSparseMatrixTest, ConversionFromSparseEigen_RowMajor_SparseMatrixTest
ASSERT_NO_THROW(ssm->finalize());
ASSERT_EQ(ssm->getState(), storm::storage::SquareSparseMatrix<int>::MatrixStatus::ReadReady);
+
+ const std::vector<uint_fast64_t> rowP = ssm->getRowIndicationsPointer();
+ const std::vector<uint_fast64_t> colP = ssm->getColumnIndicationsPointer();
+ const std::vector<int> valP = ssm->getStoragePointer();
int target = -1;
-
for (auto &coeff: tripletList) {
- ASSERT_TRUE(ssm->getValue(coeff.row(), coeff.col(), &target));
+ bool retVal = ssm->getValue(coeff.row(), coeff.col(), &target);
+ ASSERT_TRUE(retVal);
ASSERT_EQ(target, coeff.value());
}
@@ -300,7 +300,7 @@ TEST(SquareSparseMatrixTest, ConversionToSparseEigen_RowMajor_SparseMatrixTest)
values[i] = i;
}
- ASSERT_NO_THROW(ssm->initialize(100 - 10));
+ ASSERT_NO_THROW(ssm->initialize(100));
ASSERT_EQ(ssm->getState(), storm::storage::SquareSparseMatrix<int>::MatrixStatus::Initialized);
for (uint_fast32_t row = 0; row < 10; ++row) {