sp
/
tempest
1
0
Fork 0
Code Issues Pull Requests Projects Releases Wiki Activity
Browse Source

Renamed and refactored the helpers a little 

Former-commit-id: e2cd1d76eb
main
TimQu 9 years ago
parent
d5f1e6d5e7
commit
c86c6953b5
18 changed files with 621 additions and 1288 deletions
Whitespace
Show all changes Ignore whitespace when comparing lines Ignore changes in amount of whitespace Ignore changes in whitespace at EOL
Split View
Diff Options
Show Stats Download Patch File Download Diff File
  1. 160
      examples/multi-objective/mdp/power/power-timed.nm
  2. 11
      examples/multi-objective/mdp/power/power-timed.pctl
  3. 169
      examples/multi-objective/mdp/zeroconf_time/zeroconf_host_multi2_time.nm
  4. 174
      examples/multi-objective/mdp/zeroconf_time/zeroconf_host_multi4_time.nm
  5. 13
      examples/multi-objective/mdp/zeroconf_time/zeroconf_host_multi_time.pctl
  6. 18
      src/modelchecker/multiobjective/SparseMdpMultiObjectiveModelChecker.cpp
  7. 71
      src/modelchecker/multiobjective/helper/SparseMdpMultiObjectivePreprocessingHelper.h
  8. 327
      src/modelchecker/multiobjective/helper/SparseMultiObjectiveHelper.cpp
  9. 67
      src/modelchecker/multiobjective/helper/SparseMultiObjectiveHelper.h
  10. 56
      src/modelchecker/multiobjective/helper/SparseMultiObjectiveHelperReturnType.h
  11. 490
      src/modelchecker/multiobjective/helper/SparseMultiObjectiveModelCheckerHelper.cpp
  12. 80
      src/modelchecker/multiobjective/helper/SparseMultiObjectiveModelCheckerHelper.h
  13. 154
      src/modelchecker/multiobjective/helper/SparseMultiObjectivePreprocessor.cpp
  14. 67
      src/modelchecker/multiobjective/helper/SparseMultiObjectivePreprocessor.h
  15. 26
      src/modelchecker/multiobjective/helper/SparseMultiObjectivePreprocessorReturnType.h
  16. 2
      src/modelchecker/multiobjective/helper/SparseMultiObjectiveWeightVectorChecker.cpp
  17. 12
      src/modelchecker/multiobjective/helper/SparseMultiObjectiveWeightVectorChecker.h
  18. 12
      src/utility/vector.h

160
examples/multi-objective/mdp/power/power-timed.nm
View File

