From 153339c5be7ceae07420a16102e666125cc742d9 Mon Sep 17 00:00:00 2001 From: dehnert Date: Tue, 28 Mar 2017 20:57:17 +0200 Subject: [PATCH] first draft of policy iteration using DDs --- resources/3rdparty/CMakeLists.txt | 3 + resources/3rdparty/cudd-3.0.0/Makefile.in | 1 - resources/3rdparty/cudd-3.0.0/configure | 14 +- resources/3rdparty/sylvan/src/sylvan_int.h | 2 + .../3rdparty/sylvan/src/sylvan_mtbdd_storm.c | 96 +++-- .../sylvan/src/sylvan_obj_mtbdd_storm.hpp | 3 + .../3rdparty/sylvan/src/sylvan_obj_storm.cpp | 21 +- .../sylvan/src/sylvan_storm_rational_number.c | 391 ++++++++++++++++++ .../sylvan/src/sylvan_storm_rational_number.h | 12 +- .../prctl/SymbolicMdpPrctlModelChecker.cpp | 4 +- .../prctl/SymbolicMdpPrctlModelChecker.h | 6 +- .../prctl/helper/SymbolicMdpPrctlHelper.cpp | 12 +- .../prctl/helper/SymbolicMdpPrctlHelper.h | 14 +- .../SymbolicQuantitativeCheckResult.cpp | 68 ++- ...ymbolicEliminationLinearEquationSolver.cpp | 35 +- .../SymbolicEliminationLinearEquationSolver.h | 12 + .../solver/SymbolicLinearEquationSolver.cpp | 4 + .../solver/SymbolicLinearEquationSolver.h | 2 + .../SymbolicMinMaxLinearEquationSolver.cpp | 183 +++++++- .../SymbolicMinMaxLinearEquationSolver.h | 87 ++-- .../SymbolicNativeLinearEquationSolver.cpp | 4 +- .../storage/dd/sylvan/InternalSylvanAdd.cpp | 10 + src/storm/utility/solver.cpp | 10 - src/storm/utility/solver.h | 6 - .../SymbolicMdpPrctlModelCheckerTest.cpp | 8 +- 25 files changed, 833 insertions(+), 175 deletions(-) diff --git a/resources/3rdparty/CMakeLists.txt b/resources/3rdparty/CMakeLists.txt index b62cfb3d0..3e6630c75 100644 --- a/resources/3rdparty/CMakeLists.txt +++ b/resources/3rdparty/CMakeLists.txt @@ -383,6 +383,9 @@ ExternalProject_Add( LOG_CONFIGURE ON LOG_BUILD ON BUILD_BYPRODUCTS ${STORM_3RDPARTY_BINARY_DIR}/sylvan/src/libsylvan${STATIC_EXT} + + # TODO: remove + BUILD_ALWAYS 1 ) ExternalProject_Get_Property(sylvan source_dir) ExternalProject_Get_Property(sylvan binary_dir) diff --git a/resources/3rdparty/cudd-3.0.0/Makefile.in b/resources/3rdparty/cudd-3.0.0/Makefile.in index c18a9a946..75846bc4c 100644 --- a/resources/3rdparty/cudd-3.0.0/Makefile.in +++ b/resources/3rdparty/cudd-3.0.0/Makefile.in @@ -795,7 +795,6 @@ pdfdir = @pdfdir@ prefix = @prefix@ program_transform_name = @program_transform_name@ psdir = @psdir@ -runstatedir = @runstatedir@ sbindir = @sbindir@ sharedstatedir = @sharedstatedir@ srcdir = @srcdir@ diff --git a/resources/3rdparty/cudd-3.0.0/configure b/resources/3rdparty/cudd-3.0.0/configure index ef3832be0..1494802d5 100755 --- a/resources/3rdparty/cudd-3.0.0/configure +++ b/resources/3rdparty/cudd-3.0.0/configure @@ -754,7 +754,6 @@ infodir docdir oldincludedir includedir -runstatedir localstatedir sharedstatedir sysconfdir @@ -841,7 +840,6 @@ datadir='${datarootdir}' sysconfdir='${prefix}/etc' sharedstatedir='${prefix}/com' localstatedir='${prefix}/var' -runstatedir='${localstatedir}/run' includedir='${prefix}/include' oldincludedir='/usr/include' docdir='${datarootdir}/doc/${PACKAGE_TARNAME}' @@ -1094,15 +1092,6 @@ do | -silent | --silent | --silen | --sile | --sil) silent=yes ;; - -runstatedir | --runstatedir | --runstatedi | --runstated \ - | --runstate | --runstat | --runsta | --runst | --runs \ - | --run | --ru | --r) - ac_prev=runstatedir ;; - -runstatedir=* | --runstatedir=* | --runstatedi=* | --runstated=* \ - | --runstate=* | --runstat=* | --runsta=* | --runst=* | --runs=* \ - | --run=* | --ru=* | --r=*) - runstatedir=$ac_optarg ;; - -sbindir | --sbindir | --sbindi | --sbind | --sbin | --sbi | --sb) ac_prev=sbindir ;; -sbindir=* | --sbindir=* | --sbindi=* | --sbind=* | --sbin=* \ @@ -1240,7 +1229,7 @@ fi for ac_var in exec_prefix prefix bindir sbindir libexecdir datarootdir \ datadir sysconfdir sharedstatedir localstatedir includedir \ oldincludedir docdir infodir htmldir dvidir pdfdir psdir \ - libdir localedir mandir runstatedir + libdir localedir mandir do eval ac_val=\$$ac_var # Remove trailing slashes. @@ -1393,7 +1382,6 @@ Fine tuning of the installation directories: --sysconfdir=DIR read-only single-machine data [PREFIX/etc] --sharedstatedir=DIR modifiable architecture-independent data [PREFIX/com] --localstatedir=DIR modifiable single-machine data [PREFIX/var] - --runstatedir=DIR modifiable per-process data [LOCALSTATEDIR/run] --libdir=DIR object code libraries [EPREFIX/lib] --includedir=DIR C header files [PREFIX/include] --oldincludedir=DIR C header files for non-gcc [/usr/include] diff --git a/resources/3rdparty/sylvan/src/sylvan_int.h b/resources/3rdparty/sylvan/src/sylvan_int.h index df12c14f2..e8a96df81 100755 --- a/resources/3rdparty/sylvan/src/sylvan_int.h +++ b/resources/3rdparty/sylvan/src/sylvan_int.h @@ -103,6 +103,8 @@ extern llmsset_t nodes; #define CACHE_MTBDD_EQUAL_NORM_RF (66LL<<40) #define CACHE_MTBDD_EQUAL_NORM_REL_RF (67LL<<40) +#define CACHE_MTBDD_ABSTRACT_REPRESENTATIVE (68LL<<40) + #ifdef __cplusplus } #endif /* __cplusplus */ diff --git a/resources/3rdparty/sylvan/src/sylvan_mtbdd_storm.c b/resources/3rdparty/sylvan/src/sylvan_mtbdd_storm.c index 28d8e4e11..1915886f1 100644 --- a/resources/3rdparty/sylvan/src/sylvan_mtbdd_storm.c +++ b/resources/3rdparty/sylvan/src/sylvan_mtbdd_storm.c @@ -586,24 +586,24 @@ TASK_IMPL_2(MTBDD, mtbdd_op_complement, MTBDD, a, size_t, k) (void)k; // unused variable } -TASK_IMPL_3(BDD, mtbdd_min_abstract_representative, MTBDD, a, BDD, variables, BDDVAR, prev_level) { +TASK_IMPL_3(BDD, mtbdd_min_abstract_representative, MTBDD, a, BDD, v, BDDVAR, prev_level) { /* Maybe perform garbage collection */ sylvan_gc_test(); - if (sylvan_set_isempty(variables)) { + if (sylvan_set_isempty(v)) { return sylvan_true; } /* Cube is guaranteed to be a cube at this point. */ if (mtbdd_isleaf(a)) { - BDD _v = sylvan_set_next(variables); + BDD _v = sylvan_set_next(v); BDD res = CALL(mtbdd_min_abstract_representative, a, _v, prev_level); if (res == sylvan_invalid) { return sylvan_invalid; } sylvan_ref(res); - BDD res1 = sylvan_not(sylvan_ite(sylvan_ithvar(bddnode_getvariable(MTBDD_GETNODE(variables))), sylvan_true, sylvan_not(res))); + BDD res1 = sylvan_not(sylvan_ite(sylvan_ithvar(bddnode_getvariable(MTBDD_GETNODE(v))), sylvan_true, sylvan_not(res))); if (res1 == sylvan_invalid) { sylvan_deref(res); return sylvan_invalid; @@ -614,12 +614,12 @@ TASK_IMPL_3(BDD, mtbdd_min_abstract_representative, MTBDD, a, BDD, variables, BD mtbddnode_t na = MTBDD_GETNODE(a); uint32_t va = mtbddnode_getvariable(na); - bddnode_t nv = MTBDD_GETNODE(variables); + bddnode_t nv = MTBDD_GETNODE(v); BDDVAR vv = bddnode_getvariable(nv); /* Abstract a variable that does not appear in a. */ if (va > vv) { - BDD _v = sylvan_set_next(variables); + BDD _v = sylvan_set_next(v); BDD res = CALL(mtbdd_min_abstract_representative, a, _v, va); if (res == sylvan_invalid) { return sylvan_invalid; @@ -636,18 +636,19 @@ TASK_IMPL_3(BDD, mtbdd_min_abstract_representative, MTBDD, a, BDD, variables, BD return res1; } - /* TODO: Caching here. */ - /*if ((res = cuddCacheLookup2(manager, Cudd_addMinAbstractRepresentative, f, cube)) != NULL) { - return(res); - }*/ - + /* Check cache */ + MTBDD result; + if (cache_get3(CACHE_MTBDD_ABSTRACT_REPRESENTATIVE, a, v, (size_t)1, &result)) { + sylvan_stats_count(MTBDD_ABSTRACT_CACHED); + return result; + } MTBDD E = mtbdd_getlow(a); MTBDD T = mtbdd_gethigh(a); /* If the two indices are the same, so are their levels. */ if (va == vv) { - BDD _v = sylvan_set_next(variables); + BDD _v = sylvan_set_next(v); BDD res1 = CALL(mtbdd_min_abstract_representative, E, _v, va); if (res1 == sylvan_invalid) { return sylvan_invalid; @@ -723,19 +724,21 @@ TASK_IMPL_3(BDD, mtbdd_min_abstract_representative, MTBDD, a, BDD, variables, BD sylvan_deref(res1Inf); sylvan_deref(res2Inf); - /* TODO: Caching here. */ - //cuddCacheInsert2(manager, Cudd_addMinAbstractRepresentative, f, cube, res); - + /* Store in cache */ + if (cache_put3(CACHE_MTBDD_ABSTRACT_REPRESENTATIVE, a, v, (size_t)1, res)) { + sylvan_stats_count(MTBDD_ABSTRACT_CACHEDPUT); + } + sylvan_deref(res); return res; } else { /* if (va < vv) */ - BDD res1 = CALL(mtbdd_min_abstract_representative, E, variables, va); + BDD res1 = CALL(mtbdd_min_abstract_representative, E, v, va); if (res1 == sylvan_invalid) { return sylvan_invalid; } sylvan_ref(res1); - BDD res2 = CALL(mtbdd_min_abstract_representative, T, variables, va); + BDD res2 = CALL(mtbdd_min_abstract_representative, T, v, va); if (res2 == sylvan_invalid) { sylvan_deref(res1); return sylvan_invalid; @@ -750,46 +753,54 @@ TASK_IMPL_3(BDD, mtbdd_min_abstract_representative, MTBDD, a, BDD, variables, BD } sylvan_deref(res1); sylvan_deref(res2); - /* TODO: Caching here. */ - //cuddCacheInsert2(manager, Cudd_addMinAbstractRepresentative, f, cube, res); + + /* Store in cache */ + if (cache_put3(CACHE_MTBDD_ABSTRACT_REPRESENTATIVE, a, v, (size_t)1, res)) { + sylvan_stats_count(MTBDD_ABSTRACT_CACHEDPUT); + } return res; } } -TASK_IMPL_3(BDD, mtbdd_max_abstract_representative, MTBDD, a, MTBDD, variables, uint32_t, prev_level) { +TASK_IMPL_3(BDD, mtbdd_max_abstract_representative, MTBDD, a, MTBDD, v, uint32_t, prev_level) { /* Maybe perform garbage collection */ sylvan_gc_test(); - if (sylvan_set_isempty(variables)) { + if (sylvan_set_isempty(v)) { return sylvan_true; } + /* Count operation */ + sylvan_stats_count(MTBDD_ABSTRACT); + /* Cube is guaranteed to be a cube at this point. */ if (mtbdd_isleaf(a)) { - BDD _v = sylvan_set_next(variables); + /* Compute result */ + BDD _v = sylvan_set_next(v); BDD res = CALL(mtbdd_max_abstract_representative, a, _v, prev_level); if (res == sylvan_invalid) { return sylvan_invalid; } sylvan_ref(res); - BDD res1 = sylvan_not(sylvan_ite(sylvan_ithvar(bddnode_getvariable(MTBDD_GETNODE(variables))), sylvan_true, sylvan_not(res))); + BDD res1 = sylvan_not(sylvan_ite(sylvan_ithvar(bddnode_getvariable(MTBDD_GETNODE(v))), sylvan_true, sylvan_not(res))); if (res1 == sylvan_invalid) { sylvan_deref(res); return sylvan_invalid; } sylvan_deref(res); + return res1; } mtbddnode_t na = MTBDD_GETNODE(a); uint32_t va = mtbddnode_getvariable(na); - bddnode_t nv = MTBDD_GETNODE(variables); + bddnode_t nv = MTBDD_GETNODE(v); BDDVAR vv = bddnode_getvariable(nv); /* Abstract a variable that does not appear in a. */ - if (va > vv) { - BDD _v = sylvan_set_next(variables); + if (vv < va) { + BDD _v = sylvan_set_next(v); BDD res = CALL(mtbdd_max_abstract_representative, a, _v, va); if (res == sylvan_invalid) { return sylvan_invalid; @@ -806,18 +817,19 @@ TASK_IMPL_3(BDD, mtbdd_max_abstract_representative, MTBDD, a, MTBDD, variables, return res1; } - /* TODO: Caching here. */ - /*if ((res = cuddCacheLookup2(manager, Cudd_addMinAbstractRepresentative, f, cube)) != NULL) { - return(res); - }*/ - + /* Check cache */ + MTBDD result; + if (cache_get3(CACHE_MTBDD_ABSTRACT_REPRESENTATIVE, a, v, (size_t)0, &result)) { + sylvan_stats_count(MTBDD_ABSTRACT_CACHED); + return result; + } MTBDD E = mtbdd_getlow(a); MTBDD T = mtbdd_gethigh(a); /* If the two indices are the same, so are their levels. */ if (va == vv) { - BDD _v = sylvan_set_next(variables); + BDD _v = sylvan_set_next(v); BDD res1 = CALL(mtbdd_max_abstract_representative, E, _v, va); if (res1 == sylvan_invalid) { return sylvan_invalid; @@ -893,19 +905,21 @@ TASK_IMPL_3(BDD, mtbdd_max_abstract_representative, MTBDD, a, MTBDD, variables, sylvan_deref(res1Inf); sylvan_deref(res2Inf); - /* TODO: Caching here. */ - //cuddCacheInsert2(manager, Cudd_addMinAbstractRepresentative, f, cube, res); - + /* Store in cache */ + if (cache_put3(CACHE_MTBDD_ABSTRACT_REPRESENTATIVE, a, v, (size_t)0, res)) { + sylvan_stats_count(MTBDD_ABSTRACT_CACHEDPUT); + } + sylvan_deref(res); return res; } else { /* if (va < vv) */ - BDD res1 = CALL(mtbdd_max_abstract_representative, E, variables, va); + BDD res1 = CALL(mtbdd_max_abstract_representative, E, v, va); if (res1 == sylvan_invalid) { return sylvan_invalid; } sylvan_ref(res1); - BDD res2 = CALL(mtbdd_max_abstract_representative, T, variables, va); + BDD res2 = CALL(mtbdd_max_abstract_representative, T, v, va); if (res2 == sylvan_invalid) { sylvan_deref(res1); return sylvan_invalid; @@ -920,8 +934,12 @@ TASK_IMPL_3(BDD, mtbdd_max_abstract_representative, MTBDD, a, MTBDD, variables, } sylvan_deref(res1); sylvan_deref(res2); - /* TODO: Caching here. */ - //cuddCacheInsert2(manager, Cudd_addMinAbstractRepresentative, f, cube, res); + + /* Store in cache */ + if (cache_put3(CACHE_MTBDD_ABSTRACT_REPRESENTATIVE, a, v, (size_t)0, res)) { + sylvan_stats_count(MTBDD_ABSTRACT_CACHEDPUT); + } + return res; } } diff --git a/resources/3rdparty/sylvan/src/sylvan_obj_mtbdd_storm.hpp b/resources/3rdparty/sylvan/src/sylvan_obj_mtbdd_storm.hpp index 59dda64fd..23ae08b1e 100644 --- a/resources/3rdparty/sylvan/src/sylvan_obj_mtbdd_storm.hpp +++ b/resources/3rdparty/sylvan/src/sylvan_obj_mtbdd_storm.hpp @@ -56,6 +56,9 @@ Mtbdd AbstractPlusRN(const BddSet &variables) const; Mtbdd AbstractMinRN(const BddSet &variables) const; Mtbdd AbstractMaxRN(const BddSet &variables) const; +Bdd AbstractMinRepresentativeRN(const BddSet &variables) const; +Bdd AbstractMaxRepresentativeRN(const BddSet &variables) const; + Bdd BddThresholdRN(storm::RationalNumber const& rn) const; Bdd BddStrictThresholdRN(storm::RationalNumber const& rn) const; diff --git a/resources/3rdparty/sylvan/src/sylvan_obj_storm.cpp b/resources/3rdparty/sylvan/src/sylvan_obj_storm.cpp index 6e3b524c7..86d2afc85 100644 --- a/resources/3rdparty/sylvan/src/sylvan_obj_storm.cpp +++ b/resources/3rdparty/sylvan/src/sylvan_obj_storm.cpp @@ -298,21 +298,36 @@ Mtbdd::AndExistsRN(const Mtbdd &other, const BddSet &variables) const { return sylvan_storm_rational_number_and_exists(mtbdd, other.mtbdd, variables.set.bdd); } -Mtbdd Mtbdd::AbstractPlusRN(const BddSet &variables) const { +Mtbdd +Mtbdd::AbstractPlusRN(const BddSet &variables) const { LACE_ME; return sylvan_storm_rational_number_abstract_plus(mtbdd, variables.set.bdd); } -Mtbdd Mtbdd::AbstractMinRN(const BddSet &variables) const { +Mtbdd +Mtbdd::AbstractMinRN(const BddSet &variables) const { LACE_ME; return sylvan_storm_rational_number_abstract_min(mtbdd, variables.set.bdd); } -Mtbdd Mtbdd::AbstractMaxRN(const BddSet &variables) const { +Mtbdd +Mtbdd::AbstractMaxRN(const BddSet &variables) const { LACE_ME; return sylvan_storm_rational_number_abstract_max(mtbdd, variables.set.bdd); } +Bdd +Mtbdd::AbstractMinRepresentativeRN(const BddSet &variables) const { + LACE_ME; + return sylvan_storm_rational_number_min_abstract_representative(mtbdd, variables.set.bdd); +} + +Bdd +Mtbdd::AbstractMaxRepresentativeRN(const BddSet &variables) const { + LACE_ME; + return sylvan_storm_rational_number_max_abstract_representative(mtbdd, variables.set.bdd); +} + Bdd Mtbdd::BddThresholdRN(storm::RationalNumber const& rn) const { LACE_ME; diff --git a/resources/3rdparty/sylvan/src/sylvan_storm_rational_number.c b/resources/3rdparty/sylvan/src/sylvan_storm_rational_number.c index ef48ccab3..5008607bf 100644 --- a/resources/3rdparty/sylvan/src/sylvan_storm_rational_number.c +++ b/resources/3rdparty/sylvan/src/sylvan_storm_rational_number.c @@ -227,6 +227,23 @@ TASK_IMPL_2(MTBDD, sylvan_storm_rational_number_op_less, MTBDD*, pa, MTBDD*, pb) return mtbdd_invalid; } +TASK_IMPL_2(MTBDD, sylvan_storm_rational_number_op_greater, MTBDD*, pa, MTBDD*, pb) { + MTBDD a = *pa, b = *pb; + + /* Check for partial functions and for Boolean (filter) */ + if (a == mtbdd_false || b == mtbdd_false) return mtbdd_false; + + /* If both leaves, compute less */ + if (mtbdd_isleaf(a) && mtbdd_isleaf(b)) { + storm_rational_number_ptr ma = mtbdd_getstorm_rational_number_ptr(a); + storm_rational_number_ptr mb = mtbdd_getstorm_rational_number_ptr(b); + + return storm_rational_number_less_or_equal(ma, mb) ? mtbdd_false : mtbdd_true; + } + + return mtbdd_invalid; +} + TASK_IMPL_2(MTBDD, sylvan_storm_rational_number_op_less_or_equal, MTBDD*, pa, MTBDD*, pb) { MTBDD a = *pa, b = *pb; @@ -244,6 +261,23 @@ TASK_IMPL_2(MTBDD, sylvan_storm_rational_number_op_less_or_equal, MTBDD*, pa, MT return mtbdd_invalid; } +TASK_IMPL_2(MTBDD, sylvan_storm_rational_number_op_greater_or_equal, MTBDD*, pa, MTBDD*, pb) { + MTBDD a = *pa, b = *pb; + + /* Check for partial functions and for Boolean (filter) */ + if (a == mtbdd_false || b == mtbdd_false) return mtbdd_false; + + /* If both leaves, compute less or equal */ + if (mtbdd_isleaf(a) && mtbdd_isleaf(b)) { + storm_rational_number_ptr ma = mtbdd_getstorm_rational_number_ptr(a); + storm_rational_number_ptr mb = mtbdd_getstorm_rational_number_ptr(b); + + return storm_rational_number_less(ma, mb) ? mtbdd_false : mtbdd_true; + } + + return mtbdd_invalid; +} + TASK_IMPL_2(MTBDD, sylvan_storm_rational_number_op_mod, MTBDD*, pa, MTBDD*, pb) { MTBDD a = *pa, b = *pb; @@ -788,3 +822,360 @@ TASK_IMPL_3(MTBDD, sylvan_storm_rational_number_equal_norm_rel_d, MTBDD, a, MTBD return CALL(sylvan_storm_rational_number_equal_norm_rel_d2, a, b, d, &shortcircuit); } +TASK_IMPL_3(BDD, sylvan_storm_rational_number_min_abstract_representative, MTBDD, a, BDD, v, BDDVAR, prev_level) { + /* Maybe perform garbage collection */ + sylvan_gc_test(); + + if (sylvan_set_isempty(v)) { + return sylvan_true; + } + + /* Cube is guaranteed to be a cube at this point. */ + if (mtbdd_isleaf(a)) { + BDD _v = sylvan_set_next(v); + BDD res = CALL(sylvan_storm_rational_number_min_abstract_representative, a, _v, prev_level); + if (res == sylvan_invalid) { + return sylvan_invalid; + } + sylvan_ref(res); + + BDD res1 = sylvan_not(sylvan_ite(sylvan_ithvar(bddnode_getvariable(MTBDD_GETNODE(v))), sylvan_true, sylvan_not(res))); + if (res1 == sylvan_invalid) { + sylvan_deref(res); + return sylvan_invalid; + } + sylvan_deref(res); + return res1; + } + + mtbddnode_t na = MTBDD_GETNODE(a); + uint32_t va = mtbddnode_getvariable(na); + bddnode_t nv = MTBDD_GETNODE(v); + BDDVAR vv = bddnode_getvariable(nv); + + /* Abstract a variable that does not appear in a. */ + if (va > vv) { + BDD _v = sylvan_set_next(v); + BDD res = CALL(sylvan_storm_rational_number_min_abstract_representative, a, _v, va); + if (res == sylvan_invalid) { + return sylvan_invalid; + } + + // Fill in the missing variables to make representative unique. + sylvan_ref(res); + BDD res1 = sylvan_ite(sylvan_ithvar(vv), sylvan_false, res); + if (res1 == sylvan_invalid) { + sylvan_deref(res); + return sylvan_invalid; + } + sylvan_deref(res); + return res1; + } + + /* Check cache */ + MTBDD result; + if (cache_get3(CACHE_MTBDD_ABSTRACT_REPRESENTATIVE, a, v, (size_t)1, &result)) { + sylvan_stats_count(MTBDD_ABSTRACT_CACHED); + return result; + } + + MTBDD E = mtbdd_getlow(a); + MTBDD T = mtbdd_gethigh(a); + + /* If the two indices are the same, so are their levels. */ + if (va == vv) { + BDD _v = sylvan_set_next(v); + BDD res1 = CALL(sylvan_storm_rational_number_min_abstract_representative, E, _v, va); + if (res1 == sylvan_invalid) { + return sylvan_invalid; + } + sylvan_ref(res1); + + BDD res2 = CALL(sylvan_storm_rational_number_min_abstract_representative, T, _v, va); + if (res2 == sylvan_invalid) { + sylvan_deref(res1); + return sylvan_invalid; + } + sylvan_ref(res2); + + MTBDD left = sylvan_storm_rational_number_abstract_min(E, _v); + if (left == mtbdd_invalid) { + sylvan_deref(res1); + sylvan_deref(res2); + return sylvan_invalid; + } + mtbdd_ref(left); + + MTBDD right = sylvan_storm_rational_number_abstract_min(T, _v); + if (right == mtbdd_invalid) { + sylvan_deref(res1); + sylvan_deref(res2); + mtbdd_deref(left); + return sylvan_invalid; + } + mtbdd_ref(right); + + BDD tmp = sylvan_storm_rational_number_less_or_equal(left, right); + if (tmp == sylvan_invalid) { + sylvan_deref(res1); + sylvan_deref(res2); + mtbdd_deref(left); + mtbdd_deref(right); + return sylvan_invalid; + } + sylvan_ref(tmp); + + mtbdd_deref(left); + mtbdd_deref(right); + + BDD res1Inf = sylvan_ite(tmp, res1, sylvan_false); + if (res1Inf == sylvan_invalid) { + sylvan_deref(res1); + sylvan_deref(res2); + sylvan_deref(tmp); + return sylvan_invalid; + } + sylvan_ref(res1Inf); + sylvan_deref(res1); + + BDD res2Inf = sylvan_ite(tmp, sylvan_false, res2); + if (res2Inf == sylvan_invalid) { + sylvan_deref(res2); + sylvan_deref(res1Inf); + sylvan_deref(tmp); + return sylvan_invalid; + } + sylvan_ref(res2Inf); + sylvan_deref(res2); + sylvan_deref(tmp); + + BDD res = (res1Inf == res2Inf) ? sylvan_ite(sylvan_ithvar(va), sylvan_false, res1Inf) : sylvan_ite(sylvan_ithvar(va), res2Inf, res1Inf); + + if (res == sylvan_invalid) { + sylvan_deref(res1Inf); + sylvan_deref(res2Inf); + return sylvan_invalid; + } + sylvan_ref(res); + sylvan_deref(res1Inf); + sylvan_deref(res2Inf); + + /* Store in cache */ + if (cache_put3(CACHE_MTBDD_ABSTRACT_REPRESENTATIVE, a, v, (size_t)1, res)) { + sylvan_stats_count(MTBDD_ABSTRACT_CACHEDPUT); + } + + sylvan_deref(res); + return res; + } + else { /* if (va < vv) */ + BDD res1 = CALL(sylvan_storm_rational_number_min_abstract_representative, E, v, va); + if (res1 == sylvan_invalid) { + return sylvan_invalid; + } + sylvan_ref(res1); + BDD res2 = CALL(sylvan_storm_rational_number_min_abstract_representative, T, v, va); + if (res2 == sylvan_invalid) { + sylvan_deref(res1); + return sylvan_invalid; + } + sylvan_ref(res2); + + BDD res = (res1 == res2) ? res1 : sylvan_ite(sylvan_ithvar(va), res2, res1); + if (res == sylvan_invalid) { + sylvan_deref(res1); + sylvan_deref(res2); + return sylvan_invalid; + } + sylvan_deref(res1); + sylvan_deref(res2); + + /* Store in cache */ + if (cache_put3(CACHE_MTBDD_ABSTRACT_REPRESENTATIVE, a, v, (size_t)1, res)) { + sylvan_stats_count(MTBDD_ABSTRACT_CACHEDPUT); + } + return res; + } +} + +TASK_IMPL_3(BDD, sylvan_storm_rational_number_max_abstract_representative, MTBDD, a, MTBDD, v, uint32_t, prev_level) { + /* Maybe perform garbage collection */ + sylvan_gc_test(); + + if (sylvan_set_isempty(v)) { + return sylvan_true; + } + + /* Count operation */ + sylvan_stats_count(MTBDD_ABSTRACT); + + /* Cube is guaranteed to be a cube at this point. */ + if (mtbdd_isleaf(a)) { + /* Compute result */ + BDD _v = sylvan_set_next(v); + BDD res = CALL(sylvan_storm_rational_number_max_abstract_representative, a, _v, prev_level); + if (res == sylvan_invalid) { + return sylvan_invalid; + } + sylvan_ref(res); + + BDD res1 = sylvan_not(sylvan_ite(sylvan_ithvar(bddnode_getvariable(MTBDD_GETNODE(v))), sylvan_true, sylvan_not(res))); + if (res1 == sylvan_invalid) { + sylvan_deref(res); + return sylvan_invalid; + } + sylvan_deref(res); + + return res1; + } + + mtbddnode_t na = MTBDD_GETNODE(a); + uint32_t va = mtbddnode_getvariable(na); + bddnode_t nv = MTBDD_GETNODE(v); + BDDVAR vv = bddnode_getvariable(nv); + + /* Abstract a variable that does not appear in a. */ + if (vv < va) { + BDD _v = sylvan_set_next(v); + BDD res = CALL(sylvan_storm_rational_number_max_abstract_representative, a, _v, va); + if (res == sylvan_invalid) { + return sylvan_invalid; + } + + // Fill in the missing variables to make representative unique. + sylvan_ref(res); + BDD res1 = sylvan_ite(sylvan_ithvar(vv), sylvan_false, res); + if (res1 == sylvan_invalid) { + sylvan_deref(res); + return sylvan_invalid; + } + sylvan_deref(res); + return res1; + } + + /* Check cache */ + MTBDD result; + if (cache_get3(CACHE_MTBDD_ABSTRACT_REPRESENTATIVE, a, v, (size_t)0, &result)) { + sylvan_stats_count(MTBDD_ABSTRACT_CACHED); + return result; + } + + MTBDD E = mtbdd_getlow(a); + MTBDD T = mtbdd_gethigh(a); + + /* If the two indices are the same, so are their levels. */ + if (va == vv) { + BDD _v = sylvan_set_next(v); + BDD res1 = CALL(sylvan_storm_rational_number_max_abstract_representative, E, _v, va); + if (res1 == sylvan_invalid) { + return sylvan_invalid; + } + sylvan_ref(res1); + + BDD res2 = CALL(sylvan_storm_rational_number_max_abstract_representative, T, _v, va); + if (res2 == sylvan_invalid) { + sylvan_deref(res1); + return sylvan_invalid; + } + sylvan_ref(res2); + + MTBDD left = sylvan_storm_rational_number_abstract_max(E, _v); + if (left == mtbdd_invalid) { + sylvan_deref(res1); + sylvan_deref(res2); + return sylvan_invalid; + } + mtbdd_ref(left); + + MTBDD right = sylvan_storm_rational_number_abstract_max(T, _v); + if (right == mtbdd_invalid) { + sylvan_deref(res1); + sylvan_deref(res2); + mtbdd_deref(left); + return sylvan_invalid; + } + mtbdd_ref(right); + + BDD tmp = sylvan_storm_rational_number_greater_or_equal(left, right); + if (tmp == sylvan_invalid) { + sylvan_deref(res1); + sylvan_deref(res2); + mtbdd_deref(left); + mtbdd_deref(right); + return sylvan_invalid; + } + sylvan_ref(tmp); + + mtbdd_deref(left); + mtbdd_deref(right); + + BDD res1Inf = sylvan_ite(tmp, res1, sylvan_false); + if (res1Inf == sylvan_invalid) { + sylvan_deref(res1); + sylvan_deref(res2); + sylvan_deref(tmp); + return sylvan_invalid; + } + sylvan_ref(res1Inf); + sylvan_deref(res1); + + BDD res2Inf = sylvan_ite(tmp, sylvan_false, res2); + if (res2Inf == sylvan_invalid) { + sylvan_deref(res2); + sylvan_deref(res1Inf); + sylvan_deref(tmp); + return sylvan_invalid; + } + sylvan_ref(res2Inf); + sylvan_deref(res2); + sylvan_deref(tmp); + + BDD res = (res1Inf == res2Inf) ? sylvan_ite(sylvan_ithvar(va), sylvan_false, res1Inf) : sylvan_ite(sylvan_ithvar(va), res2Inf, res1Inf); + + if (res == sylvan_invalid) { + sylvan_deref(res1Inf); + sylvan_deref(res2Inf); + return sylvan_invalid; + } + sylvan_ref(res); + sylvan_deref(res1Inf); + sylvan_deref(res2Inf); + + /* Store in cache */ + if (cache_put3(CACHE_MTBDD_ABSTRACT_REPRESENTATIVE, a, v, (size_t)0, res)) { + sylvan_stats_count(MTBDD_ABSTRACT_CACHEDPUT); + } + + sylvan_deref(res); + return res; + } + else { /* if (va < vv) */ + BDD res1 = CALL(sylvan_storm_rational_number_max_abstract_representative, E, v, va); + if (res1 == sylvan_invalid) { + return sylvan_invalid; + } + sylvan_ref(res1); + BDD res2 = CALL(sylvan_storm_rational_number_max_abstract_representative, T, v, va); + if (res2 == sylvan_invalid) { + sylvan_deref(res1); + return sylvan_invalid; + } + sylvan_ref(res2); + + BDD res = (res1 == res2) ? res1 : sylvan_ite(sylvan_ithvar(va), res2, res1); + if (res == sylvan_invalid) { + sylvan_deref(res1); + sylvan_deref(res2); + return sylvan_invalid; + } + sylvan_deref(res1); + sylvan_deref(res2); + + /* Store in cache */ + if (cache_put3(CACHE_MTBDD_ABSTRACT_REPRESENTATIVE, a, v, (size_t)0, res)) { + sylvan_stats_count(MTBDD_ABSTRACT_CACHEDPUT); + } + + return res; + } +} diff --git a/resources/3rdparty/sylvan/src/sylvan_storm_rational_number.h b/resources/3rdparty/sylvan/src/sylvan_storm_rational_number.h index 1425e1bfc..d7accd375 100644 --- a/resources/3rdparty/sylvan/src/sylvan_storm_rational_number.h +++ b/resources/3rdparty/sylvan/src/sylvan_storm_rational_number.h @@ -31,7 +31,9 @@ TASK_DECL_3(MTBDD, sylvan_storm_rational_number_abstract_op_min, MTBDD, MTBDD, i TASK_DECL_3(MTBDD, sylvan_storm_rational_number_abstract_op_max, MTBDD, MTBDD, int) TASK_DECL_2(MTBDD, sylvan_storm_rational_number_op_less, MTBDD*, MTBDD*) +TASK_DECL_2(MTBDD, sylvan_storm_rational_number_op_greater, MTBDD*, MTBDD*) TASK_DECL_2(MTBDD, sylvan_storm_rational_number_op_less_or_equal, MTBDD*, MTBDD*) +TASK_DECL_2(MTBDD, sylvan_storm_rational_number_op_greater_or_equal, MTBDD*, MTBDD*) TASK_DECL_2(MTBDD, sylvan_storm_rational_number_op_mod, MTBDD*, MTBDD*) TASK_DECL_2(MTBDD, sylvan_storm_rational_number_op_pow, MTBDD*, MTBDD*) @@ -47,8 +49,10 @@ TASK_DECL_2(MTBDD, sylvan_storm_rational_number_op_ceil, MTBDD, size_t) #define sylvan_storm_rational_number_times(a, b) mtbdd_apply(a, b, TASK(sylvan_storm_rational_number_op_times)) #define sylvan_storm_rational_number_divide(a, b) mtbdd_apply(a, b, TASK(sylvan_storm_rational_number_op_divide)) #define sylvan_storm_rational_number_less(a, b) mtbdd_apply(a, b, TASK(sylvan_storm_rational_number_op_less)) +#define sylvan_storm_rational_number_greater(a, b) mtbdd_apply(a, b, TASK(sylvan_storm_rational_number_op_greater)) #define sylvan_storm_rational_number_less_or_equal(a, b) mtbdd_apply(a, b, TASK(sylvan_storm_rational_number_op_less_or_equal)) -#define sylvan_storm_rational_number_mod(a, b) mtbdd_apply(a, b, TASK(sylvan_storm_rational_number_op_less_or_equal)) +#define sylvan_storm_rational_number_greater_or_equal(a, b) mtbdd_apply(a, b, TASK(sylvan_storm_rational_number_op_greater_or_equal)) +#define sylvan_storm_rational_number_mod(a, b) mtbdd_apply(a, b, TASK(sylvan_storm_rational_number_op_mod)) #define sylvan_storm_rational_number_min(a, b) mtbdd_apply(a, b, TASK(sylvan_storm_rational_number_op_min)) #define sylvan_storm_rational_number_max(a, b) mtbdd_apply(a, b, TASK(sylvan_storm_rational_number_op_max)) #define sylvan_storm_rational_number_neg(a) mtbdd_uapply(a, TASK(sylvan_storm_rational_number_op_neg), 0) @@ -87,6 +91,12 @@ TASK_DECL_3(MTBDD, sylvan_storm_rational_number_equal_norm_d, MTBDD, MTBDD, stor TASK_DECL_3(MTBDD, sylvan_storm_rational_number_equal_norm_rel_d, MTBDD, MTBDD, storm_rational_number_ptr); #define sylvan_storm_rational_number_equal_norm_rel_d(a, b, epsilon) CALL(sylvan_storm_rational_number_equal_norm_rel_d, a, b, epsilon) +TASK_DECL_3(BDD, sylvan_storm_rational_number_min_abstract_representative, MTBDD, MTBDD, uint32_t); +#define sylvan_storm_rational_number_min_abstract_representative(a, vars) (CALL(sylvan_storm_rational_number_min_abstract_representative, a, vars, 0)) + +TASK_DECL_3(BDD, sylvan_storm_rational_number_max_abstract_representative, MTBDD, MTBDD, uint32_t); +#define sylvan_storm_rational_number_max_abstract_representative(a, vars) (CALL(sylvan_storm_rational_number_max_abstract_representative, a, vars, 0)) + #ifdef __cplusplus } #endif /* __cplusplus */ diff --git a/src/storm/modelchecker/prctl/SymbolicMdpPrctlModelChecker.cpp b/src/storm/modelchecker/prctl/SymbolicMdpPrctlModelChecker.cpp index cf88dd848..1eebe4842 100644 --- a/src/storm/modelchecker/prctl/SymbolicMdpPrctlModelChecker.cpp +++ b/src/storm/modelchecker/prctl/SymbolicMdpPrctlModelChecker.cpp @@ -22,12 +22,12 @@ namespace storm { namespace modelchecker { template - SymbolicMdpPrctlModelChecker::SymbolicMdpPrctlModelChecker(ModelType const& model, std::unique_ptr>&& linearEquationSolverFactory) : SymbolicPropositionalModelChecker(model), linearEquationSolverFactory(std::move(linearEquationSolverFactory)) { + SymbolicMdpPrctlModelChecker::SymbolicMdpPrctlModelChecker(ModelType const& model, std::unique_ptr>&& linearEquationSolverFactory) : SymbolicPropositionalModelChecker(model), linearEquationSolverFactory(std::move(linearEquationSolverFactory)) { // Intentionally left empty. } template - SymbolicMdpPrctlModelChecker::SymbolicMdpPrctlModelChecker(ModelType const& model) : SymbolicPropositionalModelChecker(model), linearEquationSolverFactory(new storm::utility::solver::SymbolicMinMaxLinearEquationSolverFactory()) { + SymbolicMdpPrctlModelChecker::SymbolicMdpPrctlModelChecker(ModelType const& model) : SymbolicPropositionalModelChecker(model), linearEquationSolverFactory(new storm::solver::SymbolicGeneralMinMaxLinearEquationSolverFactory()) { // Intentionally left empty. } diff --git a/src/storm/modelchecker/prctl/SymbolicMdpPrctlModelChecker.h b/src/storm/modelchecker/prctl/SymbolicMdpPrctlModelChecker.h index 832a6b959..c45544d98 100644 --- a/src/storm/modelchecker/prctl/SymbolicMdpPrctlModelChecker.h +++ b/src/storm/modelchecker/prctl/SymbolicMdpPrctlModelChecker.h @@ -5,7 +5,7 @@ #include "storm/models/symbolic/Mdp.h" -#include "storm/utility/solver.h" +#include "storm/solver/SymbolicMinMaxLinearEquationSolver.h" namespace storm { namespace modelchecker { @@ -16,7 +16,7 @@ namespace storm { static const storm::dd::DdType DdType = ModelType::DdType; explicit SymbolicMdpPrctlModelChecker(ModelType const& model); - explicit SymbolicMdpPrctlModelChecker(ModelType const& model, std::unique_ptr>&& linearEquationSolverFactory); + explicit SymbolicMdpPrctlModelChecker(ModelType const& model, std::unique_ptr>&& linearEquationSolverFactory); // The implemented methods of the AbstractModelChecker interface. virtual bool canHandle(CheckTask const& checkTask) const override; @@ -31,7 +31,7 @@ namespace storm { private: // An object that is used for retrieving linear equation solvers. - std::unique_ptr> linearEquationSolverFactory; + std::unique_ptr> linearEquationSolverFactory; }; } // namespace modelchecker diff --git a/src/storm/modelchecker/prctl/helper/SymbolicMdpPrctlHelper.cpp b/src/storm/modelchecker/prctl/helper/SymbolicMdpPrctlHelper.cpp index c2d86f663..74a52cf90 100644 --- a/src/storm/modelchecker/prctl/helper/SymbolicMdpPrctlHelper.cpp +++ b/src/storm/modelchecker/prctl/helper/SymbolicMdpPrctlHelper.cpp @@ -23,7 +23,7 @@ namespace storm { namespace helper { template - std::unique_ptr SymbolicMdpPrctlHelper::computeUntilProbabilities(OptimizationDirection dir, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& phiStates, storm::dd::Bdd const& psiStates, bool qualitative, storm::utility::solver::SymbolicMinMaxLinearEquationSolverFactory const& linearEquationSolverFactory) { + std::unique_ptr SymbolicMdpPrctlHelper::computeUntilProbabilities(OptimizationDirection dir, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& phiStates, storm::dd::Bdd const& psiStates, bool qualitative, storm::solver::SymbolicGeneralMinMaxLinearEquationSolverFactory const& linearEquationSolverFactory) { // We need to identify the states which have to be taken out of the matrix, i.e. all states that have // probability 0 and 1 of satisfying the until-formula. std::pair, storm::dd::Bdd> statesWithProbability01; @@ -76,7 +76,7 @@ namespace storm { } template - std::unique_ptr SymbolicMdpPrctlHelper::computeGloballyProbabilities(OptimizationDirection dir, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& psiStates, bool qualitative, storm::utility::solver::SymbolicMinMaxLinearEquationSolverFactory const& linearEquationSolverFactory) { + std::unique_ptr SymbolicMdpPrctlHelper::computeGloballyProbabilities(OptimizationDirection dir, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& psiStates, bool qualitative, storm::solver::SymbolicGeneralMinMaxLinearEquationSolverFactory const& linearEquationSolverFactory) { std::unique_ptr result = computeUntilProbabilities(dir == OptimizationDirection::Minimize ? OptimizationDirection::Maximize : OptimizationDirection::Minimize, model, transitionMatrix, model.getReachableStates(), !psiStates && model.getReachableStates(), qualitative, linearEquationSolverFactory); result->asQuantitativeCheckResult().oneMinus(); return result; @@ -94,7 +94,7 @@ namespace storm { } template - std::unique_ptr SymbolicMdpPrctlHelper::computeBoundedUntilProbabilities(OptimizationDirection dir, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& phiStates, storm::dd::Bdd const& psiStates, uint_fast64_t stepBound, storm::utility::solver::SymbolicMinMaxLinearEquationSolverFactory const& linearEquationSolverFactory) { + std::unique_ptr SymbolicMdpPrctlHelper::computeBoundedUntilProbabilities(OptimizationDirection dir, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& phiStates, storm::dd::Bdd const& psiStates, uint_fast64_t stepBound, storm::solver::SymbolicGeneralMinMaxLinearEquationSolverFactory const& linearEquationSolverFactory) { // We need to identify the states which have to be taken out of the matrix, i.e. all states that have // probability 0 or 1 of satisfying the until-formula. storm::dd::Bdd statesWithProbabilityGreater0; @@ -132,7 +132,7 @@ namespace storm { } template - std::unique_ptr SymbolicMdpPrctlHelper::computeInstantaneousRewards(OptimizationDirection dir, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, RewardModelType const& rewardModel, uint_fast64_t stepBound, storm::utility::solver::SymbolicMinMaxLinearEquationSolverFactory const& linearEquationSolverFactory) { + std::unique_ptr SymbolicMdpPrctlHelper::computeInstantaneousRewards(OptimizationDirection dir, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, RewardModelType const& rewardModel, uint_fast64_t stepBound, storm::solver::SymbolicGeneralMinMaxLinearEquationSolverFactory const& linearEquationSolverFactory) { // Only compute the result if the model has at least one reward this->getModel(). STORM_LOG_THROW(rewardModel.hasStateRewards(), storm::exceptions::InvalidPropertyException, "Missing reward model for formula. Skipping formula."); @@ -144,7 +144,7 @@ namespace storm { } template - std::unique_ptr SymbolicMdpPrctlHelper::computeCumulativeRewards(OptimizationDirection dir, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, RewardModelType const& rewardModel, uint_fast64_t stepBound, storm::utility::solver::SymbolicMinMaxLinearEquationSolverFactory const& linearEquationSolverFactory) { + std::unique_ptr SymbolicMdpPrctlHelper::computeCumulativeRewards(OptimizationDirection dir, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, RewardModelType const& rewardModel, uint_fast64_t stepBound, storm::solver::SymbolicGeneralMinMaxLinearEquationSolverFactory const& linearEquationSolverFactory) { // Only compute the result if the model has at least one reward this->getModel(). STORM_LOG_THROW(!rewardModel.empty(), storm::exceptions::InvalidPropertyException, "Missing reward model for formula. Skipping formula."); @@ -159,7 +159,7 @@ namespace storm { } template - std::unique_ptr SymbolicMdpPrctlHelper::computeReachabilityRewards(OptimizationDirection dir, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, RewardModelType const& rewardModel, storm::dd::Bdd const& targetStates, bool qualitative, storm::utility::solver::SymbolicMinMaxLinearEquationSolverFactory const& linearEquationSolverFactory) { + std::unique_ptr SymbolicMdpPrctlHelper::computeReachabilityRewards(OptimizationDirection dir, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, RewardModelType const& rewardModel, storm::dd::Bdd const& targetStates, bool qualitative, storm::solver::SymbolicGeneralMinMaxLinearEquationSolverFactory const& linearEquationSolverFactory) { // Only compute the result if there is at least one reward model. STORM_LOG_THROW(!rewardModel.empty(), storm::exceptions::InvalidPropertyException, "Missing reward model for formula. Skipping formula."); diff --git a/src/storm/modelchecker/prctl/helper/SymbolicMdpPrctlHelper.h b/src/storm/modelchecker/prctl/helper/SymbolicMdpPrctlHelper.h index 1d422c416..0dd68f188 100644 --- a/src/storm/modelchecker/prctl/helper/SymbolicMdpPrctlHelper.h +++ b/src/storm/modelchecker/prctl/helper/SymbolicMdpPrctlHelper.h @@ -6,7 +6,7 @@ #include "storm/storage/dd/Add.h" #include "storm/storage/dd/Bdd.h" -#include "storm/utility/solver.h" +#include "storm/solver/SymbolicMinMaxLinearEquationSolver.h" #include "storm/solver/SolveGoal.h" namespace storm { @@ -21,19 +21,19 @@ namespace storm { public: typedef typename storm::models::symbolic::NondeterministicModel::RewardModelType RewardModelType; - static std::unique_ptr computeBoundedUntilProbabilities(OptimizationDirection dir, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& phiStates, storm::dd::Bdd const& psiStates, uint_fast64_t stepBound, storm::utility::solver::SymbolicMinMaxLinearEquationSolverFactory const& linearEquationSolverFactory); + static std::unique_ptr computeBoundedUntilProbabilities(OptimizationDirection dir, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& phiStates, storm::dd::Bdd const& psiStates, uint_fast64_t stepBound, storm::solver::SymbolicGeneralMinMaxLinearEquationSolverFactory const& linearEquationSolverFactory); static std::unique_ptr computeNextProbabilities(OptimizationDirection dir, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& nextStates); - static std::unique_ptr computeUntilProbabilities(OptimizationDirection dir, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& phiStates, storm::dd::Bdd const& psiStates, bool qualitative, storm::utility::solver::SymbolicMinMaxLinearEquationSolverFactory const& linearEquationSolverFactory); + static std::unique_ptr computeUntilProbabilities(OptimizationDirection dir, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& phiStates, storm::dd::Bdd const& psiStates, bool qualitative, storm::solver::SymbolicGeneralMinMaxLinearEquationSolverFactory const& linearEquationSolverFactory); - static std::unique_ptr computeGloballyProbabilities(OptimizationDirection dir, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& psiStates, bool qualitative, storm::utility::solver::SymbolicMinMaxLinearEquationSolverFactory const& linearEquationSolverFactory); + static std::unique_ptr computeGloballyProbabilities(OptimizationDirection dir, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, storm::dd::Bdd const& psiStates, bool qualitative, storm::solver::SymbolicGeneralMinMaxLinearEquationSolverFactory const& linearEquationSolverFactory); - static std::unique_ptr computeCumulativeRewards(OptimizationDirection dir, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, RewardModelType const& rewardModel, uint_fast64_t stepBound, storm::utility::solver::SymbolicMinMaxLinearEquationSolverFactory const& linearEquationSolverFactory); + static std::unique_ptr computeCumulativeRewards(OptimizationDirection dir, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, RewardModelType const& rewardModel, uint_fast64_t stepBound, storm::solver::SymbolicGeneralMinMaxLinearEquationSolverFactory const& linearEquationSolverFactory); - static std::unique_ptr computeInstantaneousRewards(OptimizationDirection dir, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, RewardModelType const& rewardModel, uint_fast64_t stepBound, storm::utility::solver::SymbolicMinMaxLinearEquationSolverFactory const& linearEquationSolverFactory); + static std::unique_ptr computeInstantaneousRewards(OptimizationDirection dir, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, RewardModelType const& rewardModel, uint_fast64_t stepBound, storm::solver::SymbolicGeneralMinMaxLinearEquationSolverFactory const& linearEquationSolverFactory); - static std::unique_ptr computeReachabilityRewards(OptimizationDirection dir, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, RewardModelType const& rewardModel, storm::dd::Bdd const& targetStates, bool qualitative, storm::utility::solver::SymbolicMinMaxLinearEquationSolverFactory const& linearEquationSolverFactory); + static std::unique_ptr computeReachabilityRewards(OptimizationDirection dir, storm::models::symbolic::NondeterministicModel const& model, storm::dd::Add const& transitionMatrix, RewardModelType const& rewardModel, storm::dd::Bdd const& targetStates, bool qualitative, storm::solver::SymbolicGeneralMinMaxLinearEquationSolverFactory const& linearEquationSolverFactory); }; } diff --git a/src/storm/modelchecker/results/SymbolicQuantitativeCheckResult.cpp b/src/storm/modelchecker/results/SymbolicQuantitativeCheckResult.cpp index 87e8ff5d3..0b47a0533 100644 --- a/src/storm/modelchecker/results/SymbolicQuantitativeCheckResult.cpp +++ b/src/storm/modelchecker/results/SymbolicQuantitativeCheckResult.cpp @@ -57,11 +57,62 @@ namespace storm { return values; } + template + void print(std::ostream& out, ValueType const& value) { + if (value == storm::utility::infinity()) { + out << "inf"; + } else { + out << value; + if (std::is_same::value) { + out << " (approx. " << storm::utility::convertNumber(value) << ")"; + } + } + } + + template + void printRange(std::ostream& out, ValueType const& min, ValueType const& max) { + out << "["; + if (min == storm::utility::infinity()) { + out << "inf"; + } else { + out << min; + } + out << ", "; + if (max == storm::utility::infinity()) { + out << "inf"; + } else { + out << max; + } + out << "]"; + if (std::is_same::value) { + out << " (approx. ["; + if (min == storm::utility::infinity()) { + out << "inf"; + } else { + out << storm::utility::convertNumber(min); + } + out << ", "; + if (max == storm::utility::infinity()) { + out << "inf"; + } else { + out << storm::utility::convertNumber(max); + } + out << "])"; + } + out << " (range)"; + } + template std::ostream& SymbolicQuantitativeCheckResult::writeToStream(std::ostream& out) const { - if (states.getNonZeroCount() == 1) { - out << this->values.sumAbstract(this->values.getContainedMetaVariables()).getValue(); - } else if (states.getNonZeroCount() < 10 || std::is_same::value) { + uint64_t totalNumberOfStates = states.getNonZeroCount(); + bool minMaxSupported = std::is_same::value || std::is_same::value; + bool printAsRange = false; + + if (totalNumberOfStates == 1) { + print(out, this->values.sumAbstract(this->values.getContainedMetaVariables()).getValue()); + } else if (states.getNonZeroCount() >= 10 || !minMaxSupported) { + printAsRange = true; + } else { out << "{"; if (this->values.isZero()) { out << "0"; @@ -73,18 +124,17 @@ namespace storm { } else { first = false; } - out << valuationValuePair.second; + print(out, valuationValuePair.second); } if (states.getNonZeroCount() != this->values.getNonZeroCount()) { out << ", 0"; } } out << "}"; - } else { - ValueType min = this->getMin(); - ValueType max = this->getMax(); - - out << "[" << min << ", " << max << "] (range)"; + } + + if (printAsRange) { + printRange(out, this->getMin(), this->getMax()); } return out; } diff --git a/src/storm/solver/SymbolicEliminationLinearEquationSolver.cpp b/src/storm/solver/SymbolicEliminationLinearEquationSolver.cpp index 0b515451a..84a343cda 100644 --- a/src/storm/solver/SymbolicEliminationLinearEquationSolver.cpp +++ b/src/storm/solver/SymbolicEliminationLinearEquationSolver.cpp @@ -12,27 +12,14 @@ namespace storm { template SymbolicEliminationLinearEquationSolver::SymbolicEliminationLinearEquationSolver(storm::dd::Add const& A, storm::dd::Bdd const& allRows, std::set const& rowMetaVariables, std::set const& columnMetaVariables, std::vector> const& rowColumnMetaVariablePairs) : SymbolicLinearEquationSolver(A, allRows, rowMetaVariables, columnMetaVariables, rowColumnMetaVariablePairs) { - // Intentionally left empty. - } - - template - storm::dd::Add SymbolicEliminationLinearEquationSolver::solveEquations(storm::dd::Add const& x, storm::dd::Add const& b) const { - storm::dd::DdManager& ddManager = x.getDdManager(); + storm::dd::DdManager& ddManager = A.getDdManager(); - // We start by creating triple-layered meta variables for all original meta variables. We will use them later in the elimination process. - std::vector> oldToNewMapping; - std::set newRowVariables; - std::set newColumnVariables; - std::set helperVariables; - std::vector> newRowColumnMetaVariablePairs; - std::vector> columnHelperMetaVariablePairs; - std::vector> rowRowMetaVariablePairs; - std::vector> columnColumnMetaVariablePairs; + // Create triple-layered meta variables for all original meta variables. We will use them later in the elimination process. for (auto const& metaVariablePair : this->rowColumnMetaVariablePairs) { auto rowVariable = metaVariablePair.first; storm::dd::DdMetaVariable const& metaVariable = ddManager.getMetaVariable(rowVariable); - + std::vector newMetaVariables; if (metaVariable.getType() == storm::dd::MetaVariableType::Bool) { newMetaVariables = ddManager.addMetaVariable("tmp_" + metaVariable.getName(), 3); @@ -43,7 +30,7 @@ namespace storm { newRowVariables.insert(newMetaVariables[0]); newColumnVariables.insert(newMetaVariables[1]); helperVariables.insert(newMetaVariables[2]); - + newRowColumnMetaVariablePairs.emplace_back(newMetaVariables[0], newMetaVariables[1]); columnHelperMetaVariablePairs.emplace_back(newMetaVariables[1], newMetaVariables[2]); @@ -53,21 +40,25 @@ namespace storm { oldToNewMapping.emplace_back(std::move(newMetaVariables)); } - std::vector> oldNewMetaVariablePairs = rowRowMetaVariablePairs; + oldNewMetaVariablePairs = rowRowMetaVariablePairs; for (auto const& entry : columnColumnMetaVariablePairs) { oldNewMetaVariablePairs.emplace_back(entry.first, entry.second); } - std::vector> shiftMetaVariablePairs = newRowColumnMetaVariablePairs; + shiftMetaVariablePairs = newRowColumnMetaVariablePairs; for (auto const& entry : columnHelperMetaVariablePairs) { shiftMetaVariablePairs.emplace_back(entry.first, entry.second); } - + } + + template + storm::dd::Add SymbolicEliminationLinearEquationSolver::solveEquations(storm::dd::Add const& x, storm::dd::Add const& b) const { + storm::dd::DdManager& ddManager = x.getDdManager(); + // Build diagonal BDD over new meta variables. storm::dd::Bdd diagonal = storm::utility::dd::getRowColumnDiagonal(x.getDdManager(), this->rowColumnMetaVariablePairs); diagonal &= this->allRows; - diagonal = diagonal.swapVariables(oldNewMetaVariablePairs); - + diagonal = diagonal.swapVariables(this->oldNewMetaVariablePairs); storm::dd::Add rowsAdd = this->allRows.swapVariables(rowRowMetaVariablePairs).template toAdd(); storm::dd::Add diagonalAdd = diagonal.template toAdd(); diff --git a/src/storm/solver/SymbolicEliminationLinearEquationSolver.h b/src/storm/solver/SymbolicEliminationLinearEquationSolver.h index 0d96a3cd6..be8164d43 100644 --- a/src/storm/solver/SymbolicEliminationLinearEquationSolver.h +++ b/src/storm/solver/SymbolicEliminationLinearEquationSolver.h @@ -17,6 +17,18 @@ namespace storm { SymbolicEliminationLinearEquationSolver(storm::dd::Add const& A, storm::dd::Bdd const& allRows, std::set const& rowMetaVariables, std::set const& columnMetaVariables, std::vector> const& rowColumnMetaVariablePairs); virtual storm::dd::Add solveEquations(storm::dd::Add const& x, storm::dd::Add const& b) const override; + + private: + std::vector> oldToNewMapping; + std::set newRowVariables; + std::set newColumnVariables; + std::set helperVariables; + std::vector> newRowColumnMetaVariablePairs; + std::vector> columnHelperMetaVariablePairs; + std::vector> rowRowMetaVariablePairs; + std::vector> columnColumnMetaVariablePairs; + std::vector> oldNewMetaVariablePairs = rowRowMetaVariablePairs; + std::vector> shiftMetaVariablePairs = newRowColumnMetaVariablePairs; }; template diff --git a/src/storm/solver/SymbolicLinearEquationSolver.cpp b/src/storm/solver/SymbolicLinearEquationSolver.cpp index e5afcb08a..13c220418 100644 --- a/src/storm/solver/SymbolicLinearEquationSolver.cpp +++ b/src/storm/solver/SymbolicLinearEquationSolver.cpp @@ -40,6 +40,10 @@ namespace storm { return xCopy; } + template + void SymbolicLinearEquationSolver::setMatrix(storm::dd::Add const& newA) { + this->A = newA; + } template std::unique_ptr> GeneralSymbolicLinearEquationSolverFactory::create(storm::dd::Add const& A, storm::dd::Bdd const& allRows, std::set const& rowMetaVariables, std::set const& columnMetaVariables, std::vector> const& rowColumnMetaVariablePairs) const { diff --git a/src/storm/solver/SymbolicLinearEquationSolver.h b/src/storm/solver/SymbolicLinearEquationSolver.h index dddc80596..a83054d5c 100644 --- a/src/storm/solver/SymbolicLinearEquationSolver.h +++ b/src/storm/solver/SymbolicLinearEquationSolver.h @@ -71,6 +71,8 @@ namespace storm { */ virtual storm::dd::Add multiply(storm::dd::Add const& x, storm::dd::Add const* b = nullptr, uint_fast64_t n = 1) const; + void setMatrix(storm::dd::Add const& newA); + protected: // The matrix defining the coefficients of the linear equation system. storm::dd::Add A; diff --git a/src/storm/solver/SymbolicMinMaxLinearEquationSolver.cpp b/src/storm/solver/SymbolicMinMaxLinearEquationSolver.cpp index 9eccde6d4..9769da00f 100644 --- a/src/storm/solver/SymbolicMinMaxLinearEquationSolver.cpp +++ b/src/storm/solver/SymbolicMinMaxLinearEquationSolver.cpp @@ -1,5 +1,7 @@ #include "storm/solver/SymbolicMinMaxLinearEquationSolver.h" + + #include "storm/storage/dd/DdManager.h" #include "storm/storage/dd/Add.h" @@ -8,35 +10,98 @@ #include "storm/utility/constants.h" #include "storm/settings/SettingsManager.h" -#include "storm/settings/modules/NativeEquationSolverSettings.h" +#include "storm/settings/modules/MinMaxEquationSolverSettings.h" + +#include "storm/utility/dd.h" +#include "storm/utility/macros.h" +#include "storm/exceptions/InvalidSettingsException.h" namespace storm { namespace solver { + template + SymbolicMinMaxLinearEquationSolverSettings::SymbolicMinMaxLinearEquationSolverSettings() { + // Get the settings object to customize linear solving. + storm::settings::modules::MinMaxEquationSolverSettings const& settings = storm::settings::getModule(); + + maximalNumberOfIterations = settings.getMaximalIterationCount(); + precision = storm::utility::convertNumber(settings.getPrecision()); + relative = settings.getConvergenceCriterion() == storm::settings::modules::MinMaxEquationSolverSettings::ConvergenceCriterion::Relative; + + auto method = settings.getMinMaxEquationSolvingMethod(); + switch (method) { + case MinMaxMethod::ValueIteration: this->solutionMethod = SolutionMethod::ValueIteration; break; + case MinMaxMethod::PolicyIteration: this->solutionMethod = SolutionMethod::PolicyIteration; break; + default: + STORM_LOG_THROW(false, storm::exceptions::InvalidSettingsException, "Unsupported technique."); + } + } + + template + void SymbolicMinMaxLinearEquationSolverSettings::setSolutionMethod(SolutionMethod const& solutionMethod) { + this->solutionMethod = solutionMethod; + } + + template + void SymbolicMinMaxLinearEquationSolverSettings::setMaximalNumberOfIterations(uint64_t maximalNumberOfIterations) { + this->maximalNumberOfIterations = maximalNumberOfIterations; + } + + template + void SymbolicMinMaxLinearEquationSolverSettings::setRelativeTerminationCriterion(bool value) { + this->relative = value; + } + + template + void SymbolicMinMaxLinearEquationSolverSettings::setPrecision(ValueType precision) { + this->precision = precision; + } + + template + typename SymbolicMinMaxLinearEquationSolverSettings::SolutionMethod const& SymbolicMinMaxLinearEquationSolverSettings::getSolutionMethod() const { + return solutionMethod; + } + + template + uint64_t SymbolicMinMaxLinearEquationSolverSettings::getMaximalNumberOfIterations() const { + return maximalNumberOfIterations; + } + + template + ValueType SymbolicMinMaxLinearEquationSolverSettings::getPrecision() const { + return precision; + } + + template + bool SymbolicMinMaxLinearEquationSolverSettings::getRelativeTerminationCriterion() const { + return relative; + } + template - SymbolicMinMaxLinearEquationSolver::SymbolicMinMaxLinearEquationSolver(storm::dd::Add const& A, storm::dd::Bdd const& allRows, storm::dd::Bdd const& illegalMask, std::set const& rowMetaVariables, std::set const& columnMetaVariables, std::set const& choiceVariables, std::vector> const& rowColumnMetaVariablePairs, double precision, uint_fast64_t maximalNumberOfIterations, bool relative) : A(A), allRows(allRows), illegalMaskAdd(illegalMask.ite(A.getDdManager().getConstant(storm::utility::infinity()), A.getDdManager().template getAddZero())), rowMetaVariables(rowMetaVariables), columnMetaVariables(columnMetaVariables), choiceVariables(choiceVariables), rowColumnMetaVariablePairs(rowColumnMetaVariablePairs), precision(precision), maximalNumberOfIterations(maximalNumberOfIterations), relative(relative) { + SymbolicMinMaxLinearEquationSolver::SymbolicMinMaxLinearEquationSolver(storm::dd::Add const& A, storm::dd::Bdd const& allRows, storm::dd::Bdd const& illegalMask, std::set const& rowMetaVariables, std::set const& columnMetaVariables, std::set const& choiceVariables, std::vector> const& rowColumnMetaVariablePairs, std::unique_ptr>&& linearEquationSolverFactory, SymbolicMinMaxLinearEquationSolverSettings const& settings) : A(A), allRows(allRows), illegalMaskAdd(illegalMask.ite(A.getDdManager().getConstant(storm::utility::infinity()), A.getDdManager().template getAddZero())), rowMetaVariables(rowMetaVariables), columnMetaVariables(columnMetaVariables), choiceVariables(choiceVariables), rowColumnMetaVariablePairs(rowColumnMetaVariablePairs), linearEquationSolverFactory(std::move(linearEquationSolverFactory)), settings(settings) { // Intentionally left empty. } template - SymbolicMinMaxLinearEquationSolver::SymbolicMinMaxLinearEquationSolver(storm::dd::Add const& A, storm::dd::Bdd const& allRows, storm::dd::Bdd const& illegalMask, std::set const& rowMetaVariables, std::set const& columnMetaVariables, std::set const& choiceVariables, std::vector> const& rowColumnMetaVariablePairs) : A(A), allRows(allRows), illegalMaskAdd(illegalMask.ite(A.getDdManager().getConstant(storm::utility::infinity()), A.getDdManager().template getAddZero())), rowMetaVariables(rowMetaVariables), columnMetaVariables(columnMetaVariables), choiceVariables(choiceVariables), rowColumnMetaVariablePairs(rowColumnMetaVariablePairs) { - // Get the settings object to customize solving. - storm::settings::modules::NativeEquationSolverSettings const& settings = storm::settings::getModule(); - - // Get appropriate settings. - maximalNumberOfIterations = settings.getMaximalIterationCount(); - precision = settings.getPrecision(); - relative = settings.getConvergenceCriterion() == storm::settings::modules::NativeEquationSolverSettings::ConvergenceCriterion::Relative; + storm::dd::Add SymbolicMinMaxLinearEquationSolver::solveEquations(bool minimize, storm::dd::Add const& x, storm::dd::Add const& b) const { + switch (this->getSettings().getSolutionMethod()) { + case SymbolicMinMaxLinearEquationSolverSettings::SolutionMethod::ValueIteration: + return solveEquationsValueIteration(minimize, x, b); + break; + case SymbolicMinMaxLinearEquationSolverSettings::SolutionMethod::PolicyIteration: + return solveEquationsPolicyIteration(minimize, x, b); + break; + } } template - storm::dd::Add SymbolicMinMaxLinearEquationSolver::solveEquations(bool minimize, storm::dd::Add const& x, storm::dd::Add const& b) const { + storm::dd::Add SymbolicMinMaxLinearEquationSolver::solveEquationsValueIteration(bool minimize, storm::dd::Add const& x, storm::dd::Add const& b) const { // Set up the environment. storm::dd::Add xCopy = x; uint_fast64_t iterations = 0; bool converged = false; - while (!converged && iterations < maximalNumberOfIterations) { + while (!converged && iterations < this->settings.getMaximalNumberOfIterations()) { // Compute tmp = A * x + b storm::dd::Add xCopyAsColumn = xCopy.swapVariables(this->rowColumnMetaVariablePairs); storm::dd::Add tmp = this->A.multiplyMatrix(xCopyAsColumn, this->columnMetaVariables); @@ -50,22 +115,79 @@ namespace storm { } // Now check if the process already converged within our precision. - converged = xCopy.equalModuloPrecision(tmp, precision, relative); + converged = xCopy.equalModuloPrecision(tmp, this->settings.getPrecision(), this->settings.getRelativeTerminationCriterion()); xCopy = tmp; ++iterations; } - + if (converged) { - STORM_LOG_TRACE("Iterative solver converged in " << iterations << " iterations."); + STORM_LOG_TRACE("Iterative solver (value iteration) converged in " << iterations << " iterations."); } else { - STORM_LOG_WARN("Iterative solver did not converge in " << iterations << " iterations."); + STORM_LOG_WARN("Iterative solver (value iteration) did not converge in " << iterations << " iterations."); } - + return xCopy; } + template + storm::dd::Add SymbolicMinMaxLinearEquationSolver::solveEquationsPolicyIteration(bool minimize, storm::dd::Add const& x, storm::dd::Add const& b) const { + // Set up the environment. + storm::dd::Add currentSolution = x; + storm::dd::Add diagonal = (storm::utility::dd::getRowColumnDiagonal(x.getDdManager(), this->rowColumnMetaVariablePairs) && this->allRows).template toAdd(); + uint_fast64_t iterations = 0; + bool converged = false; + + // Pick arbitrary initial scheduler. + storm::dd::Bdd scheduler = this->A.sumAbstract(this->columnMetaVariables).maxAbstractRepresentative(this->choiceVariables); + + // And apply it to the matrix and vector. + storm::dd::Add schedulerA = diagonal - scheduler.ite(this->A, scheduler.getDdManager().template getAddZero()).sumAbstract(this->choiceVariables); + storm::dd::Add schedulerB = scheduler.ite(b, scheduler.getDdManager().template getAddZero()).sumAbstract(this->choiceVariables); + + // Initialize linear equation solver. + std::unique_ptr> linearEquationSolver = linearEquationSolverFactory->create(schedulerA, this->allRows, this->rowMetaVariables, this->columnMetaVariables, this->rowColumnMetaVariablePairs); + + // Iteratively solve and improve the scheduler. + while (!converged && iterations < this->settings.getMaximalNumberOfIterations()) { + // Solve for the value of the scheduler. + storm::dd::Add schedulerX = linearEquationSolver->solveEquations(currentSolution, schedulerB); + + // Policy improvement step. + storm::dd::Add tmp = this->A.multiplyMatrix(schedulerX.swapVariables(this->rowColumnMetaVariablePairs), this->columnMetaVariables); + storm::dd::Bdd nextScheduler; + if (minimize) { + tmp += illegalMaskAdd; + nextScheduler = tmp.minAbstractRepresentative(this->choiceVariables); + } else { + nextScheduler = tmp.maxAbstractRepresentative(this->choiceVariables); + } + + // Check for convergence. + converged = nextScheduler == scheduler; + + // Set up next iteration. + if (!converged) { + scheduler = nextScheduler; + schedulerA = diagonal - scheduler.ite(this->A, scheduler.getDdManager().template getAddZero()).sumAbstract(this->choiceVariables); + linearEquationSolver->setMatrix(schedulerA); + schedulerB = scheduler.ite(b, scheduler.getDdManager().template getAddZero()).sumAbstract(this->choiceVariables); + } + + currentSolution = schedulerX; + ++iterations; + } + + if (converged) { + STORM_LOG_TRACE("Iterative solver (policy iteration) converged in " << iterations << " iterations."); + } else { + STORM_LOG_WARN("Iterative solver (policy iteration) did not converge in " << iterations << " iterations."); + } + + return currentSolution; + } + template storm::dd::Add SymbolicMinMaxLinearEquationSolver::multiply(bool minimize, storm::dd::Add const& x, storm::dd::Add const* b, uint_fast64_t n) const { storm::dd::Add xCopy = x; @@ -90,16 +212,37 @@ namespace storm { return xCopy; } + + template + SymbolicMinMaxLinearEquationSolverSettings const& SymbolicMinMaxLinearEquationSolver::getSettings() const { + return settings; + } template - ValueType const& SymbolicMinMaxLinearEquationSolver::getPrecision() const { - return precision; + std::unique_ptr> SymbolicGeneralMinMaxLinearEquationSolverFactory::create(storm::dd::Add const& A, storm::dd::Bdd const& allRows, storm::dd::Bdd const& illegalMask, std::set const& rowMetaVariables, std::set const& columnMetaVariables, std::set const& choiceVariables, std::vector> const& rowColumnMetaVariablePairs) const { + return std::make_unique>(A, allRows, illegalMask, rowMetaVariables, columnMetaVariables, choiceVariables, rowColumnMetaVariablePairs, std::make_unique>(), settings); } + + template + SymbolicMinMaxLinearEquationSolverSettings& SymbolicGeneralMinMaxLinearEquationSolverFactory::getSettings() { + return settings; + } + + template + SymbolicMinMaxLinearEquationSolverSettings const& SymbolicGeneralMinMaxLinearEquationSolverFactory::getSettings() const { + return settings; + } + + template class SymbolicMinMaxLinearEquationSolverSettings; + template class SymbolicMinMaxLinearEquationSolverSettings; template class SymbolicMinMaxLinearEquationSolver; template class SymbolicMinMaxLinearEquationSolver; - template class SymbolicMinMaxLinearEquationSolver; + template class SymbolicGeneralMinMaxLinearEquationSolverFactory; + template class SymbolicGeneralMinMaxLinearEquationSolverFactory; + template class SymbolicGeneralMinMaxLinearEquationSolverFactory; + } } diff --git a/src/storm/solver/SymbolicMinMaxLinearEquationSolver.h b/src/storm/solver/SymbolicMinMaxLinearEquationSolver.h index c343b281b..650934211 100644 --- a/src/storm/solver/SymbolicMinMaxLinearEquationSolver.h +++ b/src/storm/solver/SymbolicMinMaxLinearEquationSolver.h @@ -6,6 +6,8 @@ #include #include +#include "storm/solver/SymbolicLinearEquationSolver.h" + #include "storm/storage/expressions/Variable.h" #include "storm/storage/dd/DdType.h" @@ -19,6 +21,32 @@ namespace storm { } namespace solver { + template + class SymbolicMinMaxLinearEquationSolverSettings { + public: + SymbolicMinMaxLinearEquationSolverSettings(); + + enum class SolutionMethod { + ValueIteration, PolicyIteration + }; + + void setSolutionMethod(SolutionMethod const& solutionMethod); + void setMaximalNumberOfIterations(uint64_t maximalNumberOfIterations); + void setRelativeTerminationCriterion(bool value); + void setPrecision(ValueType precision); + + SolutionMethod const& getSolutionMethod() const; + uint64_t getMaximalNumberOfIterations() const; + ValueType getPrecision() const; + bool getRelativeTerminationCriterion() const; + + private: + SolutionMethod solutionMethod; + uint64_t maximalNumberOfIterations; + ValueType precision; + bool relative; + }; + /*! * An interface that represents an abstract symbolic linear equation solver. In addition to solving a system of * linear equations, the functionality to repeatedly multiply a matrix with a given vector is provided. @@ -38,27 +66,10 @@ namespace storm { * @param choiceVariables The variables encoding the choices of each row group. * @param rowColumnMetaVariablePairs The pairs of row meta variables and the corresponding column meta * variables. + * @param settings The settings to use. */ - SymbolicMinMaxLinearEquationSolver(storm::dd::Add const& A, storm::dd::Bdd const& allRows, storm::dd::Bdd const& illegalMask, std::set const& rowMetaVariables, std::set const& columnMetaVariables, std::set const& choiceVariables, std::vector> const& rowColumnMetaVariablePairs); - - /*! - * Constructs a symbolic linear equation solver with the given meta variable sets and pairs. - * - * @param A The matrix defining the coefficients of the linear equation system. - * @param allRows A BDD characterizing all rows of the equation system. - * @param illegalMask A mask that characterizes all illegal choices (that are therefore not to be taken). - * @param rowMetaVariables The meta variables used to encode the rows of the matrix. - * @param columnMetaVariables The meta variables used to encode the columns of the matrix. - * @param choiceVariables The variables encoding the choices of each row group. - * @param rowColumnMetaVariablePairs The pairs of row meta variables and the corresponding column meta - * variables. - * @param precision The precision to achieve. - * @param maximalNumberOfIterations The maximal number of iterations to perform when solving a linear - * equation system iteratively. - * @param relative Sets whether or not to use a relativ stopping criterion rather than an absolute one. - */ - SymbolicMinMaxLinearEquationSolver(storm::dd::Add const& A, storm::dd::Bdd const& allRows, storm::dd::Bdd const& illegalMask, std::set const& rowMetaVariables, std::set const& columnMetaVariables, std::set const& choiceVariables, std::vector> const& rowColumnMetaVariablePairs, double precision, uint_fast64_t maximalNumberOfIterations, bool relative); - + SymbolicMinMaxLinearEquationSolver(storm::dd::Add const& A, storm::dd::Bdd const& allRows, storm::dd::Bdd const& illegalMask, std::set const& rowMetaVariables, std::set const& columnMetaVariables, std::set const& choiceVariables, std::vector> const& rowColumnMetaVariablePairs, std::unique_ptr>&& linearEquationSolverFactory, SymbolicMinMaxLinearEquationSolverSettings const& settings = SymbolicMinMaxLinearEquationSolverSettings()); + /*! * Solves the equation system A*x = b. The matrix A is required to be square and have a unique solution. * The solution of the set of linear equations will be written to the vector x. Note that the matrix A has @@ -89,7 +100,12 @@ namespace storm { */ virtual storm::dd::Add multiply(bool minimize, storm::dd::Add const& x, storm::dd::Add const* b = nullptr, uint_fast64_t n = 1) const; - virtual ValueType const& getPrecision() const;// override; + SymbolicMinMaxLinearEquationSolverSettings const& getSettings() const; + + private: + storm::dd::Add solveEquationsValueIteration(bool minimize, storm::dd::Add const& x, storm::dd::Add const& b) const; + storm::dd::Add solveEquationsPolicyIteration(bool minimize, storm::dd::Add const& x, storm::dd::Add const& b) const; + protected: // The matrix defining the coefficients of the linear equation system. storm::dd::Add A; @@ -106,20 +122,35 @@ namespace storm { // The column variables. std::set columnMetaVariables; - // The choice variables + // The choice variables. std::set choiceVariables; // The pairs of meta variables used for renaming. std::vector> rowColumnMetaVariablePairs; - // The precision to achieve. - ValueType precision; + // A factory for creating linear equation solvers when needed. + std::unique_ptr> linearEquationSolverFactory; + + // The settings to use. + SymbolicMinMaxLinearEquationSolverSettings settings; + }; + + template + class SymbolicMinMaxLinearEquationSolverFactory { + public: + virtual std::unique_ptr> create(storm::dd::Add const& A, storm::dd::Bdd const& allRows, storm::dd::Bdd const& illegalMask, std::set const& rowMetaVariables, std::set const& columnMetaVariables, std::set const& choiceVariables, std::vector> const& rowColumnMetaVariablePairs) const = 0; + }; + + template + class SymbolicGeneralMinMaxLinearEquationSolverFactory : public SymbolicMinMaxLinearEquationSolverFactory { + public: + virtual std::unique_ptr> create(storm::dd::Add const& A, storm::dd::Bdd const& allRows, storm::dd::Bdd const& illegalMask, std::set const& rowMetaVariables, std::set const& columnMetaVariables, std::set const& choiceVariables, std::vector> const& rowColumnMetaVariablePairs) const; - // The maximal number of iterations to perform. - uint_fast64_t maximalNumberOfIterations; + SymbolicMinMaxLinearEquationSolverSettings& getSettings(); + SymbolicMinMaxLinearEquationSolverSettings const& getSettings() const; - // A flag indicating whether a relative or an absolute stopping criterion is to be used. - bool relative; + private: + SymbolicMinMaxLinearEquationSolverSettings settings; }; } // namespace solver diff --git a/src/storm/solver/SymbolicNativeLinearEquationSolver.cpp b/src/storm/solver/SymbolicNativeLinearEquationSolver.cpp index 612d29427..4bd28f38b 100644 --- a/src/storm/solver/SymbolicNativeLinearEquationSolver.cpp +++ b/src/storm/solver/SymbolicNativeLinearEquationSolver.cpp @@ -60,11 +60,13 @@ namespace storm { template storm::dd::Add SymbolicNativeLinearEquationSolver::solveEquations(storm::dd::Add const& x, storm::dd::Add const& b) const { + storm::dd::DdManager& manager = this->A.getDdManager(); + // Start by computing the Jacobi decomposition of the matrix A. storm::dd::Bdd diagonal = storm::utility::dd::getRowColumnDiagonal(x.getDdManager(), this->rowColumnMetaVariablePairs); diagonal &= this->allRows; - storm::dd::Add lu = diagonal.ite(this->A.getDdManager().template getAddZero(), this->A); + storm::dd::Add lu = diagonal.ite(manager.template getAddZero(), this->A); storm::dd::Add diagonalAdd = diagonal.template toAdd(); storm::dd::Add diag = diagonalAdd.multiplyMatrix(this->A, this->columnMetaVariables); diff --git a/src/storm/storage/dd/sylvan/InternalSylvanAdd.cpp b/src/storm/storage/dd/sylvan/InternalSylvanAdd.cpp index e44433c3f..dfa171eef 100644 --- a/src/storm/storage/dd/sylvan/InternalSylvanAdd.cpp +++ b/src/storm/storage/dd/sylvan/InternalSylvanAdd.cpp @@ -475,6 +475,11 @@ namespace storm { return InternalBdd(ddManager, this->sylvanMtbdd.AbstractMinRepresentative(cube.sylvanBdd)); } + template<> + InternalBdd InternalAdd::minAbstractRepresentative(InternalBdd const& cube) const { + return InternalBdd(ddManager, this->sylvanMtbdd.AbstractMinRepresentativeRN(cube.sylvanBdd)); + } + template<> InternalAdd InternalAdd::minAbstract(InternalBdd const& cube) const { return InternalAdd(ddManager, this->sylvanMtbdd.AbstractMinRN(cube.sylvanBdd)); @@ -497,6 +502,11 @@ namespace storm { return InternalBdd(ddManager, this->sylvanMtbdd.AbstractMaxRepresentative(cube.sylvanBdd)); } + template<> + InternalBdd InternalAdd::maxAbstractRepresentative(InternalBdd const& cube) const { + return InternalBdd(ddManager, this->sylvanMtbdd.AbstractMaxRepresentativeRN(cube.sylvanBdd)); + } + template<> InternalAdd InternalAdd::maxAbstract(InternalBdd const& cube) const { return InternalAdd(ddManager, this->sylvanMtbdd.AbstractMaxRN(cube.sylvanBdd)); diff --git a/src/storm/utility/solver.cpp b/src/storm/utility/solver.cpp index 80a1db702..6453a6668 100644 --- a/src/storm/utility/solver.cpp +++ b/src/storm/utility/solver.cpp @@ -25,13 +25,6 @@ namespace storm { namespace utility { namespace solver { - - - template - std::unique_ptr> SymbolicMinMaxLinearEquationSolverFactory::create(storm::dd::Add const& A, storm::dd::Bdd const& allRows, storm::dd::Bdd const& illegalMask, std::set const& rowMetaVariables, std::set const& columnMetaVariables, std::set const& choiceVariables, std::vector> const& rowColumnMetaVariablePairs) const { - return std::unique_ptr>(new storm::solver::SymbolicMinMaxLinearEquationSolver(A, allRows, illegalMask, rowMetaVariables, columnMetaVariables, choiceVariables, rowColumnMetaVariablePairs)); - } - template std::unique_ptr> SymbolicGameSolverFactory::create(storm::dd::Add const& A, storm::dd::Bdd const& allRows, storm::dd::Bdd const& illegalPlayer1Mask, storm::dd::Bdd const& illegalPlayer2Mask, std::set const& rowMetaVariables, std::set const& columnMetaVariables, std::vector> const& rowColumnMetaVariablePairs, std::set const& player1Variables, std::set const& player2Variables) const { return std::unique_ptr>(new storm::solver::SymbolicGameSolver(A, allRows, illegalPlayer1Mask, illegalPlayer2Mask, rowMetaVariables, columnMetaVariables, rowColumnMetaVariablePairs, player1Variables, player2Variables)); @@ -100,9 +93,6 @@ namespace storm { return factory->create(manager); } - template class SymbolicMinMaxLinearEquationSolverFactory; - template class SymbolicMinMaxLinearEquationSolverFactory; - template class SymbolicMinMaxLinearEquationSolverFactory; template class SymbolicGameSolverFactory; template class SymbolicGameSolverFactory; template class GameSolverFactory; diff --git a/src/storm/utility/solver.h b/src/storm/utility/solver.h index 1d3465e8c..074d2c9bd 100644 --- a/src/storm/utility/solver.h +++ b/src/storm/utility/solver.h @@ -57,12 +57,6 @@ namespace storm { namespace utility { namespace solver { - template - class SymbolicMinMaxLinearEquationSolverFactory { - public: - virtual std::unique_ptr> create(storm::dd::Add const& A, storm::dd::Bdd const& allRows, storm::dd::Bdd const& illegalMask, std::set const& rowMetaVariables, std::set const& columnMetaVariables, std::set const& choiceVariables, std::vector> const& rowColumnMetaVariablePairs) const; - }; - template class SymbolicGameSolverFactory { public: diff --git a/src/test/modelchecker/SymbolicMdpPrctlModelCheckerTest.cpp b/src/test/modelchecker/SymbolicMdpPrctlModelCheckerTest.cpp index 872a8dbdf..d5906064e 100644 --- a/src/test/modelchecker/SymbolicMdpPrctlModelCheckerTest.cpp +++ b/src/test/modelchecker/SymbolicMdpPrctlModelCheckerTest.cpp @@ -41,7 +41,7 @@ TEST(SymbolicMdpPrctlModelCheckerTest, Dice_Cudd) { std::shared_ptr> mdp = model->as>(); - storm::modelchecker::SymbolicMdpPrctlModelChecker> checker(*mdp, std::unique_ptr>(new storm::utility::solver::SymbolicMinMaxLinearEquationSolverFactory())); + storm::modelchecker::SymbolicMdpPrctlModelChecker> checker(*mdp, std::unique_ptr>(new storm::solver::SymbolicGeneralMinMaxLinearEquationSolverFactory())); std::shared_ptr formula = formulaParser.parseSingleFormulaFromString("Pmin=? [F \"two\"]"); @@ -139,7 +139,7 @@ TEST(SymbolicMdpPrctlModelCheckerTest, Dice_Sylvan) { std::shared_ptr> mdp = model->as>(); - storm::modelchecker::SymbolicMdpPrctlModelChecker> checker(*mdp, std::unique_ptr>(new storm::utility::solver::SymbolicMinMaxLinearEquationSolverFactory())); + storm::modelchecker::SymbolicMdpPrctlModelChecker> checker(*mdp, std::unique_ptr>(new storm::solver::SymbolicGeneralMinMaxLinearEquationSolverFactory())); std::shared_ptr formula = formulaParser.parseSingleFormulaFromString("Pmin=? [F \"two\"]"); @@ -237,7 +237,7 @@ TEST(SymbolicMdpPrctlModelCheckerTest, AsynchronousLeader_Cudd) { std::shared_ptr> mdp = model->as>(); - storm::modelchecker::SymbolicMdpPrctlModelChecker> checker(*mdp, std::unique_ptr>(new storm::utility::solver::SymbolicMinMaxLinearEquationSolverFactory())); + storm::modelchecker::SymbolicMdpPrctlModelChecker> checker(*mdp, std::unique_ptr>(new storm::solver::SymbolicGeneralMinMaxLinearEquationSolverFactory())); std::shared_ptr formula = formulaParser.parseSingleFormulaFromString("Pmin=? [F \"elected\"]"); @@ -317,7 +317,7 @@ TEST(SymbolicMdpPrctlModelCheckerTest, AsynchronousLeader_Sylvan) { std::shared_ptr> mdp = model->as>(); - storm::modelchecker::SymbolicMdpPrctlModelChecker> checker(*mdp, std::unique_ptr>(new storm::utility::solver::SymbolicMinMaxLinearEquationSolverFactory())); + storm::modelchecker::SymbolicMdpPrctlModelChecker> checker(*mdp, std::unique_ptr>(new storm::solver::SymbolicGeneralMinMaxLinearEquationSolverFactory())); std::shared_ptr formula = formulaParser.parseSingleFormulaFromString("Pmin=? [F \"elected\"]");