diff --git a/CHANGELOG.md b/CHANGELOG.md index 6c81481e8d..d9b1aa3918 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -9,6 +9,7 @@ Version 1.7.x ------------- ## Version 1.7.1 +- Experimental support for compiling on Apple Silicon - Added SoPlex as a possible LP solver - Upgraded shipped version of sylvan - Upgraded repo / version for carl (for polynomials) diff --git a/CMakeLists.txt b/CMakeLists.txt index 674f7018ba..6320878a49 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -191,7 +191,7 @@ set(SHIPPED_CARL_USE_GINAC ON) # Detect Apple Silicon and adjust settings if(APPLE AND ${CMAKE_SYSTEM_PROCESSOR} MATCHES arm64) message(STATUS "Storm - Detected that target system uses Apple Silicon.") - message(WARNING "Storm currently does not compile natively on Apple Silicon! Please consider building Storm using Rosetta 2. For more information visit https://www.stormchecker.org/documentation/obtain-storm/apple-silicon.html") + message(WARNING "Compiling natively on Apple Silicon is experimental. Pleare report issues to support@stormchecker.org. For more information visit https://www.stormchecker.org/documentation/obtain-storm/apple-silicon.html") set(APPLE_SILICON 1) if(STORM_USE_CLN_EA OR STORM_USE_CLN_RF) message(WARNING "CLN and GiNaC are currently not supported on Apple Silicon-based architectures. Disabling Storm and carl usage of the libraries.") diff --git a/resources/3rdparty/CMakeLists.txt b/resources/3rdparty/CMakeLists.txt index 4174549e50..12307ae806 100644 --- a/resources/3rdparty/CMakeLists.txt +++ b/resources/3rdparty/CMakeLists.txt @@ -46,22 +46,42 @@ list(APPEND STORM_DEP_TARGETS gmm) ## ############################################################# -# Checkout Eigen version 3.3.9 -message (STATUS "Storm - Including Eigen 3.3.9.") -ExternalProject_Add( - eigen_src - GIT_REPOSITORY https://gitlab.com/libeigen/eigen.git - GIT_SHALLOW 1 - GIT_TAG 0fd6b4f71dd85b2009ee4d1aeb296e2c11fc9d68 - SOURCE_DIR ${STORM_3RDPARTY_INCLUDE_DIR}/StormEigen - PREFIX ${STORM_3RDPARTY_BINARY_DIR}/StormEigen-3.3.9 - PATCH_COMMAND git apply ${STORM_3RDPARTY_SOURCE_DIR}/patches/eigen.patch - UPDATE_COMMAND "" - CONFIGURE_COMMAND "" - BUILD_COMMAND "" - INSTALL_COMMAND "" - LOG_INSTALL ON -) +# Checkout Eigen +# Version 3.3.9 has a known issue when compiled on Apple Silicon (FullPivLU segfaults, relevant for multi-obj model checking). +# Version 3.4.0 has issues with gcc https://github.com/moves-rwth/storm/issues/162 +# We use this as a workaround: +if (APPLE_SILICON) + message (STATUS "Storm - Including Eigen 3.4.0 commit b0eded878d5d162d61583a286c0d8a45406ad1bc.") + ExternalProject_Add( + eigen_src + GIT_REPOSITORY https://gitlab.com/libeigen/eigen.git + GIT_TAG b0eded878d5d162d61583a286c0d8a45406ad1bc + SOURCE_DIR ${STORM_3RDPARTY_INCLUDE_DIR}/StormEigen + PREFIX ${STORM_3RDPARTY_BINARY_DIR}/StormEigen-3.4.0 + PATCH_COMMAND git apply ${STORM_3RDPARTY_SOURCE_DIR}/patches/eigen340.patch + UPDATE_COMMAND "" + CONFIGURE_COMMAND "" + BUILD_COMMAND "" + INSTALL_COMMAND "" + LOG_INSTALL ON + ) +else() + message (STATUS "Storm - Including Eigen 3.3.9.") + ExternalProject_Add( + eigen_src + GIT_REPOSITORY https://gitlab.com/libeigen/eigen.git + GIT_SHALLOW 1 + GIT_TAG 0fd6b4f71dd85b2009ee4d1aeb296e2c11fc9d68 + SOURCE_DIR ${STORM_3RDPARTY_INCLUDE_DIR}/StormEigen + PREFIX ${STORM_3RDPARTY_BINARY_DIR}/StormEigen-3.3.9 + PATCH_COMMAND git apply ${STORM_3RDPARTY_SOURCE_DIR}/patches/eigen339.patch + UPDATE_COMMAND "" + CONFIGURE_COMMAND "" + BUILD_COMMAND "" + INSTALL_COMMAND "" + LOG_INSTALL ON + ) +endif() add_imported_library_interface(StormEigen "${STORM_3RDPARTY_INCLUDE_DIR}/StormEigen") list(APPEND STORM_DEP_TARGETS StormEigen) add_dependencies(StormEigen eigen_src) diff --git a/resources/3rdparty/patches/eigen.patch b/resources/3rdparty/patches/eigen339.patch similarity index 100% rename from resources/3rdparty/patches/eigen.patch rename to resources/3rdparty/patches/eigen339.patch diff --git a/resources/3rdparty/patches/eigen340.patch b/resources/3rdparty/patches/eigen340.patch new file mode 100644 index 0000000000..1480655f03 --- /dev/null +++ b/resources/3rdparty/patches/eigen340.patch @@ -0,0 +1,226 @@ +From 938211b9e47c1650fb057e7ff88957cb2d25af7e Mon Sep 17 00:00:00 2001 +From: Tim Quatmann <tim.quatmann@cs.rwth-aachen.de> +Date: Wed, 12 Apr 2023 21:05:54 +0200 +Subject: [PATCH] stormeigen + +--- + Eigen/src/SparseLU/SparseLU.h | 8 +-- + Eigen/src/SparseLU/SparseLU_column_bmod.h | 2 +- + Eigen/src/SparseLU/SparseLU_copy_to_ucol.h | 2 +- + Eigen/src/SparseLU/SparseLU_panel_bmod.h | 6 +- + Eigen/src/SparseLU/SparseLU_pivotL.h | 79 +++++++++++++++------- + 5 files changed, 62 insertions(+), 35 deletions(-) + +diff --git a/Eigen/src/SparseLU/SparseLU.h b/Eigen/src/SparseLU/SparseLU.h +index f70aab197..33ba7e1d7 100644 +--- a/Eigen/src/SparseLU/SparseLU.h ++++ b/Eigen/src/SparseLU/SparseLU.h +@@ -158,12 +158,12 @@ class SparseLU : public SparseSolverBase<SparseLU<MatrixType_,OrderingType_> >, + + public: + +- SparseLU():m_lastError(""),m_Ustore(0,0,0,0,0,0),m_symmetricmode(false),m_diagpivotthresh(1.0),m_detPermR(1) ++ SparseLU():m_lastError(""),m_Ustore(0,0,0,0,0,0),m_symmetricmode(false),m_diagpivotthresh(1),m_detPermR(1) + { + initperfvalues(); + } + explicit SparseLU(const MatrixType& matrix) +- : m_lastError(""),m_Ustore(0,0,0,0,0,0),m_symmetricmode(false),m_diagpivotthresh(1.0),m_detPermR(1) ++ : m_lastError(""),m_Ustore(0,0,0,0,0,0),m_symmetricmode(false),m_diagpivotthresh(1),m_detPermR(1) + { + initperfvalues(); + compute(matrix); +@@ -355,7 +355,7 @@ class SparseLU : public SparseSolverBase<SparseLU<MatrixType_,OrderingType_> >, + using std::abs; + eigen_assert(m_factorizationIsOk && "The matrix should be factorized first."); + // Initialize with the determinant of the row matrix +- Scalar det = Scalar(1.); ++ Scalar det = Scalar(1); + // Note that the diagonal blocks of U are stored in supernodes, + // which are available in the L part :) + for (Index j = 0; j < this->cols(); ++j) +@@ -438,7 +438,7 @@ class SparseLU : public SparseSolverBase<SparseLU<MatrixType_,OrderingType_> >, + { + eigen_assert(m_factorizationIsOk && "The matrix should be factorized first."); + // Initialize with the determinant of the row matrix +- Scalar det = Scalar(1.); ++ Scalar det = Scalar(1); + // Note that the diagonal blocks of U are stored in supernodes, + // which are available in the L part :) + for (Index j = 0; j < this->cols(); ++j) +diff --git a/Eigen/src/SparseLU/SparseLU_column_bmod.h b/Eigen/src/SparseLU/SparseLU_column_bmod.h +index d5c29b353..09abdc402 100644 +--- a/Eigen/src/SparseLU/SparseLU_column_bmod.h ++++ b/Eigen/src/SparseLU/SparseLU_column_bmod.h +@@ -131,7 +131,7 @@ Index SparseLUImpl<Scalar,StorageIndex>::column_bmod(const Index jcol, const Ind + { + irow = glu.lsub(isub); + glu.lusup(nextlu) = dense(irow); +- dense(irow) = Scalar(0.0); ++ dense(irow) = Scalar(0); + ++nextlu; + } + +diff --git a/Eigen/src/SparseLU/SparseLU_copy_to_ucol.h b/Eigen/src/SparseLU/SparseLU_copy_to_ucol.h +index e06b2a03a..cf91f49df 100644 +--- a/Eigen/src/SparseLU/SparseLU_copy_to_ucol.h ++++ b/Eigen/src/SparseLU/SparseLU_copy_to_ucol.h +@@ -89,7 +89,7 @@ Index SparseLUImpl<Scalar,StorageIndex>::copy_to_ucol(const Index jcol, const In + irow = glu.lsub(isub); + glu.usub(nextu) = perm_r(irow); // Unlike the L part, the U part is stored in its final order + glu.ucol(nextu) = dense(irow); +- dense(irow) = Scalar(0.0); ++ dense(irow) = Scalar(0); + nextu++; + isub++; + } +diff --git a/Eigen/src/SparseLU/SparseLU_panel_bmod.h b/Eigen/src/SparseLU/SparseLU_panel_bmod.h +index fd1243d72..a34cbad6d 100644 +--- a/Eigen/src/SparseLU/SparseLU_panel_bmod.h ++++ b/Eigen/src/SparseLU/SparseLU_panel_bmod.h +@@ -124,7 +124,7 @@ void SparseLUImpl<Scalar,StorageIndex>::panel_bmod(const Index m, const Index w, + + Index isub = lptr + no_zeros; + Index off = u_rows-segsize; +- for (Index i = 0; i < off; i++) U(i,u_col) = 0; ++ for (Index i = 0; i < off; i++) U(i,u_col) = Scalar(0); + for (Index i = 0; i < segsize; i++) + { + Index irow = glu.lsub(isub); +@@ -173,7 +173,7 @@ void SparseLUImpl<Scalar,StorageIndex>::panel_bmod(const Index m, const Index w, + { + Index irow = glu.lsub(isub++); + dense_col(irow) = U.coeff(i+off,u_col); +- U.coeffRef(i+off,u_col) = 0; ++ U.coeffRef(i+off,u_col) = Scalar(0); + } + + // Scatter l into SPA dense[] +@@ -181,7 +181,7 @@ void SparseLUImpl<Scalar,StorageIndex>::panel_bmod(const Index m, const Index w, + { + Index irow = glu.lsub(isub++); + dense_col(irow) -= L.coeff(i,u_col); +- L.coeffRef(i,u_col) = 0; ++ L.coeffRef(i,u_col) = Scalar(0); + } + u_col++; + } +diff --git a/Eigen/src/SparseLU/SparseLU_pivotL.h b/Eigen/src/SparseLU/SparseLU_pivotL.h +index 6daed9133..78da06f62 100644 +--- a/Eigen/src/SparseLU/SparseLU_pivotL.h ++++ b/Eigen/src/SparseLU/SparseLU_pivotL.h +@@ -34,7 +34,46 @@ + + namespace Eigen { + namespace internal { +- ++ template<typename Scalar, typename StorageIndex, typename std::enable_if<!std::is_same<Scalar, storm::RationalFunction>::value, int>::type = 0> ++ void findLargestAbsolutePivotElement(const Index& nsupc, const Index& nsupr, Scalar const* lu_col_ptr, StorageIndex const* lsub_ptr, Index const& diagind, Index& diag, Index& pivptr, Scalar& pivmax) { ++ Scalar rtemp = 0; ++ for (Index isub = nsupc; isub < nsupr; ++isub) { ++ rtemp = storm::utility::abs(lu_col_ptr[isub]); ++ if (rtemp > pivmax) { ++ pivmax = rtemp; ++ pivptr = isub; ++ } ++ if (lsub_ptr[isub] == diagind) diag = isub; ++ } ++ } ++ ++ template<typename Scalar, typename StorageIndex, typename std::enable_if<std::is_same<Scalar, storm::RationalFunction>::value, int>::type = 0> ++ void findLargestAbsolutePivotElement(const Index& nsupc, const Index& nsupr, Scalar const* lu_col_ptr, StorageIndex const* lsub_ptr, Index const& diagind, Index& diag, Index& pivptr, Scalar& pivmax) { ++ bool foundPivotElement = false; ++ for (Index isub = nsupc; isub < nsupr; ++isub) { ++ Scalar const& rtemp = lu_col_ptr[isub]; ++ if (!foundPivotElement && !storm::utility::isZero(rtemp)) { ++ foundPivotElement = true; ++ pivmax = rtemp; ++ pivptr = isub; ++ } ++ if (lsub_ptr[isub] == diagind) diag = isub; ++ } ++ } ++ ++ template<typename Scalar, typename std::enable_if<std::is_same<Scalar, double>::value, int>::type = 0> ++ bool diagonalElementCanBePivot(Scalar const* lu_col_ptr, Index const& diag, Scalar const& diagpivotthresh, Scalar const& pivmax) { ++ Scalar thresh = diagpivotthresh * pivmax; ++ double rtemp = std::abs(lu_col_ptr[diag]); ++ if (!storm::utility::isZero(rtemp) && rtemp >= thresh) return true; ++ return false; ++ } ++ ++ template<typename Scalar, typename std::enable_if<!std::is_same<Scalar, double>::value, int>::type = 0> ++ bool diagonalElementCanBePivot(Scalar const* lu_col_ptr, Index const& diag, Scalar const& diagpivotthresh, Scalar const& pivmax) { ++ if (!storm::utility::isZero(lu_col_ptr[diag])) return true; ++ return false; ++ } + /** + * \brief Performs the numerical pivotin on the current column of L, and the CDIV operation. + * +@@ -72,42 +111,30 @@ Index SparseLUImpl<Scalar,StorageIndex>::pivotL(const Index jcol, const RealScal + StorageIndex* lsub_ptr = &(glu.lsub.data()[lptr]); // Start of row indices of the supernode + + // Determine the largest abs numerical value for partial pivoting +- Index diagind = iperm_c(jcol); // diagonal index +- RealScalar pivmax(-1.0); +- Index pivptr = nsupc; +- Index diag = emptyIdxLU; +- RealScalar rtemp; +- Index isub, icol, itemp, k; +- for (isub = nsupc; isub < nsupr; ++isub) { +- using std::abs; +- rtemp = abs(lu_col_ptr[isub]); +- if (rtemp > pivmax) { +- pivmax = rtemp; +- pivptr = isub; +- } +- if (lsub_ptr[isub] == diagind) diag = isub; +- } ++ Index diagind = iperm_c(jcol); // diagonal index ++ RealScalar pivmax(-1); ++ Index pivptr = nsupc; ++ Index diag = emptyIdxLU; ++ Index icol, itemp, k; ++ ++ findLargestAbsolutePivotElement(nsupc, nsupr, lu_col_ptr, lsub_ptr, diagind, diag, pivptr, pivmax); + + // Test for singularity +- if ( pivmax <= RealScalar(0.0) ) { ++ bool columnStructurallyEmpty = nsupr <= nsupc; ++ if ( columnStructurallyEmpty || storm::utility::isZero(pivmax)) { + // if pivmax == -1, the column is structurally empty, otherwise it is only numerically zero +- pivrow = pivmax < RealScalar(0.0) ? diagind : lsub_ptr[pivptr]; ++ pivrow = columnStructurallyEmpty ? diagind : lsub_ptr[pivptr]; + perm_r(pivrow) = StorageIndex(jcol); + return (jcol+1); + } + +- RealScalar thresh = diagpivotthresh * pivmax; +- +- // Choose appropriate pivotal element +- ++ // Choose appropriate pivotal element + { + // Test if the diagonal element can be used as a pivot (given the threshold value) + if (diag >= 0 ) + { + // Diagonal element exists +- using std::abs; +- rtemp = abs(lu_col_ptr[diag]); +- if (rtemp != RealScalar(0.0) && rtemp >= thresh) pivptr = diag; ++ if (diagonalElementCanBePivot(lu_col_ptr, diag, diagpivotthresh, pivmax)) pivptr = diag; + } + pivrow = lsub_ptr[pivptr]; + } +@@ -127,7 +154,7 @@ Index SparseLUImpl<Scalar,StorageIndex>::pivotL(const Index jcol, const RealScal + } + } + // cdiv operations +- Scalar temp = Scalar(1.0) / lu_col_ptr[nsupc]; ++ Scalar temp = Scalar(1) / lu_col_ptr[nsupc]; + for (k = nsupc+1; k < nsupr; k++) + lu_col_ptr[k] *= temp; + return 0; +-- +2.32.1 (Apple Git-133) +