@ -1,160 +0,0 @@
// power manager example
mdp
const int QMAX; // max queue size
// to model the pm making a choice and then a move being made we need
// two clock ticks for each transition
// first the pm decides tick1 and then the system moves tick2
module timer
c : [0..1];
[tick1] c=0 -> (c'=1);
[tick2] c=1 -> (c'=0);
endmodule
//-------------------------------------------------------------------------
// POWER MANAGER
module PM
pm : [0..4] init 4;
// 0 - go to active
// 1 - go to idle
// 2 - go to idlelp
// 3 - go to stby
// 4 - go to sleep
[tick1] true -> (pm'=0);
[tick1] true -> (pm'=1);
[tick1] true -> (pm'=2);
[tick1] true -> (pm'=3);
[tick1] true -> (pm'=4);
endmodule
//-------------------------------------------------------------------------
// SERVICE REQUESTER
module SR
sr : [0..1] init 0;
// 0 idle
// 1 1req
[tick2] sr=0 -> 0.898: (sr'=0) + 0.102: (sr'=1);
[tick2] sr=1 -> 0.454: (sr'=0) + 0.546: (sr'=1);
endmodule
//-------------------------------------------------------------------------
// SERVICE PROVIDER
module SP
sp : [0..10] init 9;
// 0 active
// 1 idle
// 2 active_idlelp
// 3 idlelp
// 4 idlelp_active
// 5 active_stby
// 6 stby
// 7 stby_active
// 8 active_sleep
// 9 sleep
// 10 sleep_active
// states where PM has no control (transient states)
[tick2] sp=2 -> 0.75 : (sp'=2) + 0.25 : (sp'=3); // active_idlelp
[tick2] sp=4 -> 0.25 : (sp'=0) + 0.75 : (sp'=4); // idlelp_active
[tick2] sp=5 -> 0.995 : (sp'=5) + 0.005 : (sp'=6); // active_stby
[tick2] sp=7 -> 0.005 : (sp'=0) + 0.995 : (sp'=7); // stby_active
[tick2] sp=8 -> 0.9983 : (sp'=8) + 0.0017 : (sp'=9); // active_sleep
[tick2] sp=10 -> 0.0017 : (sp'=0) + 0.9983 : (sp'=10); // sleep_active
// states where PM has control
// goto_active
[tick2] sp=0 & pm=0 -> (sp'=0); // active
[tick2] sp=1 & pm=0 -> (sp'=0); // idle
[tick2] sp=3 & pm=0 -> (sp'=4); // idlelp
[tick2] sp=6 & pm=0 -> (sp'=7); // stby
[tick2] sp=9 & pm=0 -> (sp'=10); // sleep
// goto_idle
[tick2] sp=0 & pm=1 -> (sp'=1); // active
[tick2] sp=1 & pm=1 -> (sp'=1); // idle
[tick2] sp=3 & pm=1 -> (sp'=3); // idlelp
[tick2] sp=6 & pm=1 -> (sp'=6); // stby
[tick2] sp=9 & pm=1 -> (sp'=9); // sleep
// goto_idlelp
[tick2] sp=0 & pm=2 -> (sp'=2); // active
[tick2] sp=1 & pm=2 -> (sp'=2); // idle
[tick2] sp=3 & pm=2 -> (sp'=3); // idlelp
[tick2] sp=6 & pm=2 -> (sp'=6); // stby
[tick2] sp=9 & pm=2 -> (sp'=9); // sleep
// goto_stby
[tick2] sp=0 & pm=3 -> (sp'=5); // active
[tick2] sp=1 & pm=3 -> (sp'=5); // idle
[tick2] sp=3 & pm=3 -> (sp'=5); // idlelp
[tick2] sp=6 & pm=3 -> (sp'=6); // stby
[tick2] sp=9 & pm=3 -> (sp'=9); // sleep
// goto_sleep
[tick2] sp=0 & pm=4 -> (sp'=8); // active
[tick2] sp=1 & pm=4 -> (sp'=8); // idle
[tick2] sp=3 & pm=4 -> (sp'=8); // idlelp
[tick2] sp=6 & pm=4 -> (sp'=8); // stby
[tick2] sp=9 & pm=4 -> (sp'=9); // sleep
endmodule
//-------------------------------------------------------------------------
// SQ
module SQ
q : [0..QMAX] init 0;
// serve if busy
[tick2] sr=0 & sp=0 -> (q'=max(q-1,0));
[tick2] sr=1 & sp=0 -> (q'=q);
// otherwise do nothing
[tick2] sr=0 & sp>0 -> (q'=q);
[tick2] sr=1 & sp>0 -> (q'=min(q+1,QMAX));
endmodule
//-------------------------------------------------------------------------
//rewards "time"
// [tick2] bat=1 : 1;
//endrewards
rewards "power"
[tick2] sp=0 & c=1 : 2.5;
[tick2] sp=1 & c=1 : 1.5;
[tick2] sp=2 & c=1 : 2.5;
[tick2] sp=3 & c=1 : 0.8;
[tick2] sp=4 & c=1 : 2.5;
[tick2] sp=5 & c=1 : 2.5;
[tick2] sp=6 & c=1 : 0.3;
[tick2] sp=7 & c=1 : 2.5;
[tick2] sp=8 & c=1 : 2.5;
[tick2] sp=9 & c=1 : 0.1;
[tick2] sp=10 & c=1 : 2.5;
endrewards
// is an instantaneous property but I suppose we can look at average size
// i.e. divide by the expected number of time steps
rewards "queue"
[tick2] c=1 : q;
endrewards
rewards "lost"
[tick2] sr=1 & sp>0 & q=2 : 1;
endrewards

11
examples/multi-objective/mdp/power/power-timed.pctl
View File

@ -1,11 +0,0 @@
// Average queue size
const double Q;
// Time bound
const int k;
// Minimum energy usage over k time-steps, such that average queue size remains below Q
"num_energy": multi(R{"power"}min=? [ C<=k ], R{"queue"}<=Q*k [ C<=k ])
// Pareto query: minimum energy usage vs minimum average queue size
"pareto": multi(R{"power"}min=? [ C<=k ], R{"queue"}min=? [ C<=k ])

169
examples/multi-objective/mdp/zeroconf_time/zeroconf_host_multi2_time.nm
View File

@ -1,169 +0,0 @@
// IPv4: PTA model with digitial clocks
// multi-objective model of the host
// gxn/dxp 28/09/09
mdp
//-------------------------------------------------------------
// VARIABLES
const int N=20; // number of abstract hosts
const int K=2; // number of probes to send
// PROBABILITIES
const double old = N/65024; // probability pick an ip address being used
//const double old = 1/2; // probability pick an ip address being used
const double new = (1-old); // probability pick a new ip address
// TIMING CONSTANTS
const int CONSEC = 2; // time interval between sending consecutive probles
const int TRANSTIME = 1; // upper bound on transmission time delay
const int LONGWAIT = 60; // minimum time delay after a high number of address collisions
const int DEFEND = 10;
const int TIME_MAX_X = 60; // max value of clock x
const int TIME_MAX_Y = 10; // max value of clock y
const int TIME_MAX_Z = 1; // max value of clock z
// OTHER CONSTANTS
const int MAXCOLL = 10; // maximum number of collisions before long wait
//-------------------------------------------------------------
// CONCRETE HOST
module host0
x : [0..TIME_MAX_X]; // first clock of the host
y : [0..TIME_MAX_Y]; // second clock of the host
coll : [0..MAXCOLL]; // number of address collisions
probes : [0..K]; // counter (number of probes sent)
mess : [0..1]; // need to send a message or not
defend : [0..1]; // defend (if =1, try to defend IP address)
ip : [1..2]; // ip address (1 - in use & 2 - fresh)
l : [0..4] init 1; // location
// 0 : RECONFIGURE
// 1 : RANDOM
// 2 : WAITSP
// 3 : WAITSG
// 4 : USE
// RECONFIGURE
[reset] l=0 -> (l'=1);
// RANDOM (choose IP address)
[rec0] (l=1) -> true; // get message (ignore since have no ip address)
[rec1] (l=1) -> true; // get message (ignore since have no ip address)
// small number of collisions (choose straight away)
[] l=1 & coll<MAXCOLL -> 1/3*old : (l'=2) & (ip'=1) & (x'=0)
+ 1/3*old : (l'=2) & (ip'=1) & (x'=1)
+ 1/3*old : (l'=2) & (ip'=1) & (x'=2)
+ 1/3*new : (l'=2) & (ip'=2) & (x'=0)
+ 1/3*new : (l'=2) & (ip'=2) & (x'=1)
+ 1/3*new : (l'=2) & (ip'=2) & (x'=2);
// large number of collisions: (wait for LONGWAIT)
[time] l=1 & coll=MAXCOLL & x<LONGWAIT -> (x'=min(x+1,TIME_MAX_X));
[] l=1 & coll=MAXCOLL & x=LONGWAIT -> 1/3*old : (l'=2) & (ip'=1) & (x'=0)
+ 1/3*old : (l'=2) & (ip'=1) & (x'=1)
+ 1/3*old : (l'=2) & (ip'=1) & (x'=2)
+ 1/3*new : (l'=2) & (ip'=2) & (x'=0)
+ 1/3*new : (l'=2) & (ip'=2) & (x'=1)
+ 1/3*new : (l'=2) & (ip'=2) & (x'=2);
// WAITSP
// let time pass
[time] l=2 & x<2 -> (x'=min(x+1,2));
// send probe
[send1] l=2 & ip=1 & x=2 & probes<K -> (x'=0) & (probes'=probes+1);
[send2] l=2 & ip=2 & x=2 & probes<K -> (x'=0) & (probes'=probes+1);
// sent K probes and waited 2 seconds
[configured] l=2 & x=2 & probes=K -> (l'=3) & (probes'=0) & (coll'=0) & (x'=0);
// get message and ip does not match: ignore
[rec0] l=2 & ip!=0 -> (l'=l);
[rec1] l=2 & ip!=1 -> (l'=l);
// get a message with matching ip: reconfigure
[rec1] l=2 & ip=1 -> (l'=0) & (coll'=min(coll+1,MAXCOLL)) & (x'=0) & (probes'=0);
// WAITSG (sends two gratuitious arp probes)
// time passage
[time] l=3 & mess=0 & defend=0 & x<CONSEC -> (x'=min(x+1,TIME_MAX_X));
[time] l=3 & mess=0 & defend=1 & x<CONSEC -> (x'=min(x+1,TIME_MAX_X)) & (y'=min(y+1,DEFEND));
// receive message and same ip: defend
[rec1] l=3 & mess=0 & ip=1 & (defend=0 | y>=DEFEND) -> (defend'=1) & (mess'=1) & (y'=0);
// receive message and same ip: defer
[rec1] l=3 & mess=0 & ip=1 & (defend=0 | y<DEFEND) -> (l'=0) & (probes'=0) & (defend'=0) & (x'=0) & (y'=0);
// receive message and different ip
[rec0] l=3 & mess=0 & ip!=0 -> (l'=l);
[rec1] l=3 & mess=0 & ip!=1 -> (l'=l);
// send probe reply or message for defence
[send1] l=3 & ip=1 & mess=1 -> (mess'=0);
[send2] l=3 & ip=2 & mess=1 -> (mess'=0);
// send first gratuitous arp message
[send1] l=3 & ip=1 & mess=0 & x=CONSEC & probes<1 -> (x'=0) & (probes'=probes+1);
[send2] l=3 & ip=2 & mess=0 & x=CONSEC & probes<1 -> (x'=0) & (probes'=probes+1);
// send second gratuitous arp message (move to use)
[send1] l=3 & ip=1 & mess=0 & x=CONSEC & probes=1 -> (l'=4) & (x'=0) & (y'=0) & (probes'=0);
[send2] l=3 & ip=2 & mess=0 & x=CONSEC & probes=1 -> (l'=4) & (x'=0) & (y'=0) & (probes'=0);
// USE (only interested in reaching this state so do not need to add anything here)
[] l=4 -> true;
endmodule
//-------------------------------------------------------------
// error automaton for the environment assumption
// do not get a reply when K probes are sent
const int M; // time between sending and receiving a message
module env_error2
env : [0..1]; // 0 active and 1 done
k : [0..2]; // counts the number of messages sent
c1 : [0..M+1]; // time since first message
c2 : [0..M+1]; // time since second message
error : [0..1];
// message with new ip address arrives so done
[send2] error=0 & env=0 -> (env'=1);
// message with old ip address arrives so count
[send1] error=0 & env=0 -> (k'=min(k+1,K));
// time passgae so update relevant clocks
[time] error=0 & env=0 & k=0 -> true;
[time] error=0 & env=0 & k=1 & min(c1,c2)<M -> (c1'=min(c1+1,M+1));
[time] error=0 & env=0 & k=2 & min(c1,c2)<M -> (c1'=min(c1+1,M+1)) & (c2'=min(c2+1,M+1));
// all clocks reached their bound so an error
[time] error=0 & env=0 & min(c1,c2)=M -> (error'=1);
// send a reply (then done)
[rec1] error=0 & env=0 & k>0 & min(c1,c2)<=M -> (env'=1);
// finished so any action can be performed
[time] error=1 | env=1 -> true;
[send1] error=1 | env=1 -> true;
[send2] error=1 | env=1 -> true;
[rec1] error=1 | env=1 -> true;
endmodule
//-------------------------------------------------------------
// error automaton for the time bounded assumption
// host does not send configured signal within T seconds
const int T;
module time_error
time_error : [0..1];
done : [0..1];
t : [0..T];
[time] t<T-1 & done=0 & time_error=0 -> (t'=t+1); // time passes and bound not reached
[time] t=T-1 & done=0 & time_error=0 -> (time_error'=1); // bound reached so error
[configured] time_error=0 -> (done'=1); // configured within the time bound
// when in error or done state can loop with either action
[configured] time_error=1 | done=1 -> true;
[time] time_error=1 | done=1 -> true;
endmodule

174
examples/multi-objective/mdp/zeroconf_time/zeroconf_host_multi4_time.nm
View File

@ -1,174 +0,0 @@
// IPv4: PTA model with digitial clocks
// multi-objective model of the host
// gxn/dxp 28/09/09
mdp
//-------------------------------------------------------------
// VARIABLES
const int N=20; // number of abstract hosts
const int K=4; // number of probes to send
// PROBABILITIES
const double old = N/65024; // probability pick an ip address being used
//const double old = 0.5; // probability pick an ip address being used
const double new = (1-old); // probability pick a new ip address
// TIMING CONSTANTS
const int CONSEC = 2; // time interval between sending consecutive probles
const int TRANSTIME = 1; // upper bound on transmission time delay
const int LONGWAIT = 60; // minimum time delay after a high number of address collisions
const int DEFEND = 10;
const int TIME_MAX_X = 60; // max value of clock x
const int TIME_MAX_Y = 10; // max value of clock y
const int TIME_MAX_Z = 1; // max value of clock z
// OTHER CONSTANTS
const int MAXCOLL = 10; // maximum number of collisions before long wait
//-------------------------------------------------------------
// CONCRETE HOST
module host0
x : [0..TIME_MAX_X]; // first clock of the host
y : [0..TIME_MAX_Y]; // second clock of the host
coll : [0..MAXCOLL]; // number of address collisions
probes : [0..K]; // counter (number of probes sent)
mess : [0..1]; // need to send a message or not
defend : [0..1]; // defend (if =1, try to defend IP address)
ip : [1..2]; // ip address (1 - in use & 2 - fresh)
l : [0..4] init 1; // location
// 0 : RECONFIGURE
// 1 : RANDOM
// 2 : WAITSP
// 3 : WAITSG
// 4 : USE
// RECONFIGURE
[reset] l=0 -> (l'=1);
// RANDOM (choose IP address)
[rec0] (l=1) -> true; // get message (ignore since have no ip address)
[rec1] (l=1) -> true; // get message (ignore since have no ip address)
// small number of collisions (choose straight away)
[] l=1 & coll<MAXCOLL -> 1/3*old : (l'=2) & (ip'=1) & (x'=0)
+ 1/3*old : (l'=2) & (ip'=1) & (x'=1)
+ 1/3*old : (l'=2) & (ip'=1) & (x'=2)
+ 1/3*new : (l'=2) & (ip'=2) & (x'=0)
+ 1/3*new : (l'=2) & (ip'=2) & (x'=1)
+ 1/3*new : (l'=2) & (ip'=2) & (x'=2);
// large number of collisions: (wait for LONGWAIT)
[time] l=1 & coll=MAXCOLL & x<LONGWAIT -> (x'=min(x+1,TIME_MAX_X));
[] l=1 & coll=MAXCOLL & x=LONGWAIT -> 1/3*old : (l'=2) & (ip'=1) & (x'=0)
+ 1/3*old : (l'=2) & (ip'=1) & (x'=1)
+ 1/3*old : (l'=2) & (ip'=1) & (x'=2)
+ 1/3*new : (l'=2) & (ip'=2) & (x'=0)
+ 1/3*new : (l'=2) & (ip'=2) & (x'=1)
+ 1/3*new : (l'=2) & (ip'=2) & (x'=2);
// WAITSP
// let time pass
[time] l=2 & x<2 -> (x'=min(x+1,2));
// send probe
[send1] l=2 & ip=1 & x=2 & probes<K -> (x'=0) & (probes'=probes+1);
[send2] l=2 & ip=2 & x=2 & probes<K -> (x'=0) & (probes'=probes+1);
// sent K probes and waited 2 seconds
[configured] l=2 & x=2 & probes=K -> (l'=3) & (probes'=0) & (coll'=0) & (x'=0);
// get message and ip does not match: ignore
[rec0] l=2 & ip!=0 -> (l'=l);
[rec1] l=2 & ip!=1 -> (l'=l);
// get a message with matching ip: reconfigure
[rec1] l=2 & ip=1 -> (l'=0) & (coll'=min(coll+1,MAXCOLL)) & (x'=0) & (probes'=0);
// WAITSG (sends two gratuitious arp probes)
// time passage
[time] l=3 & mess=0 & defend=0 & x<CONSEC -> (x'=min(x+1,TIME_MAX_X));
[time] l=3 & mess=0 & defend=1 & x<CONSEC -> (x'=min(x+1,TIME_MAX_X)) & (y'=min(y+1,DEFEND));
// receive message and same ip: defend
[rec1] l=3 & mess=0 & ip=1 & (defend=0 | y>=DEFEND) -> (defend'=1) & (mess'=1) & (y'=0);
// receive message and same ip: defer
[rec1] l=3 & mess=0 & ip=1 & (defend=0 | y<DEFEND) -> (l'=0) & (probes'=0) & (defend'=0) & (x'=0) & (y'=0);
// receive message and different ip
[rec0] l=3 & mess=0 & ip!=0 -> (l'=l);
[rec1] l=3 & mess=0 & ip!=1 -> (l'=l);
// send probe reply or message for defence
[send1] l=3 & ip=1 & mess=1 -> (mess'=0);
[send2] l=3 & ip=2 & mess=1 -> (mess'=0);
// send first gratuitous arp message
[send1] l=3 & ip=1 & mess=0 & x=CONSEC & probes<1 -> (x'=0) & (probes'=probes+1);
[send2] l=3 & ip=2 & mess=0 & x=CONSEC & probes<1 -> (x'=0) & (probes'=probes+1);
// send second gratuitous arp message (move to use)
[send1] l=3 & ip=1 & mess=0 & x=CONSEC & probes=1 -> (l'=4) & (x'=0) & (y'=0) & (probes'=0);
[send2] l=3 & ip=2 & mess=0 & x=CONSEC & probes=1 -> (l'=4) & (x'=0) & (y'=0) & (probes'=0);
// USE (only interested in reaching this state so do not need to add anything here)
[] l=4 -> true;
endmodule
//-------------------------------------------------------------
// error automaton for the environment assumption
// do not get a reply when K probes are sent
const int M; // time between sending and receiving a message
module env_error4
env : [0..1]; // 0 active and 1 done
k : [0..4]; // counts the number of messages sent
c1 : [0..M+1]; // time since first message
c2 : [0..M+1]; // time since second message
c3 : [0..M+1]; // time since third message
c4 : [0..M+1]; // time since fourth message
error : [0..1];
// message with new ip address arrives so done
[send2] error=0 & env=0 -> (env'=1);
// message with old ip address arrives so count
[send1] error=0 & env=0 -> (k'=min(k+1,K));
// time passgae so update relevant clocks
[time] error=0 & env=0 & k=0 -> true;
[time] error=0 & env=0 & k=1 & min(c1,c2,c3,c4)<M -> (c1'=min(c1+1,M+1));
[time] error=0 & env=0 & k=2 & min(c1,c2,c3,c4)<M -> (c1'=min(c1+1,M+1)) & (c2'=min(c2+1,M+1));
[time] error=0 & env=0 & k=3 & min(c1,c2,c3,c4)<M -> (c1'=min(c1+1,M+1)) & (c2'=min(c2+1,M+1)) & (c3'=min(c3+1,M+1));
[time] error=0 & env=0 & k=4 & min(c1,c2,c3,c4)<M -> (c1'=min(c1+1,M+1)) & (c2'=min(c2+1,M+1)) & (c3'=min(c3+1,M+1)) & (c4'=min(c4+1,M+1));
// all clocks reached their bound so an error
[time] error=0 & env=0 & min(c1,c2,c3,c4)=M -> (error'=1);
// send a reply (then done)
[rec1] error=0 & env=0 & k>0 & min(c1,c2,c3,c4)<=M -> (env'=1);
// finished so any action can be performed
[time] error=1 | env=1 -> true;
[send1] error=1 | env=1 -> true;
[send2] error=1 | env=1 -> true;
[send2] error=1 | env=1 -> true;
[rec1] error=1 | env=1 -> true;
endmodule
//-------------------------------------------------------------
// error automaton for the time bounded assumption
// host does not send configured signal within T seconds
const int T;
module time_error
time_error : [0..1];
done : [0..1];
t : [0..T];
[time] t<T-1 & done=0 & time_error=0 -> (t'=t+1); // time passes and bound not reached
[time] t=T-1 & done=0 & time_error=0 -> (time_error'=1); // bound reached so error
[configured] time_error=0 -> (done'=1); // configured within the time bound
// when in error or done state can loop with either action
[configured] time_error=1 | done=1 -> true;
[time] time_error=1 | done=1 -> true;
endmodule

13
examples/multi-objective/mdp/zeroconf_time/zeroconf_host_multi_time.pctl
View File

@ -1,13 +0,0 @@
// Max probability of component violating assumption property (checked separately)
const double p_fail =
K=2 ? 0.19 :
K=4 ? 0.006859000000000001 :
K=6 ? 2.476099000000001E-4 :
K=8 ? 8.938717390000006E-6 :
0;
// Assume-guarantee check via multi-objective
"num_ag": multi(Pmax=? [ F time_error=1 ] , P>=1-p_fail [ G (error=0) ])
// Pareto query for assume-guarantee check
"pareto": multi(Pmax=? [ F time_error=1 ] , Pmax=? [ G (error=0) ])

18
src/modelchecker/multiobjective/SparseMdpMultiObjectiveModelChecker.cpp
View File

@ -9,9 +9,8 @@
#include "src/modelchecker/results/ExplicitQualitativeCheckResult.h"
#include "src/modelchecker/results/ExplicitQuantitativeCheckResult.h"
#include "src/modelchecker/multiobjective/helper/SparseMultiObjectiveModelCheckerInformation.h"
#include "src/modelchecker/multiobjective/helper/SparseMdpMultiObjectivePreprocessingHelper.h"
#include "src/modelchecker/multiobjective/helper/SparseMultiObjectiveModelCheckerHelper.h"
#include "src/modelchecker/multiobjective/helper/SparseMultiObjectivePreprocessor.h"
#include "src/modelchecker/multiobjective/helper/SparseMultiObjectiveHelper.h"
#include "src/exceptions/NotImplementedException.h"
@ -35,20 +34,17 @@ namespace storm {
template<typename SparseMdpModelType>
std::unique_ptr<CheckResult> SparseMdpMultiObjectiveModelChecker<SparseMdpModelType>::checkMultiObjectiveFormula(CheckTask<storm::logic::MultiObjectiveFormula> const& checkTask) {
helper::SparseMultiObjectiveModelCheckerInformation<SparseMdpModelType> info = helper::SparseMdpMultiObjectivePreprocessingHelper<SparseMdpModelType>::preprocess(checkTask.getFormula(), this->getModel());
auto preprocessedData = helper::SparseMdpMultiObjectivePreprocessingHelper<SparseMdpModelType>::preprocess(checkTask.getFormula(), this->getModel());
std::cout << std::endl;
std::cout << "Preprocessed Information:" << std::endl;
info.printInformationToStream(std::cout);
std::cout << "Preprocessing done." << std::endl;
preprocessedData.printToStream(std::cout);
#ifdef STORM_HAVE_CARL
helper::SparseMultiObjectiveModelCheckerHelper<SparseMdpModelType, storm::RationalNumber>::check(info);
helper::SparseMultiObjectiveHelper<SparseMdpModelType, storm::RationalNumber>::check(preprocessedData);
#else
STORM_LOG_THROW(false, storm::exceptions::UnexpectedException, "Multi objective model checking currently requires carl.");
STORM_LOG_THROW(false, storm::exceptions::UnexpectedException, "Multi objective model checking requires carl.");
#endif
std::cout << "Information after helper call: " << std::endl;
info.printInformationToStream(std::cout);
return std::unique_ptr<CheckResult>(new ExplicitQualitativeCheckResult());
}

71
src/modelchecker/multiobjective/helper/SparseMdpMultiObjectivePreprocessingHelper.h
View File

@ -1,71 +0,0 @@
#ifndef STORM_MODELCHECKER_MULTIOBJECTIVE_HELPER_SPARSEMDPMULTIOBJECTIVEPREPROCESSINGHELPER_H_
#define STORM_MODELCHECKER_MULTIOBJECTIVE_HELPER_SPARSEMDPMULTIOBJECTIVEPREPROCESSINGHELPER_H_
#include "src/logic/Formulas.h"
#include "src/modelchecker/multiobjective/helper/SparseMultiObjectiveModelCheckerInformation.h"
namespace storm {
namespace modelchecker {
namespace helper {
/*
* Helper Class to invoke the necessary preprocessing for multi objective model checking
*/
template <class SparseMdpModelType>
class SparseMdpMultiObjectivePreprocessingHelper {
public:
typedef typename SparseMdpModelType::ValueType ValueType;
typedef typename SparseMdpModelType::RewardModelType RewardModelType;
typedef SparseMultiObjectiveModelCheckerInformation<SparseMdpModelType> Information;
/*!
* Preprocesses the given model w.r.t. the given formulas.
* @param originalModel The considered model
* @param originalFormula the considered formula. The subformulas should only contain one OperatorFormula at top level, i.e., the formula is simple.
*/
static Information preprocess(storm::logic::MultiObjectiveFormula const& originalFormula, SparseMdpModelType const& originalModel);
private:
/*!
* Inserts the information regarding the given formula (i.e. objective) into info.objectives
*
* @param formula OperatorFormula representing the objective
* @param the information collected so far
*/
static void addObjective(std::shared_ptr<storm::logic::Formula const> const& formula, Information& info);
/*!
* Sets the timeBound for the given objective
*/
static void setStepBoundOfObjective(typename Information::ObjectiveInformation& currentObjective);
/*!
* Sets whether we should consider negated rewards
*/
static void setWhetherNegatedRewardsAreConsidered(Information& info);
/*!
* Apply the neccessary preprocessing for the given formula.
* @param formula the current (sub)formula
* @param info the information collected so far
* @param currentObjective the currently considered objective. The given formula should be a a (sub)formula of this objective
* @param optionalRewardModelName the reward model name that is considered for the formula (if available)
*/
static void preprocessFormula(storm::logic::ProbabilityOperatorFormula const& formula, Information& info, typename Information::ObjectiveInformation& currentObjective);
static void preprocessFormula(storm::logic::RewardOperatorFormula const& formula, Information& info, typename Information::ObjectiveInformation& currentObjective);
static void preprocessFormula(storm::logic::UntilFormula const& formula, Information& info, typename Information::ObjectiveInformation& currentObjective);
static void preprocessFormula(storm::logic::BoundedUntilFormula const& formula, Information& info, typename Information::ObjectiveInformation& currentObjective);
static void preprocessFormula(storm::logic::GloballyFormula const& formula, Information& info, typename Information::ObjectiveInformation& currentObjective);
static void preprocessFormula(storm::logic::EventuallyFormula const& formula, Information& info, typename Information::ObjectiveInformation& currentObjective, boost::optional<std::string> const& optionalRewardModelName = boost::none);
static void preprocessFormula(storm::logic::CumulativeRewardFormula const& formula, Information& info, typename Information::ObjectiveInformation& currentObjective, boost::optional<std::string> const& optionalRewardModelName = boost::none);
};
}
}
}
#endif /* STORM_MODELCHECKER_MULTIOBJECTIVE_HELPER_SPARSEMDPMULTIOBJECTIVEPREPROCESSINGHELPER_H_ */

327
src/modelchecker/multiobjective/helper/SparseMultiObjectiveHelper.cpp
View File

@ -0,0 +1,327 @@
#include "src/modelchecker/multiobjective/helper/SparseMultiObjectiveHelper.h"
#include "src/adapters/CarlAdapter.h"
#include "src/models/sparse/Mdp.h"
#include "src/models/sparse/StandardRewardModel.h"
#include "src/utility/constants.h"
#include "src/utility/vector.h"
#include "src/settings//SettingsManager.h"
#include "src/settings/modules/MultiObjectiveSettings.h"
#include "src/exceptions/UnexpectedException.h"
namespace storm {
namespace modelchecker {
namespace helper {
template <class SparseModelType, typename RationalNumberType>
typename SparseMultiObjectiveHelper<SparseModelType, RationalNumberType>::ReturnType SparseMultiObjectiveHelper<SparseModelType, RationalNumberType>::check(PreprocessorData const& data) {
ReturnType result;
result.overApproximation = storm::storage::geometry::Polytope<RationalNumberType>::createUniversalPolytope();
result.underApproximation = storm::storage::geometry::Polytope<RationalNumberType>::createEmptyPolytope();
uint_fast64_t numOfObjectivesWithoutThreshold = 0;
for(auto const& obj : data.objectives) {
if(!obj.threshold) {
++numOfObjectivesWithoutThreshold;
}
}
if(numOfObjectivesWithoutThreshold == 0) {
achievabilityQuery(data, result);
} else if (numOfObjectivesWithoutThreshold == 1) {
numericalQuery(data, result);
} else if (numOfObjectivesWithoutThreshold == data.objectives.size()) {
paretoQuery(data, result);
} else {
STORM_LOG_THROW(false, storm::exceptions::UnexpectedException, "The number of objecties without threshold is not valid. It should be either 0 (achievabilityQuery), 1 (numericalQuery), or " << data.objectives.size() << " (paretoQuery). Got " << numOfObjectivesWithoutThreshold << " instead.");
}
return result;
}
template <class SparseModelType, typename RationalNumberType>
void SparseMultiObjectiveHelper<SparseModelType, RationalNumberType>::achievabilityQuery(PreprocessorData const& data, ReturnType& result) {
Point thresholds;
thresholds.reserve(data.objectives.size());
storm::storage::BitVector strictThresholds(data.objectives.size());
for(uint_fast64_t objIndex = 0; objIndex < data.objectives.size(); ++objIndex) {
thresholds.push_back(storm::utility::convertNumber<RationalNumberType>(*data.objectives[objIndex].threshold));
strictThresholds.set(objIndex, data.objectives[objIndex].thresholdIsStrict);
}
storm::storage::BitVector individualObjectivesToBeChecked(data.objectives.size(), true);
bool converged=false;
SparseMultiObjectiveWeightVectorChecker<SparseModelType> weightVectorChecker(data);
do {
WeightVector separatingVector = findSeparatingVector(thresholds, result.underApproximation, individualObjectivesToBeChecked);
refineResult(separatingVector, weightVectorChecker, result);
// Check for convergence
if(!checkIfThresholdsAreSatisfied(overApproximation, thresholds, strictThresholds)){
result.thresholdsAreAchievable = false;
converged=true;
}
if(checkIfThresholdsAreSatisfied(underApproximation, thresholds, strictThresholds)){
result.thresholdsAreAchievable = true;
converged=true;
}
converged |= (storm::settings::getModule<storm::settings::modules::MultiObjectiveSettings>().isMaxIterationsSet() && result.iterations.size() >= storm::settings::getModule<storm::settings::modules::MultiObjectiveSettings>().getMaxIterations());
} while(!converged);
}
template <class SparseModelType, typename RationalNumberType>
void SparseMultiObjectiveHelper<SparseModelType, RationalNumberType>::numericalQuery(PreprocessorData const& data, ReturnType& result) {
Point thresholds;
thresholds.reserve(data.objectives.size());
storm::storage::BitVector strictThresholds(data.objectives.size());
std::vector<storm::storage::geometry::Halfspace<RationalNumberType>> thresholdConstraints;
thresholdConstraints.reserve(data.objectives.size()-1);
uint_fast64_t optimizingObjIndex = data.objectives.size(); // init with invalid value...
for(uint_fast64_t objIndex = 0; objIndex < data.objectives.size(); ++objIndex) {
if(data.objectives[objIndex].threshold) {
thresholds.push_back(storm::utility::convertNumber<RationalNumberType>(*data.objectives[objIndex].threshold));
WeightVector normalVector(data.objectives.size(), storm::utility::zero<RationalNumberType>());
normalVector[objIndex] = -storm::utility::one<RationalNumberType>();
thresholdConstraints.emplace_back(std::move(normalVector), -thresholds.back());
strictThresholds.set(objIndex, data.objectives[objIndex].thresholdIsStrict);
} else {
optimizingObjIndex = objIndex;
thresholds.push_back(storm::utility::zero<RationalNumberType>());
}
}
auto thresholdsAsPolytope = storm::storage::geometry::Polytope<RationalNumberType>::create(thresholdConstraints);
WeightVector directionOfOptimizingObjective(data.objectives.size(), storm::utility::zero<RationalNumberType>());
directionOfOptimizingObjective[optimizingObjIndex] = storm::utility::one<RationalNumberType>();
// Try to find one valid solution
storm::storage::BitVector individualObjectivesToBeChecked(data.objectives.size(), true);
individualObjectivesToBeChecked.set(optimizingObjIndex, false);
SparseMultiObjectiveWeightVectorChecker<SparseModelType> weightVectorChecker(data);
bool converged=false;
do {
WeightVector separatingVector = findSeparatingVector(thresholds, result.underApproximation, individualObjectivesToBeChecked);
refineResult(separatingVector, weightVectorChecker, result);
//Pick the threshold for the optimizing objective low enough so valid solutions are not excluded
thresholds[optimizingObjIndex] = std::min(thresholds[optimizingObjIndex], iterations.back().point[optimizingObjIndex]);
if(!checkIfThresholdsAreSatisfied(overApproximation, thresholds, strictThresholds)){
result.thresholdsAreAchievable = false;
converged=true;
}
if(checkIfThresholdsAreSatisfied(underApproximation, thresholds, strictThresholds)){
result.thresholdsAreAchievable = true;
converged=true;
}
} while(!converged);
if(*result.thresholdsAreAchievable) {
STORM_LOG_DEBUG("Found a solution that satisfies the objective thresholds.");
individualObjectivesToBeChecked.clear();
// Improve the found solution.
// Note that we do not have to care whether a threshold is strict anymore, because the resulting optimum should be
// the supremum over all strategies. Hence, one could combine a scheduler inducing the optimum value (but possibly violating strict
// thresholds) and (with very low probability) a scheduler that satisfies all (possibly strict) thresholds.
converged = false;
do {
std::pair<Point, bool> optimizationRes = underApproximation->intersection(thresholdsAsPolytope)->optimize(directionOfOptimizingObjective);
STORM_LOG_THROW(optimizationRes.second, storm::exceptions::UnexpectedException, "The underapproximation is either unbounded or empty.");
thresholds[optimizingObjIndex] = optimizationRes.first[optimizingObjIndex];
STORM_LOG_DEBUG("Best solution found so far is " << storm::utility::convertNumber<double>(thresholds[optimizingObjIndex]) << ".");
//Compute an upper bound for the optimum and check for convergence
optimizationRes = overApproximation->intersection(thresholdsAsPolytope)->optimize(directionOfOptimizingObjective);
if(optimizationRes.second) {
result.precisionOfResult = optimizationRes.first[optimizingObjIndex]; - thresholds[optimizingObjIndex];
STORM_LOG_DEBUG("Solution can be improved by at most " << storm::utility::convertNumber<double>(*result.precisionOfResult));
if(storm::utility::convertNumber<double>(*result.precisionOfResult) < storm::settings::getModule<storm::settings::modules::MultiObjectiveSettings>().getPrecision()) {
result.numericalResult = thresholds[optimizingObjIndex];
result.optimumIsAchievable = checkIfThresholdsAreSatisfied(underApproximation, thresholds, strictThresholds);
converged=true;
}
}
converged |= (storm::settings::getModule<storm::settings::modules::MultiObjectiveSettings>().isMaxIterationsSet() && result.iterations.size() >= storm::settings::getModule<storm::settings::modules::MultiObjectiveSettings>().getMaxIterations());
if(!converged) {
WeightVector separatingVector = findSeparatingVector(thresholds, result.underApproximation, individualObjectivesToBeChecked);
refineResult(separatingVector, weightVectorChecker, result);
}
} while(!converged);
}
}
template <class SparseModelType, typename RationalNumberType>
void SparseMultiObjectiveHelper<SparseModelType, RationalNumberType>::paretoQuery(PreprocessorData const& data, ReturnType& result) {
//First consider the objectives individually
for(uint_fast64_t objIndex = 0; objIndex < data.objectives.size(); ++objIndex) {
WeightVector direction(data.objectives.size(), storm::utility::zero<RationalNumberType>());
direction[objIndex] = storm::utility::one<RationalNumberType>();
refineResult(direction, weightVectorChecker, result);
}
bool converged=false;
do {
// Get the halfspace of the underApproximation with maximal distance to a vertex of the overApproximation
std::vector<storm::storage::geometry::Halfspace<RationalNumberType>> underApproxHalfspaces = result.underApproximation->getHalfspaces();
std::vector<Point> overApproxVertices = result.overApproximation->getVertices();
uint_fast64_t farestHalfspaceIndex = underApproxHalfspaces.size();
RationalNumberType farestDistance = storm::utility::zero<RationalNumberType>();
for(uint_fast64_t halfspaceIndex = 0; halfspaceIndex < underApproxHalfspaces.size(); ++halfspaceIndex) {
for(auto const& vertex : overApproxVertices) {
RationalNumberType distance = -underApproxHalfspaces[halfspaceIndex].euclideanDistance(vertex);
if(distance > farestDistance) {
farestHalfspaceIndex = halfspaceIndex;
farestDistance = distance;
}
}
}
STORM_LOG_DEBUG("Farest distance between under- and overApproximation is " << storm::utility::convertNumber<double>(farestDistance));
if(storm::utility::convertNumber<double>(farestDistance) < storm::settings::getModule<storm::settings::modules::MultiObjectiveSettings>().getPrecision()) {
result.precisionOfResult = farestDistance;
converged = true;
}
converged |= (storm::settings::getModule<storm::settings::modules::MultiObjectiveSettings>().isMaxIterationsSet() && result.iterations.size() >= storm::settings::getModule<storm::settings::modules::MultiObjectiveSettings>().getMaxIterations());
if(!converged) {
WeightVector direction = underApproxHalfspaces[farestHalfspaceIndex].normalVector();
refineResult(direction, weightVectorChecker, result);
}
} while(!converged);
}
template <class SparseModelType, typename RationalNumberType>
void SparseMultiObjectiveHelper<SparseModelType, RationalNumberType>::refineResult(WeightVector const& direction, SparseWeightedObjectivesModelCheckerHelper<SparseModelType> const& weightedObjectivesChecker, ReturnType& result) {
weightedObjectivesChecker.check(storm::utility::vector::convertNumericVector<typename SparseModelType::ValueType>(direction));
STORM_LOG_DEBUG("Result with weighted objectives is " << weightedObjectivesChecker.getInitialStateResultOfObjectives());
result.iterations.emplace_back(direction, storm::utility::vector::convertNumericVector<RationalNumberType>(weightedObjectivesChecker.getInitialStateResultOfObjectives()), weightedObjectivesChecker.getScheduler());
updateOverApproximation(result.iterations.back().point, result.iterations.back().weightVector);
updateUnderApproximation();
}
template <class SparseModelType, typename RationalNumberType>
typename SparseMultiObjectiveHelper<SparseModelType, RationalNumberType>::WeightVector SparseMultiObjectiveHelper<SparseModelType, RationalNumberType>::findSeparatingVector(Point const& pointToBeSeparated, std::shared_ptr<storm::storage::geometry::Polytope<RationalNumberType>> const& underApproximation, storm::storage::BitVector& individualObjectivesToBeChecked) const {
STORM_LOG_DEBUG("Searching a weight vector to seperate the point given by " << storm::utility::vector::convertNumericVector<double>(pointToBeSeparated) << ".");
if(underApproximation->isEmpty()) {
// In this case, every weight vector is separating
uint_fast64_t objIndex = individualObjectivesToBeChecked.getNextSetIndex(0) % pointToBeSeparated.size();
WeightVector result(pointToBeSeparated.size(), storm::utility::zero<RationalNumberType>());
result[objIndex] = storm::utility::one<RationalNumberType>();
individualObjectivesToBeChecked.set(objIndex, false);
return result;
}
// Reaching this point means that the underApproximation contains halfspaces.
// The seperating vector has to be the normal vector of one of these halfspaces.
// We pick one with maximal distance to the given point. However, weight vectors that correspond to checking individual objectives take precedence.
std::vector<storm::storage::geometry::Halfspace<RationalNumberType>> halfspaces = underApproximation->getHalfspaces();
uint_fast64_t farestHalfspaceIndex = 0;
// Note that we are looking for a halfspace that does not contain the point. Thus, the returned distances are negative.
RationalNumberType farestDistance = -halfspaces[farestHalfspaceIndex].euclideanDistance(pointToBeSeparated);
storm::storage::BitVector nonZeroVectorEntries = ~storm::utility::vector::filterZero<RationalNumberType>(halfspaces[farestHalfspaceIndex].normalVector());
bool foundSeparatingSingleObjectiveVector = nonZeroVectorEntries.getNumberOfSetBits()==1 && individualObjectivesToBeChecked.get(nonZeroVectorEntries.getNextSetIndex(0)) && farestDistance>storm::utility::zero<RationalNumberType>();
for(uint_fast64_t halfspaceIndex = 1; halfspaceIndex < halfspaces.size(); ++halfspaceIndex) {
RationalNumberType distance = -halfspaces[halfspaceIndex].euclideanDistance(pointToBeSeparated);
nonZeroVectorEntries = ~storm::utility::vector::filterZero<RationalNumberType>(halfspaces[farestHalfspaceIndex].normalVector());
bool isSeparatingSingleObjectiveVector = nonZeroVectorEntries.getNumberOfSetBits() == 1 && individualObjectivesToBeChecked.get(nonZeroVectorEntries.getNextSetIndex(0)) && distance>storm::utility::zero<RationalNumberType>();
// Check if this halfspace is 'better' than the current one
if((!foundSeparatingSingleObjectiveVector && isSeparatingSingleObjectiveVector ) || (foundSeparatingSingleObjectiveVector==isSeparatingSingleObjectiveVector && distance>farestDistance)) {
foundSeparatingSingleObjectiveVector |= isSeparatingSingleObjectiveVector;
farestHalfspaceIndex = halfspaceIndex;
farestDistance = distance;
}
}
if(foundSeparatingSingleObjectiveVector) {
nonZeroVectorEntries = ~storm::utility::vector::filterZero<RationalNumberType>(halfspaces[farestHalfspaceIndex].normalVector());
individualObjectivesToBeChecked.set(nonZeroVectorEntries.getNextSetIndex(0), false);
}
STORM_LOG_THROW(farestDistance>=storm::utility::zero<RationalNumberType>(), storm::exceptions::UnexpectedException, "There is no seperating vector.");
STORM_LOG_DEBUG("Found separating weight vector: " << storm::utility::vector::convertNumericVector<double>(halfspaces[farestHalfspaceIndex].normalVector()) << ".");
return halfspaces[farestHalfspaceIndex].normalVector();
}
template <class SparseModelType, typename RationalNumberType>
void SparseMultiObjectiveHelper<SparseModelType, RationalNumberType>::updateOverApproximation(std::vector<typename ReturnType::Iteration> const& iterations, std::shared_ptr<storm::storage::geometry::Polytope<RationalNumberType>>& overApproximation) {
storm::storage::geometry::Halfspace<RationalNumberType> h(iterations.back().weightVector, storm::utility::vector::dotProduct(iterations.back().weightVector, iterations.back().point));
// Due to numerical issues, it might be the case that the updated overapproximation does not contain the underapproximation,
// e.g., when the new point is strictly contained in the underapproximation. Check if this is the case.
RationalNumberType maximumOffset = h.offset();
for(auto const& iteration : iterations){
maximumOffset = std::max(maximumOffset, storm::utility::vector::dotProduct(h.normalVector(), iteration.point));
}
if(maximumOffset > h.offset()){
// We correct the issue by shifting the halfspace such that it contains the underapproximation
h.offset() = maximumOffset;
STORM_LOG_WARN("Numerical issues: The overapproximation would not contain the underapproximation. Hence, a halfspace is shifted by " << storm::utility::convertNumber<double>(h.euclideanDistance(iterations.back().point)) << ".");
}
overApproximation = overApproximation->intersection(h);
STORM_LOG_DEBUG("Updated OverApproximation to " << overApproximation->toString(true));
}
template <class SparseModelType, typename RationalNumberType>
void SparseMultiObjectiveHelper<SparseModelType, RationalNumberType>::updateUnderApproximation(std::vector<typename ReturnType::Iteration> const& iterations, std::shared_ptr<storm::storage::geometry::Polytope<RationalNumberType>>& underApproximation) {
std::vector<Point> paretoPointsVec;
paretoPointsVec.reserve(iterations.size());
for(auto const& iteration : iterations) {
paretoPointsVec.push_back(iterations.point);
}
STORM_LOG_WARN("REMOVE ADDING ADDITIONAL VERTICES AS SOON AS HYPRO WORKS FOR DEGENERATED POLYTOPES");
if(paretoPointsVec.front().size()==2) {
Point p1 = {-10000, -9999};
Point p2 = {-9999, -10000};
paretoPointsVec.push_back(p1);
paretoPointsVec.push_back(p2);
} else {
Point p1 = {-10000, -9999, -9999};
Point p2 = {-9999, -10000, -9999};
Point p3 = {-9999, -9999, -10000};
paretoPointsVec.push_back(p1);
paretoPointsVec.push_back(p2);
paretoPointsVec.push_back(p3);
}
boost::optional<Point> upperBounds;
if(!paretoPointsVec.empty()){
//Get the pointwise maximum of the pareto points
upperBounds = paretoPointsVec.front();
for(auto paretoPointIt = paretoPointsVec.begin()+1; paretoPointIt != paretoPointsVec.end(); ++paretoPointIt){
auto upperBoundIt = upperBounds->begin();
for(auto const& paretoPointCoordinate : *paretoPointIt){
if(paretoPointCoordinate>*upperBoundIt){
*upperBoundIt = paretoPointCoordinate;
}
++upperBoundIt;
}
}
}
underApproximation = storm::storage::geometry::Polytope<RationalNumberType>::create(paretoPointsVec)->downwardClosure(upperBounds);
STORM_LOG_DEBUG("Updated UnderApproximation to " << underApproximation->toString(true));
}
template <class SparseModelType, typename RationalNumberType>
bool SparseMultiObjectiveHelper<SparseModelType, RationalNumberType>::checkIfThresholdsAreSatisfied(std::shared_ptr<storm::storage::geometry::Polytope<RationalNumberType>> const& polytope, Point const& thresholds, storm::storage::BitVector const& strictThresholds) const {
std::vector<storm::storage::geometry::Halfspace<RationalNumberType>> halfspaces = polytope->getHalfspaces();
for(auto const& h : halfspaces) {
RationalNumberType distance = h.distance(thresholds);
if(distance < storm::utility::zero<RationalNumberType>()) {
return false;
}
if(distance == storm::utility::zero<RationalNumberType>()) {
// In this case, the thresholds point is on the boundary of the polytope.
// Check if this is problematic for the strict thresholds
for(auto strictThreshold : strictThresholds) {
if(h.normalVector()[strictThreshold] > storm::utility::zero<RationalNumberType>()) {
return false;
}
}
}
}
return true;
}
#ifdef STORM_HAVE_CARL
template class SparseMultiObjectiveHelper<storm::models::sparse::Mdp<double>, storm::RationalNumber>;
#endif
}
}
}

67
src/modelchecker/multiobjective/helper/SparseMultiObjectiveHelper.h
View File

@ -0,0 +1,67 @@
#ifndef STORM_MODELCHECKER_MULTIOBJECTIVE_HELPER_SPARSEMULTIOBJECTIVEHELPER_H_
#define STORM_MODELCHECKER_MULTIOBJECTIVE_HELPER_SPARSEMULTIOBJECTIVEHELPER_H_
#include "src/modelchecker/multiobjective/helper/SparseMultiObjectivePreprocessorReturnType.h"
#include "src/modelchecker/multiObjective/helper/SparseMultiObjectiveWeightVectorChecker.h"
#include "src/storage/geometry/Polytope.h"
#include "src/storage/TotalScheduler.h"
namespace storm {
namespace modelchecker {
namespace helper {
template <class SparseModelType, typename RationalNumberType>
class SparseMultiObjectiveHelper {
public:
typedef SpaseMultiObjectivePreprocessorReturnType<SparseModelType> PreprocessorData;
typedef SpaseMultiObjectiveHelperReturnType<RationalNumberType> ReturnType;
typedef std::vector<RationalNumberType> Point;
typedef std::vector<RationalNumberType> WeightVector;
static ReturnType check(PreprocessorData const& data);
private:
static void achievabilityQuery(PreprocessorData const& data, ReturnType& result);
static void numericalQuery(PreprocessorData const& data, ReturnType& result);
static void paretoQuery(PreprocessorData const& data, ReturnType& result);
/*
* Refines the current result w.r.t. the given direction vector
*/
static void refineResult(WeightVector const& direction, SparseWeightedObjectivesModelCheckerHelper<SparseModelType> const& weightedObjectivesChecker, ReturnType& result);
/*
* Returns a weight vector w that separates the under approximation from the given point p, i.e.,
* For each x in the under approximation, it holds that w*x <= w*p
*
* @param pointToBeSeparated the point that is to be seperated
* @param underapproximation the current under approximation
* @param objectivesToBeCheckedIndividually stores for each objective whether it still makes sense to check for this objective individually (i.e., with weight vector given by w_{objIndex} = 1 )
*/
static WeightVector findSeparatingVector(Point const& pointToBeSeparated, std::shared_ptr<storm::storage::geometry::Polytope<RationalNumberType>> const& underApproximation, storm::storage::BitVector& individualObjectivesToBeChecked) const;
/*
* Updates the over- and under approximation, respectively
*
* @param iterations the iterations performed so far. iterations.back() should be the newest iteration, whose information is not yet included in the approximation.
* @param over/underApproximation the current over/underApproximation that will be updated
*/
static void updateOverApproximation(std::vector<typename ReturnType::Iteration> const& iterations, std::shared_ptr<storm::storage::geometry::Polytope<RationalNumberType>>& overApproximation);
static void updateUnderApproximation(std::vector<typename ReturnType::Iteration> const& iterations, std::shared_ptr<storm::storage::geometry::Polytope<RationalNumberType>>& underApproximation);
/*
* Returns true iff there is a point in the given polytope that satisfies the given thresholds.
* It is assumed that the given polytope contains the downward closure of its vertices.
*/
static bool checkIfThresholdsAreSatisfied(std::shared_ptr<storm::storage::geometry::Polytope<RationalNumberType>> const& polytope, Point const& thresholds, storm::storage::BitVector const& strictThresholds) const;
};
}
}
}
#endif /* STORM_MODELCHECKER_MULTIOBJECTIVE_HELPER_SPARSEMULTIOBJECTIVEHELPER_H_ */

56
src/modelchecker/multiobjective/helper/SparseMultiObjectiveHelperReturnType.h
View File

@ -0,0 +1,56 @@
#ifndef STORM_MODELCHECKER_MULTIOBJECTIVE_HELPER_SPARSEMULTIOBJECTIVEHELPERRETURNTYPE_H_
#define STORM_MODELCHECKER_MULTIOBJECTIVE_HELPER_SPARSEMULTIOBJECTIVEHELPERRETURNTYPE_H_
#include <vector>
#include <memory>
#include <iomanip>
#include <boost/optional.hpp>
#include "src/storage/geometry/Polytope.h"
#include "src/storage/TotalScheduler.h"
namespace storm {
namespace modelchecker {
namespace helper {
template <typename RationalNumberType>
struct SparseMultiObjectiveHelperReturnType {
struct Iteration {
Iteration (std::vector<RationalNumberType> const& weightVector, std::vector<RationalNumberType> const& point, storm::storage::TotalScheduler const& scheduler) : weightVector(weightVector), point(point), scheduler(scheduler) {
//Intentionally left empty
}
Iteration (std::vector<RationalNumberType>&& weightVector, std::vector<RationalNumberType>&& point, storm::storage::TotalScheduler&& scheduler) : weightVector(weightVector), point(point), scheduler(scheduler) {
//Intentionally left empty
}
std::vector<RationalNumberType> weightVector;
std::vector<RationalNumberType> point;
storm::storage::TotalScheduler scheduler;
}
//Stores the results for the individual iterations
std::vector<Iteration> iterations;
//Stores an over/ under approximation of the set of achievable values
std::shared_ptr<storm::storage::geometry::Polytope<RationalNumberType>> overApproximation;
std::shared_ptr<storm::storage::geometry::Polytope<RationalNumberType>> underApproximation;
// Stores the result of an achievability query (if applicable)
// For a numerical query, stores whether there is one feasible solution
boost::optional<bool> thresholdsAreAchievable;
//Stores the result of a numerical query (if applicable)
boost::optional<RationalNumberType> numericalResult;
boost::optional<bool> optimumIsAchievable; //True if there is an actual scheduler that induces the computed supremum (i.e., supremum == maximum)
//Stores the achieved precision for numerical and pareto queries
boost::optional<RationalNumberType> precisionOfResult;
};
}
}
}
#endif /* STORM_MODELCHECKER_MULTIOBJECTIVE_HELPER_SPARSEMULTIOBJECTIVEHELPERRETURNTYPE_H_ */

490
src/modelchecker/multiobjective/helper/SparseMultiObjectiveModelCheckerHelper.cpp
View File

@ -1,490 +0,0 @@
#include "src/modelchecker/multiobjective/helper/SparseMultiObjectiveModelCheckerHelper.h"
#include "src/adapters/CarlAdapter.h"
#include "src/models/sparse/Mdp.h"
#include "src/models/sparse/StandardRewardModel.h"
#include "src/utility/constants.h"
#include "src/utility/vector.h"
#include "src/settings//SettingsManager.h"
#include "src/settings/modules/MultiObjectiveSettings.h"
#include "src/exceptions/UnexpectedException.h"
namespace storm {
namespace modelchecker {
namespace helper {
template <class SparseModelType, typename RationalNumberType>
SparseMultiObjectiveModelCheckerHelper<SparseModelType, RationalNumberType>::SparseMultiObjectiveModelCheckerHelper(Information& info) : info(info), weightedObjectivesChecker(info), numIterations(0) {
overApproximation = storm::storage::geometry::Polytope<RationalNumberType>::createUniversalPolytope();
underApproximation = storm::storage::geometry::Polytope<RationalNumberType>::createEmptyPolytope();
}
template <class SparseModelType, typename RationalNumberType>
SparseMultiObjectiveModelCheckerHelper<SparseModelType, RationalNumberType>::~SparseMultiObjectiveModelCheckerHelper() {
// Intentionally left empty
}
template <class SparseModelType, typename RationalNumberType>
void SparseMultiObjectiveModelCheckerHelper<SparseModelType, RationalNumberType>::check(Information & info) {
SparseMultiObjectiveModelCheckerHelper<SparseModelType, RationalNumberType> helper(info);
uint_fast64_t numOfObjectivesWithoutThreshold = 0;
for(auto const& obj : info.objectives) {
if(!obj.threshold) {
++numOfObjectivesWithoutThreshold;
}
}
if(numOfObjectivesWithoutThreshold == 0) {
helper.achievabilityQuery();
} else if (numOfObjectivesWithoutThreshold == 1) {
helper.numericalQuery();
} else if (numOfObjectivesWithoutThreshold == info.objectives.size()) {
helper.paretoQuery();
} else {
STORM_LOG_THROW(false, storm::exceptions::UnexpectedException, "The number of objecties without threshold is not valid. It should be either 0 (achievabilityQuery), 1 (numericalQuery), or " << info.objectives.size() << " (paretoQuery). Got " << numOfObjectivesWithoutThreshold << " instead.");
}
}
template <class SparseModelType, typename RationalNumberType>
void SparseMultiObjectiveModelCheckerHelper<SparseModelType, RationalNumberType>::achievabilityQuery() {
Point thresholds;
thresholds.reserve(info.objectives.size());
storm::storage::BitVector strictThresholds(info.objectives.size());
for(uint_fast64_t objIndex = 0; objIndex < info.objectives.size(); ++objIndex) {
STORM_LOG_THROW(info.objectives[objIndex].threshold, storm::exceptions::UnexpectedException, "Can not perform achievabilityQuery: No threshold given for at least one objective.");
thresholds.push_back(storm::utility::convertNumber<RationalNumberType>(*info.objectives[objIndex].threshold));
strictThresholds.set(objIndex, info.objectives[objIndex].thresholdIsStrict);
}
storm::storage::BitVector individualObjectivesToBeChecked(info.objectives.size(), true);
bool converged=false;
do {
WeightVector separatingVector = findSeparatingVector(thresholds, individualObjectivesToBeChecked);
refineSolution(separatingVector);
++numIterations;
// Check for convergence
if(!checkIfThresholdsAreSatisfied(overApproximation, thresholds, strictThresholds)){
std::cout << "PROPERTY VIOLATED" << std::endl;
converged=true;
}
if(checkIfThresholdsAreSatisfied(underApproximation, thresholds, strictThresholds)){
std::cout << "PROPERTY SATISFIED" << std::endl;
converged=true;
}
converged |= (storm::settings::getModule<storm::settings::modules::MultiObjectiveSettings>().isMaxIterationsSet() && numIterations >= storm::settings::getModule<storm::settings::modules::MultiObjectiveSettings>().getMaxIterations());
} while(!converged);
}
template <class SparseModelType, typename RationalNumberType>
void SparseMultiObjectiveModelCheckerHelper<SparseModelType, RationalNumberType>::numericalQuery() {
Point thresholds;
thresholds.reserve(info.objectives.size());
storm::storage::BitVector strictThresholds(info.objectives.size());
std::vector<storm::storage::geometry::Halfspace<RationalNumberType>> thresholdConstraints;
thresholdConstraints.reserve(info.objectives.size()-1);
uint_fast64_t optimizingObjIndex = info.objectives.size();
for(uint_fast64_t objIndex = 0; objIndex < info.objectives.size(); ++objIndex) {
if(info.objectives[objIndex].threshold) {
thresholds.push_back(storm::utility::convertNumber<RationalNumberType>(*info.objectives[objIndex].threshold));
WeightVector normalVector(info.objectives.size(), storm::utility::zero<RationalNumberType>());
normalVector[objIndex] = -storm::utility::one<RationalNumberType>();
thresholdConstraints.emplace_back(std::move(normalVector), -thresholds.back());
std::cout << "adding threshold constraint" << thresholdConstraints.back().toString(true) << std::endl;
strictThresholds.set(objIndex, info.objectives[objIndex].thresholdIsStrict);
} else {
STORM_LOG_ASSERT(optimizingObjIndex == info.objectives.size(), "Numerical Query invoked but there are multiple objectives without threshold");
optimizingObjIndex = objIndex;
thresholds.push_back(storm::utility::zero<RationalNumberType>());
}
}
STORM_LOG_ASSERT(optimizingObjIndex < info.objectives.size(), "Numerical Query invoked but there are no objectives without threshold");
auto thresholdsAsPolytope = storm::storage::geometry::Polytope<RationalNumberType>::create(thresholdConstraints);
WeightVector directionOfOptimizingObjective(info.objectives.size(), storm::utility::zero<RationalNumberType>());
directionOfOptimizingObjective[optimizingObjIndex] = storm::utility::one<RationalNumberType>();
// Try to find one valid solution
storm::storage::BitVector individualObjectivesToBeChecked(info.objectives.size(), true);
individualObjectivesToBeChecked.set(optimizingObjIndex, false);
do { //TODO introduce convergence criterion? (like max num of iterations)
WeightVector separatingVector = findSeparatingVector(thresholds, individualObjectivesToBeChecked);
refineSolution(separatingVector);
//Pick the threshold for the optimizing objective low enough so valid solutions are not excluded
for(auto const& paretoPoint : paretoOptimalPoints) {
thresholds[optimizingObjIndex] = std::min(thresholds[optimizingObjIndex], paretoPoint.first[optimizingObjIndex]);
}
if(!checkIfThresholdsAreSatisfied(overApproximation, thresholds, strictThresholds)){
std::cout << "INFEASIBLE OBJECTIVE THRESHOLDS" << std::endl;
return;
}
} while(!checkIfThresholdsAreSatisfied(underApproximation, thresholds, strictThresholds));
STORM_LOG_DEBUG("Found a solution that satisfies the objective thresholds.");
individualObjectivesToBeChecked.clear();
// Improve the found solution.
// Note that we do not have to care whether a threshold is strict anymore, because the resulting optimum should be
// the supremum over all strategies. Hence, one could combine a scheduler inducing the optimum value (but possibly violating strict
// thresholds) and (with very low probability) a scheduler that satisfies all (possibly strict) thresholds.
while (true) {
std::pair<Point, bool> optimizationRes = underApproximation->intersection(thresholdsAsPolytope)->optimize(directionOfOptimizingObjective);
STORM_LOG_THROW(optimizationRes.second, storm::exceptions::UnexpectedException, "The underapproximation is either unbounded or empty.");
thresholds[optimizingObjIndex] = optimizationRes.first[optimizingObjIndex];
STORM_LOG_DEBUG("Best solution found so far is " << storm::utility::convertNumber<double>(thresholds[optimizingObjIndex]) << ".");
//Compute an upper bound for the optimum and check for convergence
optimizationRes = overApproximation->intersection(thresholdsAsPolytope)->optimize(directionOfOptimizingObjective);
if(optimizationRes.second) {
RationalNumberType upperBoundOnOptimum = optimizationRes.first[optimizingObjIndex];
STORM_LOG_DEBUG("Solution can be improved by at most " << storm::utility::convertNumber<double>(upperBoundOnOptimum - thresholds[optimizingObjIndex]));
if(storm::utility::convertNumber<double>(upperBoundOnOptimum - thresholds[optimizingObjIndex]) < storm::settings::getModule<storm::settings::modules::MultiObjectiveSettings>().getPrecision()) { std::cout << "Found Optimum: " << storm::utility::convertNumber<double>(thresholds[optimizingObjIndex]) << "." << std::endl;
if(!checkIfThresholdsAreSatisfied(underApproximation, thresholds, strictThresholds)) {
std::cout << "Optimum is only the supremum" << std::endl;
}
return;
}
}
WeightVector separatingVector = findSeparatingVector(thresholds, individualObjectivesToBeChecked);
refineSolution(separatingVector);
}
}
template <class SparseModelType, typename RationalNumberType>
void SparseMultiObjectiveModelCheckerHelper<SparseModelType, RationalNumberType>::paretoQuery() {
//First consider the objectives individually
for(uint_fast64_t objIndex = 0; objIndex < info.objectives.size(); ++objIndex) {
WeightVector direction(info.objectives.size(), storm::utility::zero<RationalNumberType>());
direction[objIndex] = storm::utility::one<RationalNumberType>();
refineSolution(direction);
}
while(true) { //todo maxNumOfIterations ?
// Get the halfspace of the underApproximation with maximal distance to a vertex of the overApproximation
std::vector<storm::storage::geometry::Halfspace<RationalNumberType>> underApproxHalfspaces = underApproximation->getHalfspaces();
std::vector<Point> overApproxVertices = overApproximation->getVertices();
uint_fast64_t farestHalfspaceIndex = underApproxHalfspaces.size();
RationalNumberType farestDistance = storm::utility::zero<RationalNumberType>();
for(uint_fast64_t halfspaceIndex = 0; halfspaceIndex < underApproxHalfspaces.size(); ++halfspaceIndex) {
for(auto const& vertex : overApproxVertices) {
RationalNumberType distance = underApproxHalfspaces[halfspaceIndex].euclideanDistance(vertex);
// Note that the distances are negative, i.e., we have to consider the minimum distance
if(distance < farestDistance) {
farestHalfspaceIndex = halfspaceIndex;
farestDistance = distance;
}
}
}
STORM_LOG_DEBUG("Farest distance between under- and overApproximation is " << storm::utility::convertNumber<double>(farestDistance));
if(-storm::utility::convertNumber<double>(farestDistance) < storm::settings::getModule<storm::settings::modules::MultiObjectiveSettings>().getPrecision()) {
std::cout << "DONE computing the pareto curve" << std::endl;
return;
}
WeightVector direction = underApproxHalfspaces[farestHalfspaceIndex].normalVector();
// normalize the direction vector so the entries sum up to one
storm::utility::vector::scaleVectorInPlace(direction, storm::utility::one<RationalNumberType>() / std::accumulate(direction.begin(), direction.end(), storm::utility::zero<RationalNumberType>()));
refineSolution(direction);
}
}
template <class SparseModelType, typename RationalNumberType>
void SparseMultiObjectiveModelCheckerHelper<SparseModelType, RationalNumberType>::refineSolution(WeightVector const& direction) {
weightedObjectivesChecker.check(storm::utility::vector::convertNumericVector<ModelValueType>(direction));
storm::storage::TotalScheduler scheduler = weightedObjectivesChecker.getScheduler();
Point weightedObjectivesResult = storm::utility::vector::convertNumericVector<RationalNumberType>(weightedObjectivesChecker.getInitialStateResultOfObjectives());
STORM_LOG_DEBUG("Result with weighted objectives is " << storm::utility::vector::toString(weightedObjectivesResult));
// Insert the computed scheduler and check whether we have already seen it before
auto schedulerInsertRes = schedulers.insert(std::make_pair(std::move(scheduler), paretoOptimalPoints.end()));
if(schedulerInsertRes.second){
// The scheduler is new, so insert the newly computed pareto optimal point.
// For the corresponding scheduler, we safe the entry in the map.
schedulerInsertRes.first->second = paretoOptimalPoints.insert(std::make_pair(std::move(weightedObjectivesResult), std::vector<WeightVector>())).first;
}
// In the case where the scheduler is not new, we assume that the corresponding pareto optimal points for the old and new scheduler are equal.
// hence, no additional point will be inserted.
// Note that the values might not be exactly equal due to numerical issues.
Point const& paretoPoint = schedulerInsertRes.first->second->first;
std::vector<WeightVector>& weightVectors = schedulerInsertRes.first->second->second;
weightVectors.push_back(direction);
updateOverApproximation(paretoPoint, weightVectors.back());
updateUnderApproximation();
}
template <class SparseModelType, typename RationalNumberType>
typename SparseMultiObjectiveModelCheckerHelper<SparseModelType, RationalNumberType>::WeightVector SparseMultiObjectiveModelCheckerHelper<SparseModelType, RationalNumberType>::findSeparatingVector(Point const& pointToBeSeparated, storm::storage::BitVector& individualObjectivesToBeChecked) const {
STORM_LOG_DEBUG("Searching a weight vector to seperate the point given by " << storm::utility::vector::convertNumericVector<double>(pointToBeSeparated) << ".");
// First, we check 'simple' weight vectors that correspond to checking a single objective
while (!individualObjectivesToBeChecked.empty()){
uint_fast64_t objectiveIndex = individualObjectivesToBeChecked.getNextSetIndex(0);
individualObjectivesToBeChecked.set(objectiveIndex, false);
WeightVector normalVector(info.objectives.size(), storm::utility::zero<RationalNumberType>());
normalVector[objectiveIndex] = storm::utility::one<RationalNumberType>();
storm::storage::geometry::Halfspace<RationalNumberType> h(normalVector, storm::utility::vector::dotProduct(normalVector, pointToBeSeparated));
bool hIsSeparating = true;
for(auto const& paretoPoint : paretoOptimalPoints){
hIsSeparating &= h.contains(paretoPoint.first);
}
if(hIsSeparating) {
return h.normalVector();
}
}
if(underApproximation->isEmpty()) {
// In this case [1 0..0] is always separating
WeightVector result(info.objectives.size(), storm::utility::zero<RationalNumberType>());
result.front() = storm::utility::one<RationalNumberType>();
return result;
}
// Get the halfspace of the underApproximation with maximal distance to the given point
// Note that we are looking for a halfspace that does not contain the point. Thus the distances are negative and
// we are actually searching for the one with minimal distance.
std::vector<storm::storage::geometry::Halfspace<RationalNumberType>> halfspaces = underApproximation->getHalfspaces();
uint_fast64_t farestHalfspaceIndex = 0;
RationalNumberType farestDistance = halfspaces[farestHalfspaceIndex].euclideanDistance(pointToBeSeparated);
for(uint_fast64_t halfspaceIndex = 1; halfspaceIndex < halfspaces.size(); ++halfspaceIndex) {
RationalNumberType distance = halfspaces[halfspaceIndex].euclideanDistance(pointToBeSeparated);
if(distance < farestDistance) {
farestHalfspaceIndex = halfspaceIndex;
farestDistance = distance;
}
}
STORM_LOG_THROW(farestDistance<=storm::utility::zero<RationalNumberType>(), storm::exceptions::UnexpectedException, "There is no seperating halfspace.");
// normalize the result so the entries sum up to one
WeightVector result = halfspaces[farestHalfspaceIndex].normalVector();
storm::utility::vector::scaleVectorInPlace(result, storm::utility::one<RationalNumberType>() / std::accumulate(result.begin(), result.end(), storm::utility::zero<RationalNumberType>()));
return result;
/* old version using LP solving
// We now use LP solving to find a seperating halfspace.
// Note that StoRM's LPSover does not support rational numbers. Hence, we optimize a polytope instead.
// The variables of the LP are the normalVector entries w_0 ... w_{n-1} and the minimal distance d between the halfspace and the paretoOptimalPoints
uint_fast64_t numVariables = pointToBeSeparated.size() + 1;
std::vector<storm::storage::geometry::Halfspace<RationalNumberType>> constraints;
constraints.reserve((numVariables-1) + 2 + paretoOptimalPoints.size());
// w_i >= 0 and d>= 0 <==> -w_i <= 0 and -d <= 0
for(uint_fast64_t i = 0; i<numVariables; ++i) {
WeightVector constraintCoefficients(numVariables, storm::utility::zero<RationalNumberType>());
constraintCoefficients[i] = -storm::utility::one<RationalNumberType>();
RationalNumberType constraintOffset = storm::utility::zero<RationalNumberType>();
constraints.emplace_back(std::move(constraintCoefficients), std::move(constraintOffset));
}
// sum_i w_i == 1 <==> sum_i w_i <=1 and sum_i -w_i <= -1
{
WeightVector constraintCoefficients(numVariables, storm::utility::one<RationalNumberType>());
constraintCoefficients.back() = storm::utility::zero<RationalNumberType>();
RationalNumberType constraintOffset = storm::utility::one<RationalNumberType>();
constraints.emplace_back(constraintCoefficients, constraintOffset);
storm::utility::vector::scaleVectorInPlace(constraintCoefficients, -storm::utility::one<RationalNumberType>());
constraintOffset = -storm::utility::one<RationalNumberType>();
constraints.emplace_back(std::move(constraintCoefficients), std::move(constraintOffset));
}
// let q=pointToBeSeparated. For each paretoPoint x: q*w - x*w >= d <==> (x-q) * w + d <= 0
for(auto const& paretoPoint : paretoOptimalPoints){
WeightVector constraintCoefficients(numVariables-1);
storm::utility::vector::subtractVectors(paretoPoint.first, pointToBeSeparated, constraintCoefficients);
constraintCoefficients.push_back(storm::utility::one<RationalNumberType>());
RationalNumberType constraintOffset = storm::utility::zero<RationalNumberType>();
constraints.emplace_back(std::move(constraintCoefficients), std::move(constraintOffset));
}
WeightVector optimizationVector(numVariables);
optimizationVector.back() = storm::utility::one<RationalNumberType>(); // maximize d
std::pair<WeightVector, bool> optimizeResult = storm::storage::geometry::Polytope<RationalNumberType>::create(constraints)->optimize(optimizationVector);
STORM_LOG_THROW(optimizeResult.second, storm::exceptions::UnexpectedException, "There seems to be no seperating halfspace.");
optimizeResult.first.pop_back(); //drop the distance
RationalNumberType offset = storm::utility::vector::dotProduct(optimizeResult.first, pointToBeSeparated);
return storm::storage::geometry::Halfspace<RationalNumberType>(std::move(optimizeResult.first), std::move(offset));
*/
}
template <class SparseModelType, typename RationalNumberType>
void SparseMultiObjectiveModelCheckerHelper<SparseModelType, RationalNumberType>::updateOverApproximation(Point const& newPoint, WeightVector const& newWeightVector) {
storm::storage::geometry::Halfspace<RationalNumberType> h(newWeightVector, storm::utility::vector::dotProduct(newWeightVector, newPoint));
// Due to numerical issues, it might be the case that the updated overapproximation does not contain the underapproximation,
// e.g., when the new point is strictly contained in the underapproximation. Check if this is the case.
RationalNumberType maximumOffset = h.offset();
for(auto const& paretoPoint : paretoOptimalPoints){
maximumOffset = std::max(maximumOffset, storm::utility::vector::dotProduct(h.normalVector(), paretoPoint.first));
}
if(maximumOffset > h.offset()){
// We correct the issue by shifting the halfspace such that it contains the underapproximation
h.offset() = maximumOffset;
STORM_LOG_WARN("Numerical issues: The overapproximation would not contain the underapproximation. Hence, a halfspace is shifted by " << storm::utility::convertNumber<double>(h.euclideanDistance(newPoint)) << ".");
}
overApproximation = overApproximation->intersection(h);
STORM_LOG_DEBUG("Updated OverApproximation to " << overApproximation->toString(true));
}
template <class SparseModelType, typename RationalNumberType>
void SparseMultiObjectiveModelCheckerHelper<SparseModelType, RationalNumberType>::updateUnderApproximation() {
std::vector<Point> paretoPointsVec;
paretoPointsVec.reserve(paretoOptimalPoints.size());
for(auto const& paretoPoint : paretoOptimalPoints) {
paretoPointsVec.push_back(paretoPoint.first);
}
STORM_LOG_WARN("REMOVE ADDING ADDITIONAL VERTICES AS SOON AS HYPRO WORKS FOR DEGENERATED POLYTOPES");
Point p1 = {-1000, -999};
Point p2 = {-999, -1000};
paretoPointsVec.push_back(p1);
paretoPointsVec.push_back(p2);
boost::optional<Point> upperBounds;
if(!paretoPointsVec.empty()){
//Get the pointwise maximum of the pareto points
upperBounds = paretoPointsVec.front();
for(auto paretoPointIt = paretoPointsVec.begin()+1; paretoPointIt != paretoPointsVec.end(); ++paretoPointIt){
auto upperBoundIt = upperBounds->begin();
for(auto const& paretoPointCoordinate : *paretoPointIt){
if(paretoPointCoordinate>*upperBoundIt){
*upperBoundIt = paretoPointCoordinate;
}
++upperBoundIt;
}
}
}
underApproximation = storm::storage::geometry::Polytope<RationalNumberType>::create(paretoPointsVec)->downwardClosure(upperBounds);
STORM_LOG_DEBUG("Updated UnderApproximation to " << underApproximation->toString(true));
}
template <class SparseModelType, typename RationalNumberType>
bool SparseMultiObjectiveModelCheckerHelper<SparseModelType, RationalNumberType>::checkIfThresholdsAreSatisfied(std::shared_ptr<storm::storage::geometry::Polytope<RationalNumberType>> const& polytope, Point const& thresholds, storm::storage::BitVector const& strictThresholds) const {
std::vector<storm::storage::geometry::Halfspace<RationalNumberType>> halfspaces = polytope->getHalfspaces();
for(auto const& h : halfspaces) {
RationalNumberType distance = h.distance(thresholds);
if(distance < storm::utility::zero<RationalNumberType>()) {
return false;
}
if(distance == storm::utility::zero<RationalNumberType>()) {
// In this case, the thresholds point is on the boundary of the polytope.
// Check if this is problematic for the strict thresholds
for(auto strictThreshold : strictThresholds) {
if(h.normalVector()[strictThreshold] > storm::utility::zero<RationalNumberType>()) {
return false;
}
}
}
}
return true;
}
/*
bool existsObjectiveWithoutThreshold = false;
bool existsObjectiveWithStrictThreshold = false;
Point thresholdPoint;
thresholdPoint.reserve(info.objectives.size());
for(auto const& obj : info.objectives) {
if(obj.threshold) {
thresholdPoint.push_back(storm::utility::convertNumber<RationalNumber>(*obj.threshold));
if(obj.thresholdIsStrict) {
existsObjectiveWithStrictThreshold = true;
}
} else {
existsObjectiveWithoutThreshold = true;
thresholdPoint.push_back(storm::utility::one<RationalNumberType>());
}
}
if(existsObjectiveWithoutThreshold) {
// We need to find values for the objectives without thresholds in a way that the result is not affected by them.
// We build a polytope that contains valid solutions according to the thresholds and then analyze the intersection with the given polytope.
std::vector<storm::storage::geometry::Halfspace<RationalNumberType>> thresholdConstraints;
WeightVector directionForOptimization(info.objectives.size(), storm::utility::zero<RationalNumberType>());
for(uint_fast64_t objIndex = 0; objIndex < info.objectives.size(); ++objIndex) {
if(info.objectives[objIndex].threshold) {
WeightVector normalVector(info.objectives.size(), storm::utility::zero<RationalNumberType>());
normalVector[objIndex] = -storm::utility::one<RationalNumberType>();
thresholdConstraints.emplace_back(std::move(normalVector), -thresholdPoint[objIndex]);
directionForOptimization[objIndex] = -storm::utility::one<RationalNumberType>();
}
}
auto validSolutions = polytope->intersection(storm::storage::geometry::Polytope<RationalNumberType>::create(std::move(thresholdConstraints)));
if(validSolutions->isEmpty()) {
return false;
}
auto optimizationRes = validSolutions->optimize(directionForOptimization);
STORM_LOG_THROW(optimizationRes.second, storm::exceptions::UnexpectedException, "Couldn't check whether thresholds are satisfied for a given polytope: Didn't expect an infeasible or unbounded solution.");
for(uint_fast64_t objIndex = 0; objIndex < info.objectives.size(); ++objIndex) {
if(!info.objectives[objIndex].threshold) {
thresholdPoint[objIndex] = optimizationRes.first[objIndex] - storm::utility::one<RationalNumberType>();
}
}
}
// We then check whether the intersection of the thresholdPolytope and the given polytope contains any point (that also satisfies the strict thresholds);
std::cout << "Valid solutions are " << validSolutions->toString(true) << std::endl;
std::cout << "not empty" << std::endl;
if(existsObjectiveWithStrictThreshold) {
std::cout << "hasObjWithStrictThreshold" << std::endl;
std::cout << "optRes is " << storm::utility::vector::convertNumericVector<double>(optimizationRes.first) << std::endl;
for(uint_fast64_t objIndex = 0; objIndex < info.objectives.size(); ++objIndex) {
if(info.objectives[objIndex].threshold && info.objectives[objIndex].thresholdIsStrict) {
RationalNumberType threshold = storm::utility::convertNumber<RationalNumberType>( *info.objectives[objIndex].threshold);
if (optimizationRes.first[objIndex] <= threshold) {
return false;
}
}
}
}
return true;
Point thresholds;
thresholds.reserve(info.objectives.size());
for(auto const& obj : info.objectives) {
if(obj.threshold) {
thresholds.push_back(storm::utility::convertNumber<RationalNumber>(*obj.threshold));
} else {
STORM_LOG_ASSERT(numericalQueryResult, "Checking whether thresholds are satisfied for a numerical query but no threshold for the numerical objective was given.");
thresholds.push_back(*numericalQueryResult);
}
}
std::vector<storm::storage::geometry::Halfspace<RationalNumberType>> halfspaces = polytope->getHalfspaces();
for(auto const& h : halfspaces) {
RationalNumberType distance = h.distance(thresholds);
if(distance < storm::utility::zero<RationalNumberType>()) {
return false;
}
if(distance == storm::utility::zero<RationalNumberType>()) {
// In this case, the thresholds point is on the boundary of the polytope.
// Check if this is problematic for the strict thresholds
for(uint_fast64_t objIndex = 0; objIndex < info.objectives.size(); ++objIndex) {
if(info.objectives[objIndex].thresholdIsStrict && h.normalVector()[objIndex] > storm::utility::zero<RationalNumberType>()) {
return false;
}
}
}
}
return true;
}
*/
#ifdef STORM_HAVE_CARL
template class SparseMultiObjectiveModelCheckerHelper<storm::models::sparse::Mdp<double>, storm::RationalNumber>;
#endif
}
}
}

80
src/modelchecker/multiobjective/helper/SparseMultiObjectiveModelCheckerHelper.h
View File

@ -1,80 +0,0 @@
#ifndef STORM_MODELCHECKER_MULTIOBJECTIVE_HELPER_SPARSEMULTIOBJECTIVEMODELCHECKERHELPER_H_
#define STORM_MODELCHECKER_MULTIOBJECTIVE_HELPER_SPARSEMULTIOBJECTIVEMODELCHECKERHELPER_H_
#include <map>
#include <unordered_map>
#include "src/modelchecker/multiobjective/helper/SparseMultiObjectiveModelCheckerInformation.h"
#include "src/modelchecker/multiObjective/helper/SparseWeightedObjectivesModelCheckerHelper.h"
#include "src/storage/geometry/Polytope.h"
#include "src/storage/TotalScheduler.h"
namespace storm {
namespace modelchecker {
namespace helper {
template <class SparseModelType, typename RationalNumberType>
class SparseMultiObjectiveModelCheckerHelper {
public:
typedef typename SparseModelType::ValueType ModelValueType;
typedef typename SparseModelType::RewardModelType RewardModelType;
typedef SparseMultiObjectiveModelCheckerInformation<SparseModelType> Information;
typedef std::vector<RationalNumberType> Point;
typedef std::vector<RationalNumberType> WeightVector;
static void check(Information& info);
private:
SparseMultiObjectiveModelCheckerHelper(Information& info);
~SparseMultiObjectiveModelCheckerHelper();
void achievabilityQuery();
void numericalQuery();
void paretoQuery();
/*
* Refines the solution w.r.t. the given direction vector
*/
void refineSolution(WeightVector const& direction);
/*
* Returns a weight vector w that separates the underapproximation from the given point p, i.e.,
* For each x in the underApproximation, it holds that w*x <= w*p
*
* @param pointToBeSeparated the point that is to be seperated
* @param objectivesToBeCheckedIndividually stores for each objective whether it still makes sense to check for this objective individually (i.e., with weight vector given by w_{objIndex} = 1 )
*/
WeightVector findSeparatingVector(Point const& pointToBeSeparated, storm::storage::BitVector& individualObjectivesToBeChecked) const;
void updateOverApproximation(Point const& newPoint, WeightVector const& newWeightVector);
void updateUnderApproximation();
/*
* Returns true iff there is a point in the given polytope that satisfies the given thresholds.
* It is assumed that the given polytope contains the downward closure of its vertices.
*/
bool checkIfThresholdsAreSatisfied(std::shared_ptr<storm::storage::geometry::Polytope<RationalNumberType>> const& polytope, Point const& thresholds, storm::storage::BitVector const& strictThresholds) const;
Information& info;
SparseWeightedObjectivesModelCheckerHelper<SparseModelType> weightedObjectivesChecker;
//TODO: sort points as needed
std::map<Point, std::vector<WeightVector>> paretoOptimalPoints;
std::unordered_map<storm::storage::TotalScheduler, typename std::map<Point, std::vector<WeightVector>>::iterator> schedulers;
uint_fast64_t numIterations;
std::shared_ptr<storm::storage::geometry::Polytope<RationalNumberType>> overApproximation;
std::shared_ptr<storm::storage::geometry::Polytope<RationalNumberType>> underApproximation;
};
}
}
}
#endif /* STORM_MODELCHECKER_MULTIOBJECTIVE_HELPER_SPARSEMULTIOBJECTIVEMODELCHECKERHELPER_H_ */

154
src/modelchecker/multiobjective/helper/SparseMdpMultiObjectivePreprocessingHelper.cpp → src/modelchecker/multiobjective/helper/SparseMultiObjectivePreprocessor.cpp
View File

@ -1,15 +1,11 @@
#include "src/modelchecker/multiobjective/helper/SparseMdpMultiObjectivePreprocessingHelper.h"
#include "src/modelchecker/multiobjective/helper/SparseMultiObjectivePreprocessor.h"
#include "src/models/sparse/Mdp.h"
#include "src/models/sparse/StandardRewardModel.h"
#include "src/modelchecker/propositional/SparsePropositionalModelChecker.h"
#include "src/modelchecker/results/ExplicitQualitativeCheckResult.h"
#include "src/transformer/StateDuplicator.h"
#include "src/storage/BitVector.h"
#include "src/transformer/StateDuplicator.h"
#include "src/utility/macros.h"
#include "src/utility/vector.h"
#include "src/utility/graph.h"
@ -21,24 +17,22 @@ namespace storm {
namespace modelchecker {
namespace helper {
template<typename SparseMdpModelType>
SparseMultiObjectiveModelCheckerInformation<SparseMdpModelType> SparseMdpMultiObjectivePreprocessingHelper<SparseMdpModelType>::preprocess(storm::logic::MultiObjectiveFormula const& originalFormula, SparseMdpModelType const& originalModel) {
Information info(std::move(originalModel));
template<typename SparseModelType>
typename SparseMultiObjectivePreprocessor<SparseModelType>::ReturnType SparseMultiObjectivePreprocessor<SparseModelType>::preprocess(storm::logic::MultiObjectiveFormula const& originalFormula, SparseModelType const& originalModel) {
ReturnType data(std::move(originalModel));
//Initialize the state mapping.
info.newToOldStateIndexMapping = storm::utility::vector::buildVectorForRange(0, info.preprocessedModel.getNumberOfStates());
data.newToOldStateIndexMapping = storm::utility::vector::buildVectorForRange(0, data.preprocessedModel.getNumberOfStates());
//Gather information regarding the individual objectives
for(auto const& subFormula : originalFormula.getSubFormulas()){
addObjective(subFormula, info);
addObjective(subFormula, data);
}
// Find out whether negated rewards should be considered.
//setWhetherNegatedRewardsAreConsidered(info);
//Invoke preprocessing on objectives
for (auto& obj : info.objectives){
for (auto& obj : data.objectives){
STORM_LOG_DEBUG("Preprocessing objective " << *obj.originalFormula << ".");
if(obj.originalFormula->isProbabilityOperatorFormula()){
preprocessFormula(obj.originalFormula->asProbabilityOperatorFormula(), info, obj);
preprocessFormula(obj.originalFormula->asProbabilityOperatorFormula(), data, obj);
} else if(obj.originalFormula->isRewardOperatorFormula()){
preprocessFormula(obj.originalFormula->asRewardOperatorFormula(), info, obj);
preprocessFormula(obj.originalFormula->asRewardOperatorFormula(), data, obj);
} else {
STORM_LOG_THROW(false, storm::exceptions::InvalidPropertyException, "Could not preprocess the subformula " << *obj.originalFormula << " of " << originalFormula << " because it is not supported");
}
@ -46,27 +40,27 @@ namespace storm {
//We can now remove all original reward models to save some memory
std::set<std::string> origRewardModels = originalFormula.getReferencedRewardModels();
for (auto const& rewModel : origRewardModels){
info.preprocessedModel.removeRewardModel(rewModel);
data.preprocessedModel.removeRewardModel(rewModel);
}
return info;
return data;
}
template<typename SparseMdpModelType>
void SparseMdpMultiObjectivePreprocessingHelper<SparseMdpModelType>::addObjective(std::shared_ptr<storm::logic::Formula const> const& formula, Information& info) {
template<typename SparseModelType>
void SparseMultiObjectivePreprocessor<SparseModelType>::addObjective(std::shared_ptr<storm::logic::Formula const> const& formula, ReturnType& data) {
STORM_LOG_THROW(formula->isOperatorFormula(), storm::exceptions::InvalidPropertyException, "Expected an OperatorFormula as subformula of multi-objective query but got " << *formula);
typename Information::ObjectiveInformation objective;
ObjectiveInformation objective;
objective.originalFormula = formula;
objective.rewardModelName = "objective" + std::to_string(info.objectives.size());
objective.rewardModelName = "objective" + std::to_string(data.objectives.size());
// Make sure the new reward model gets a unique name
while(info.preprocessedModel.hasRewardModel(objective.rewardModelName)){
while(data.preprocessedModel.hasRewardModel(objective.rewardModelName)){
objective.rewardModelName = "_" + objective.rewardModelName;
}
info.objectives.push_back(std::move(objective));
data.objectives.push_back(std::move(objective));
}
template<typename SparseMdpModelType>
void SparseMdpMultiObjectivePreprocessingHelper<SparseMdpModelType>::setStepBoundOfObjective(typename Information::ObjectiveInformation& objective){
template<typename SparseModelType>
void SparseMultiObjectivePreprocessor<SparseModelType>::setStepBoundOfObjective(ObjectiveInformation& objective){
if(objective.originalFormula->isProbabilityOperatorFormula()){
storm::logic::Formula const& f = objective.originalFormula->asProbabilityOperatorFormula().getSubformula();
if(f.isBoundedUntilFormula()){
@ -84,24 +78,8 @@ namespace storm {
}
}
template<typename SparseMdpModelType>
void SparseMdpMultiObjectivePreprocessingHelper<SparseMdpModelType>::setWhetherNegatedRewardsAreConsidered(Information& info) {
// Negative Rewards are considered whenever there is an unbounded reward objective that requires to minimize.
// If there is no unbounded reward objective, we make a majority vote.
/* uint_fast64_t numOfMinimizingObjectives = 0;
for(auto const& obj : info.objectives){
if(obj.originalFormula->isRewardOperatorFormula() && !obj.stepBound){
info.negatedRewardsConsidered = obj.originalFormulaMinimizes;
return;
}
numOfMinimizingObjectives += obj.originalFormulaMinimizes ? 1 : 0;
}
// Reaching this point means that we make a majority vote
info.negatedRewardsConsidered = (numOfMinimizingObjectives*2) > info.objectives.size();
*/}
template<typename SparseMdpModelType>
void SparseMdpMultiObjectivePreprocessingHelper<SparseMdpModelType>::preprocessFormula(storm::logic::ProbabilityOperatorFormula const& formula, Information& info, typename Information::ObjectiveInformation& currentObjective) {
template<typename SparseModelType>
void SparseMultiObjectivePreprocessor<SparseModelType>::preprocessFormula(storm::logic::ProbabilityOperatorFormula const& formula, ReturnType& data, ObjectiveInformation& currentObjective) {
//Set the information for the objective
if(formula.hasBound()) {
currentObjective.threshold = formula.getBound().threshold;
@ -118,26 +96,26 @@ namespace storm {
// Invoke preprocessing for subformula
if(formula.getSubformula().isUntilFormula()){
preprocessFormula(formula.getSubformula().asUntilFormula(), info, currentObjective);
preprocessFormula(formula.getSubformula().asUntilFormula(), data, currentObjective);
} else if(formula.getSubformula().isBoundedUntilFormula()){
preprocessFormula(formula.getSubformula().asBoundedUntilFormula(), info, currentObjective);
preprocessFormula(formula.getSubformula().asBoundedUntilFormula(), data, currentObjective);
} else if(formula.getSubformula().isGloballyFormula()){
preprocessFormula(formula.getSubformula().asGloballyFormula(), info, currentObjective);
preprocessFormula(formula.getSubformula().asGloballyFormula(), data, currentObjective);
} else if(formula.getSubformula().isEventuallyFormula()){
preprocessFormula(formula.getSubformula().asEventuallyFormula(), info, currentObjective);
preprocessFormula(formula.getSubformula().asEventuallyFormula(), data, currentObjective);
} else {
STORM_LOG_THROW(false, storm::exceptions::InvalidPropertyException, "The subformula of " << formula << " is not supported.");
}
}
template<typename SparseMdpModelType>
void SparseMdpMultiObjectivePreprocessingHelper<SparseMdpModelType>::preprocessFormula(storm::logic::RewardOperatorFormula const& formula, Information& info, typename Information::ObjectiveInformation& currentObjective) {
template<typename SparseModelType>
void SparseMultiObjectivePreprocessor<SparseModelType>::preprocessFormula(storm::logic::RewardOperatorFormula const& formula, ReturnType& data, ObjectiveInformation& currentObjective) {
//TODO Make sure that our decision for negative rewards is consistent with the formula
// STORM_LOG_THROW(info.negatedRewardsConsidered == currentObjective.originalFormulaMinimizes, storm::exceptions::InvalidPropertyException, "Decided to consider " << (info.negatedRewardsConsidered ? "negated" : "non-negated") << "rewards, but the formula " << formula << (currentObjective.originalFormulaMinimizes ? " minimizes" : "maximizes") << ". Reward objectives should either all minimize or all maximize.");
// STORM_LOG_THROW(data.negatedRewardsConsidered == currentObjective.originalFormulaMinimizes, storm::exceptions::InvalidPropertyException, "Decided to consider " << (data.negatedRewardsConsidered ? "negated" : "non-negated") << "rewards, but the formula " << formula << (currentObjective.originalFormulaMinimizes ? " minimizes" : "maximizes") << ". Reward objectives should either all minimize or all maximize.");
// Check if the reward model is uniquely specified
STORM_LOG_THROW((formula.hasRewardModelName() && info.preprocessedModel.hasRewardModel(formula.getRewardModelName()))
|| info.preprocessedModel.hasUniqueRewardModel(), storm::exceptions::InvalidPropertyException, "The reward model is not unique and the formula " << formula << " does not specify a reward model.");
STORM_LOG_THROW((formula.hasRewardModelName() && data.preprocessedModel.hasRewardModel(formula.getRewardModelName()))
|| data.preprocessedModel.hasUniqueRewardModel(), storm::exceptions::InvalidPropertyException, "The reward model is not unique and the formula " << formula << " does not specify a reward model.");
//Set the information for the objective
if(formula.hasBound()) {
@ -155,59 +133,59 @@ namespace storm {
// Invoke preprocessing for subformula
if(formula.getSubformula().isEventuallyFormula()){
preprocessFormula(formula.getSubformula().asEventuallyFormula(), info, currentObjective, formula.getOptionalRewardModelName());
preprocessFormula(formula.getSubformula().asEventuallyFormula(), data, currentObjective, formula.getOptionalRewardModelName());
} else if(formula.getSubformula().isCumulativeRewardFormula()) {
preprocessFormula(formula.getSubformula().asCumulativeRewardFormula(), info, currentObjective, formula.getOptionalRewardModelName());
preprocessFormula(formula.getSubformula().asCumulativeRewardFormula(), data, currentObjective, formula.getOptionalRewardModelName());
} else {
STORM_LOG_THROW(false, storm::exceptions::InvalidPropertyException, "The subformula of " << formula << " is not supported.");
}
}
template<typename SparseMdpModelType>
void SparseMdpMultiObjectivePreprocessingHelper<SparseMdpModelType>::preprocessFormula(storm::logic::UntilFormula const& formula, Information& info, typename Information::ObjectiveInformation& currentObjective) {
template<typename SparseModelType>
void SparseMultiObjectivePreprocessor<SparseModelType>::preprocessFormula(storm::logic::UntilFormula const& formula, ReturnType& data, ObjectiveInformation& currentObjective) {
CheckTask<storm::logic::Formula> phiTask(formula.getLeftSubformula());
CheckTask<storm::logic::Formula> psiTask(formula.getRightSubformula());
storm::modelchecker::SparsePropositionalModelChecker<SparseMdpModelType> mc(info.preprocessedModel);
storm::modelchecker::SparsePropositionalModelChecker<SparseModelType> mc(data.preprocessedModel);
STORM_LOG_THROW(mc.canHandle(phiTask) && mc.canHandle(psiTask), storm::exceptions::InvalidPropertyException, "The subformulas of " << formula << " should be propositional.");
storm::storage::BitVector phiStates = mc.check(phiTask)->asExplicitQualitativeCheckResult().getTruthValuesVector();
storm::storage::BitVector psiStates = mc.check(psiTask)->asExplicitQualitativeCheckResult().getTruthValuesVector();
auto duplicatorResult = storm::transformer::StateDuplicator<SparseMdpModelType>::transform(info.preprocessedModel, ~phiStates | psiStates);
info.preprocessedModel = std::move(*duplicatorResult.model);
auto duplicatorResult = storm::transformer::StateDuplicator<SparseModelType>::transform(data.preprocessedModel, ~phiStates | psiStates);
data.preprocessedModel = std::move(*duplicatorResult.model);
// duplicatorResult.newToOldStateIndexMapping now reffers to the indices of the model we had before preprocessing this formula.
// This might not be the actual original model.
for(auto & originalModelStateIndex : duplicatorResult.newToOldStateIndexMapping){
originalModelStateIndex = info.newToOldStateIndexMapping[originalModelStateIndex];
originalModelStateIndex = data.newToOldStateIndexMapping[originalModelStateIndex];
}
info.newToOldStateIndexMapping = std::move(duplicatorResult.newToOldStateIndexMapping);
data.newToOldStateIndexMapping = std::move(duplicatorResult.newToOldStateIndexMapping);
// build stateAction reward vector that gives (one*transitionProbability) reward whenever a transition leads from the firstCopy to a psiState
storm::storage::BitVector newPsiStates(info.preprocessedModel.getNumberOfStates(), false);
storm::storage::BitVector newPsiStates(data.preprocessedModel.getNumberOfStates(), false);
for(auto const& oldPsiState : psiStates){
//note that psiStates are always located in the second copy
newPsiStates.set(duplicatorResult.secondCopyOldToNewStateIndexMapping[oldPsiState], true);
}
std::vector<ValueType> objectiveRewards(info.preprocessedModel.getTransitionMatrix().getRowCount(), storm::utility::zero<ValueType>());
std::vector<ValueType> objectiveRewards(data.preprocessedModel.getTransitionMatrix().getRowCount(), storm::utility::zero<ValueType>());
for (auto const& firstCopyState : duplicatorResult.firstCopy) {
for (uint_fast64_t row = info.preprocessedModel.getTransitionMatrix().getRowGroupIndices()[firstCopyState]; row < info.preprocessedModel.getTransitionMatrix().getRowGroupIndices()[firstCopyState + 1]; ++row) {
objectiveRewards[row] = info.preprocessedModel.getTransitionMatrix().getConstrainedRowSum(row, newPsiStates);
for (uint_fast64_t row = data.preprocessedModel.getTransitionMatrix().getRowGroupIndices()[firstCopyState]; row < data.preprocessedModel.getTransitionMatrix().getRowGroupIndices()[firstCopyState + 1]; ++row) {
objectiveRewards[row] = data.preprocessedModel.getTransitionMatrix().getConstrainedRowSum(row, newPsiStates);
}
}
if(currentObjective.isNegative){
storm::utility::vector::scaleVectorInPlace(objectiveRewards, -storm::utility::one<ValueType>());
}
info.preprocessedModel.addRewardModel(currentObjective.rewardModelName, RewardModelType(boost::none, objectiveRewards));
data.preprocessedModel.addRewardModel(currentObjective.rewardModelName, RewardModelType(boost::none, objectiveRewards));
}
template<typename SparseMdpModelType>
void SparseMdpMultiObjectivePreprocessingHelper<SparseMdpModelType>::preprocessFormula(storm::logic::BoundedUntilFormula const& formula, Information& info, typename Information::ObjectiveInformation& currentObjective) {
template<typename SparseModelType>
void SparseMultiObjectivePreprocessor<SparseModelType>::preprocessFormula(storm::logic::BoundedUntilFormula const& formula, ReturnType& data, ObjectiveInformation& currentObjective) {
STORM_LOG_THROW(formula.hasDiscreteTimeBound(), storm::exceptions::InvalidPropertyException, "Expected a boundedUntilFormula with a discrete time bound but got " << formula << ".");
currentObjective.stepBound = formula.getDiscreteTimeBound();
preprocessFormula(storm::logic::UntilFormula(formula.getLeftSubformula().asSharedPointer(), formula.getRightSubformula().asSharedPointer()), info, currentObjective);
preprocessFormula(storm::logic::UntilFormula(formula.getLeftSubformula().asSharedPointer(), formula.getRightSubformula().asSharedPointer()), data, currentObjective);
}
template<typename SparseMdpModelType>
void SparseMdpMultiObjectivePreprocessingHelper<SparseMdpModelType>::preprocessFormula(storm::logic::GloballyFormula const& formula, Information& info, typename Information::ObjectiveInformation& currentObjective) {
template<typename SparseModelType>
void SparseMultiObjectivePreprocessor<SparseModelType>::preprocessFormula(storm::logic::GloballyFormula const& formula, ReturnType& data, ObjectiveInformation& currentObjective) {
// The formula will be transformed to an until formula for the complementary event
currentObjective.isInverted = true;
if(currentObjective.threshold) {
@ -221,54 +199,54 @@ namespace storm {
}
currentObjective.isNegative = !currentObjective.isNegative;
auto negatedSubformula = std::make_shared<storm::logic::UnaryBooleanStateFormula>(storm::logic::UnaryBooleanStateFormula::OperatorType::Not, formula.getSubformula().asSharedPointer());
preprocessFormula(storm::logic::UntilFormula(storm::logic::Formula::getTrueFormula(), negatedSubformula), info, currentObjective);
preprocessFormula(storm::logic::UntilFormula(storm::logic::Formula::getTrueFormula(), negatedSubformula), data, currentObjective);
}
template<typename SparseMdpModelType>
void SparseMdpMultiObjectivePreprocessingHelper<SparseMdpModelType>::preprocessFormula(storm::logic::EventuallyFormula const& formula, Information& info, typename Information::ObjectiveInformation& currentObjective, boost::optional<std::string> const& optionalRewardModelName) {
template<typename SparseModelType>
void SparseMultiObjectivePreprocessor<SparseModelType>::preprocessFormula(storm::logic::EventuallyFormula const& formula, ReturnType& data, ObjectiveInformation& currentObjective, boost::optional<std::string> const& optionalRewardModelName) {
if(formula.isReachabilityProbabilityFormula()){
preprocessFormula(storm::logic::UntilFormula(storm::logic::Formula::getTrueFormula(), formula.getSubformula().asSharedPointer()), info, currentObjective);
preprocessFormula(storm::logic::UntilFormula(storm::logic::Formula::getTrueFormula(), formula.getSubformula().asSharedPointer()), data, currentObjective);
return;
}
STORM_LOG_THROW(formula.isReachabilityRewardFormula(), storm::exceptions::InvalidPropertyException, "The formula " << formula << " neither considers reachability Probabilities nor reachability rewards");
CheckTask<storm::logic::Formula> targetTask(formula.getSubformula());
storm::modelchecker::SparsePropositionalModelChecker<SparseMdpModelType> mc(info.preprocessedModel);
storm::modelchecker::SparsePropositionalModelChecker<SparseModelType> mc(data.preprocessedModel);
STORM_LOG_THROW(mc.canHandle(targetTask), storm::exceptions::InvalidPropertyException, "The subformula of " << formula << " should be propositional.");
storm::storage::BitVector targetStates = mc.check(targetTask)->asExplicitQualitativeCheckResult().getTruthValuesVector();
STORM_LOG_WARN("There is no check for reward finiteness yet"); //TODO
auto duplicatorResult = storm::transformer::StateDuplicator<SparseMdpModelType>::transform(info.preprocessedModel, targetStates);
info.preprocessedModel = std::move(*duplicatorResult.model);
auto duplicatorResult = storm::transformer::StateDuplicator<SparseModelType>::transform(data.preprocessedModel, targetStates);
data.preprocessedModel = std::move(*duplicatorResult.model);
// duplicatorResult.newToOldStateIndexMapping now reffers to the indices of the model we had before preprocessing this formula.
// This might not be the actual original model.
for(auto & originalModelStateIndex : duplicatorResult.newToOldStateIndexMapping){
originalModelStateIndex = info.newToOldStateIndexMapping[originalModelStateIndex];
originalModelStateIndex = data.newToOldStateIndexMapping[originalModelStateIndex];
}
info.newToOldStateIndexMapping = std::move(duplicatorResult.newToOldStateIndexMapping);
data.newToOldStateIndexMapping = std::move(duplicatorResult.newToOldStateIndexMapping);
// Add a reward model that gives zero reward to the states of the second copy.
std::vector<ValueType> objectiveRewards = info.preprocessedModel.getRewardModel(optionalRewardModelName ? optionalRewardModelName.get() : "").getTotalRewardVector(info.preprocessedModel.getTransitionMatrix());
std::vector<ValueType> objectiveRewards = data.preprocessedModel.getRewardModel(optionalRewardModelName ? optionalRewardModelName.get() : "").getTotalRewardVector(data.preprocessedModel.getTransitionMatrix());
storm::utility::vector::setVectorValues(objectiveRewards, duplicatorResult.secondCopy, storm::utility::zero<ValueType>());
if(currentObjective.isNegative){
storm::utility::vector::scaleVectorInPlace(objectiveRewards, -storm::utility::one<ValueType>());
}
info.preprocessedModel.addRewardModel(currentObjective.rewardModelName, RewardModelType(boost::none, objectiveRewards));
data.preprocessedModel.addRewardModel(currentObjective.rewardModelName, RewardModelType(boost::none, objectiveRewards));
}
template<typename SparseMdpModelType>
void SparseMdpMultiObjectivePreprocessingHelper<SparseMdpModelType>::preprocessFormula(storm::logic::CumulativeRewardFormula const& formula, Information& info, typename Information::ObjectiveInformation& currentObjective, boost::optional<std::string> const& optionalRewardModelName) {
template<typename SparseModelType>
void SparseMultiObjectivePreprocessor<SparseModelType>::preprocessFormula(storm::logic::CumulativeRewardFormula const& formula, ReturnType& data, ObjectiveInformation& currentObjective, boost::optional<std::string> const& optionalRewardModelName) {
STORM_LOG_THROW(formula.hasDiscreteTimeBound(), storm::exceptions::InvalidPropertyException, "Expected a cumulativeRewardFormula with a discrete time bound but got " << formula << ".");
currentObjective.stepBound = formula.getDiscreteTimeBound();
std::vector<ValueType> objectiveRewards = info.preprocessedModel.getRewardModel(optionalRewardModelName ? optionalRewardModelName.get() : "").getTotalRewardVector(info.preprocessedModel.getTransitionMatrix());
std::vector<ValueType> objectiveRewards = data.preprocessedModel.getRewardModel(optionalRewardModelName ? optionalRewardModelName.get() : "").getTotalRewardVector(data.preprocessedModel.getTransitionMatrix());
if(currentObjective.isNegative){
storm::utility::vector::scaleVectorInPlace(objectiveRewards, -storm::utility::one<ValueType>());
}
info.preprocessedModel.addRewardModel(currentObjective.rewardModelName, RewardModelType(boost::none, objectiveRewards));
data.preprocessedModel.addRewardModel(currentObjective.rewardModelName, RewardModelType(boost::none, objectiveRewards));
}
template class SparseMdpMultiObjectivePreprocessingHelper<storm::models::sparse::Mdp<double>>;
template class SparseMultiObjectivePreprocessor<storm::models::sparse::Mdp<double>>;
}

67
src/modelchecker/multiobjective/helper/SparseMultiObjectivePreprocessor.h
View File

@ -0,0 +1,67 @@
#ifndef STORM_MODELCHECKER_MULTIOBJECTIVE_HELPER_SPARSEMULTIOBJECTIVEPREPROCESSOR_H_
#define STORM_MODELCHECKER_MULTIOBJECTIVE_HELPER_SPARSEMULTIOBJECTIVEPREPROCESSOR_H_
#include "src/logic/Formulas.h"
#include "src/modelchecker/multiobjective/helper/SparseMultiObjectivePreprocessorReturnType.h"
namespace storm {
namespace modelchecker {
namespace helper {
/*
* Helper Class to invoke the necessary preprocessing for multi objective model checking
*/
template <class SparseModelType>
class SparseMultiObjectivePreprocessor {
public:
typedef typename SparseModelType::ValueType ValueType;
typedef typename SparseModelType::RewardModelType RewardModelType;
typedef SparseMultiObjectivePreprocessorReturnType<SparseModelType> ReturnType;
typedef typename SparseMultiObjectivePreprocessorReturnType<SparseModelType>::ObjectiveInformation ObjectiveInformation;
/*!
* Preprocesses the given model w.r.t. the given formulas.
* @param originalModel The considered model
* @param originalFormula the considered formula. The subformulas should only contain one OperatorFormula at top level, i.e., the formula is simple.
*/
static ReturnType preprocess(storm::logic::MultiObjectiveFormula const& originalFormula, SparseModelType const& originalModel);
private:
/*!
* Inserts the information regarding the given formula (i.e. objective) into info.objectives
*
* @param formula OperatorFormula representing the objective
* @param the information collected so far
*/
static void addObjective(std::shared_ptr<storm::logic::Formula const> const& formula, ReturnType& data);
/*!
* Sets the timeBound for the given objective
*/
static void setStepBoundOfObjective(ObjectiveInformation& currentObjective);
/*!
* Apply the neccessary preprocessing for the given formula.
* @param formula the current (sub)formula
* @param data the information collected so far
* @param currentObjective the currently considered objective. The given formula should be a a (sub)formula of this objective
* @param optionalRewardModelName the reward model name that is considered for the formula (if available)
*/
static void preprocessFormula(storm::logic::ProbabilityOperatorFormula const& formula, ReturnType& info, ObjectiveInformation& currentObjective);
static void preprocessFormula(storm::logic::RewardOperatorFormula const& formula, ReturnType& info, ObjectiveInformation& currentObjective);
static void preprocessFormula(storm::logic::UntilFormula const& formula, ReturnType& info, ObjectiveInformation& currentObjective);
static void preprocessFormula(storm::logic::BoundedUntilFormula const& formula, ReturnType& info, ObjectiveInformation& currentObjective);
static void preprocessFormula(storm::logic::GloballyFormula const& formula, ReturnType& info, ObjectiveInformation& currentObjective);
static void preprocessFormula(storm::logic::EventuallyFormula const& formula, ReturnType& info, ObjectiveInformation& currentObjective, boost::optional<std::string> const& optionalRewardModelName = boost::none);
static void preprocessFormula(storm::logic::CumulativeRewardFormula const& formula, ReturnType& info, ObjectiveInformation& currentObjective, boost::optional<std::string> const& optionalRewardModelName = boost::none);
};
}
}
}
#endif /* STORM_MODELCHECKER_MULTIOBJECTIVE_HELPER_SPARSEMULTIOBJECTIVEPREPROCESSOR_H_ */

26
src/modelchecker/multiobjective/helper/SparseMultiObjectiveModelCheckerInformation.h → src/modelchecker/multiobjective/helper/SparseMultiObjectivePreprocessorReturnType.h
View File

@ -1,9 +1,10 @@
#ifndef STORM_MODELCHECKER_MULTIOBJECTIVE_HELPER_SPARSEMULTIOBJECTIVEMODELCHECKERINFORMATION_H_
#define STORM_MODELCHECKER_MULTIOBJECTIVE_HELPER_SPARSEMULTIOBJECTIVEMODELCHECKERINFORMATION_H_
#ifndef STORM_MODELCHECKER_MULTIOBJECTIVE_HELPER_SPARSEMULTIOBJECTIVEPREPROCESSORRETURNTYPE_H_
#define STORM_MODELCHECKER_MULTIOBJECTIVE_HELPER_SPARSEMULTIOBJECTIVEPREPROCESSORRETURNTYPE_H_
#include <vector>
#include <memory>
#include <iomanip>
#include <boost/optional.hpp>
#include "src/logic/Formulas.h"
@ -28,7 +29,7 @@ namespace storm {
bool thresholdIsStrict = false;
boost::optional<uint_fast64_t> stepBound;
void printInformationToStream(std::ostream& out) const {
void printToStream(std::ostream& out) const {
out << std::setw(30) << originalFormula->toString();
out << " \t(";
out << (isNegative ? "-x, " : " ");
@ -52,26 +53,19 @@ namespace storm {
}
};
SparseMultiObjectiveModelCheckerInformation(SparseModelType const& model) : preprocessedModel(model), originalModel(model) {
//Intentionally left empty
}
SparseModelType preprocessedModel;
SparseModelType const& originalModel;
std::vector<uint_fast64_t> newToOldStateIndexMapping;
std::vector<ObjectiveInformation> objectives;
void printInformationToStream(std::ostream& out) {
void printToStream(std::ostream& out) {
out << "---------------------------------------------------------------------------------------------------------------------------------------" << std::endl;
out << " Multi-objective Model Checker Information " << std::endl;
out << " Multi-objective Preprocessor Result " << std::endl;
out << "---------------------------------------------------------------------------------------------------------------------------------------" << std::endl;
out << std::endl;
out << "Objectives:" << std::endl;
out << "--------------------------------------------------------------" << std::endl;
for (auto const& obj : objectives) {
obj.printInformationToStream(out);
obj.printToStream(out);
}
out << "--------------------------------------------------------------" << std::endl;
out << std::endl;
@ -84,10 +78,14 @@ namespace storm {
out << "---------------------------------------------------------------------------------------------------------------------------------------" << std::endl;
}
SparseModelType preprocessedModel;
SparseModelType const& originalModel;
std::vector<uint_fast64_t> newToOldStateIndexMapping;
std::vector<ObjectiveInformation> objectives;
};
}
}
}
#endif /* STORM_MODELCHECKER_MULTIOBJECTIVE_HELPER_SPARSEMULTIOBJECTIVEMODELCHECKERINFORMATION_H_ */
#endif /* STORM_MODELCHECKER_MULTIOBJECTIVE_HELPER_SPARSEMULTIOBJECTIVEPREPROCESSORRETURNTYPE_H_ */

2
src/modelchecker/multiobjective/helper/SparseWeightedObjectivesModelCheckerHelper.cpp → src/modelchecker/multiobjective/helper/SparseMultiObjectiveWeightVectorChecker.cpp
View File

@ -144,7 +144,7 @@ namespace storm {
//std::cout << "stateActionRewardVector for objective " << objIndex << " is " << storm::utility::vector::toString(info.preprocessedModel.getRewardModel(info.objectives[objIndex].rewardModelName).getStateActionRewardVector()) << std::endl;
//std::cout << "deterministic state rewards for objective " << objIndex << " are " << storm::utility::vector::toString(deterministicStateRewards) << std::endl;
storm::storage::BitVector statesWithRewards = storm::utility::vector::filter<ValueType>(deterministicStateRewards, [&] (ValueType const& value) -> bool {return !storm::utility::isZero(value);});
storm::storage::BitVector statesWithRewards = ~storm::utility::vector::filterZero(deterministicStateRewards);
// As target states, we take the states from which no reward is reachable.
STORM_LOG_WARN("TODO: target state selection is currently only valid for reachability properties...");
//TODO: we should be able to give some hint to the solver..

12
src/modelchecker/multiobjective/helper/SparseWeightedObjectivesModelCheckerHelper.h → src/modelchecker/multiobjective/helper/SparseMultiObjectiveWeightVectorChecker.h
View File

@ -1,9 +1,9 @@
#ifndef STORM_MODELCHECKER_MULTIOBJECTIVE_HELPER_SPARSEWEIGHTEDOBJECTIVESMODELCHECKERHELPER_H_
#define STORM_MODELCHECKER_MULTIOBJECTIVE_HELPER_SPARSEWEIGHTEDOBJECTIVESMODELCHECKERHELPER_H_
#ifndef STORM_MODELCHECKER_MULTIOBJECTIVE_HELPER_SPARSEMULTIOBJECTIVEWEIGHTVECTORCHECKER_H_
#define STORM_MODELCHECKER_MULTIOBJECTIVE_HELPER_SPARSEMULTIOBJECTIVEWEIGHTVECTORCHECKER_H_
#include <vector>
#include "src/modelchecker/multiobjective/helper/SparseMultiObjectiveModelCheckerInformation.h"
#include "src/modelchecker/multiobjective/helper/SparseMultiObjectivePreprocessorReturnType.h"
#include "src/storage/TotalScheduler.h"
namespace storm {
@ -11,13 +11,13 @@ namespace storm {
namespace helper {
/*!
* Helper Class that takes a MultiObjectiveInformation and a weight vector and ...
* Helper Class that takes preprocessed multi objective data and a weight vector and ...
* - computes the maximal expected reward w.r.t. the weighted sum of the rewards of the individual objectives
* - extracts the scheduler that induces this maximum
* - computes for each objective the value induced by this scheduler
*/
template <class SparseModelType>
class SparseWeightedObjectivesModelCheckerHelper {
class SparseMultiObjectiveWeightVectorChecker {
public:
typedef typename SparseModelType::ValueType ValueType;
typedef typename SparseModelType::RewardModelType RewardModelType;
@ -87,4 +87,4 @@ namespace storm {
}
}
#endif /* STORM_MODELCHECKER_MULTIOBJECTIVE_HELPER_SPARSEWEIGHTEDOBJECTIVESMODELCHECKERHELPER_H_ */
#endif /* STORM_MODELCHECKER_MULTIOBJECTIVE_HELPER_SPARSEMULTIOBJECTIVEWEIGHTEDVECTORCHECKER_H_ */

12
src/utility/vector.h
View File

@ -405,6 +405,18 @@ namespace storm {
return filter(values, fnc);
}
/*!
* Retrieves a bit vector containing all the indices for which the value at this position is equal to zero
*
* @param values The vector of values.
* @return The resulting bit vector.
*/
template<class T>
storm::storage::BitVector filterZero(std::vector<T> const& values) {
std::function<bool (T const&)> fnc = [] (T const& value) -> bool { return storm::utility::isZero(value); };
return filter(values, fnc);
}
/**
* Sum the entries from values that are set to one in the filter vector.
* @param values

Write Preview
Loading…
Cancel
Save