artur/reactics
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Semi-working version of parametric MC for rsLTL

Browse Source
This commit is contained in:
Artur Meski
2017-09-17 12:24:07 +01:00
parent f7cfeb43c1
commit 9297a7caf7
4 changed files with 188 additions and 117 deletions
Download Patch File Download Diff File Expand all files Collapse all files

21
rs/reaction_system_with_concentrations_param.py
View File

@@ -14,6 +14,8 @@ class ParameterObj(object):
def is_param(some_object): def is_param(some_object):
if isinstance(some_object, ParameterObj): if isinstance(some_object, ParameterObj):
return True return True
else:
return False
class ReactionSystemWithConcentrationsParam(ReactionSystem): class ReactionSystemWithConcentrationsParam(ReactionSystem):
@@ -21,7 +23,7 @@ class ReactionSystemWithConcentrationsParam(ReactionSystem):
self.reactions = [] self.reactions = []
self.parameters = dict() self.parameters = dict()
self.parametric_reactions = [] # self.parametric_reactions = []
self.meta_reactions = dict() self.meta_reactions = dict()
self.permanent_entities = dict() self.permanent_entities = dict()
self.background_set = [] self.background_set = []
@@ -170,9 +172,6 @@ class ReactionSystemWithConcentrationsParam(ReactionSystem):
raise RuntimeError("No products defined") raise RuntimeError("No products defined")
reaction = self.process_rip(R,I,P) reaction = self.process_rip(R,I,P)
if self.is_parametric_reaction(reaction):
self.parametric_reactions.append(reaction)
else:
self.reactions.append(reaction) self.reactions.append(reaction)
def add_reaction_without_reactants(self, R, I, P): def add_reaction_without_reactants(self, R, I, P):
@@ -255,18 +254,6 @@ class ReactionSystemWithConcentrationsParam(ReactionSystem):
else: else:
raise RuntimeError("Unknown meta-reaction type: " + repr(r_type)) raise RuntimeError("Unknown meta-reaction type: " + repr(r_type))
def show_param_reactions(self, soft=False):
print(C_MARK_INFO + " Parametric reactions:")
if soft and len(self.parametric_reactions) > 50:
print(" -> there are more than 50 reactions (" + str(len(self.parametric_reactions)) + ")")
else:
print(" "*4 + "{0: ^35}{1: ^25}{2: ^15}".format("reactants"," inhibitors"," products"))
for reactants,inhibitors,products in self.parametric_reactions:
print(" " + "- {0: ^35}{1: ^25}{2: ^15}".format(
"{ " + self.state_to_str(reactants) + " }",
" { " + self.state_to_str(inhibitors) + " }",
" { " + self.state_to_str(products) + " }"))
def show_max_concentrations(self): def show_max_concentrations(self):
print(C_MARK_INFO + " Maximal allowed concentration levels:") print(C_MARK_INFO + " Maximal allowed concentration levels:")
for e,max_conc in self.max_conc_per_ent.items(): for e,max_conc in self.max_conc_per_ent.items():
@@ -280,7 +267,7 @@ class ReactionSystemWithConcentrationsParam(ReactionSystem):
def show(self, soft=False): def show(self, soft=False):
self.show_background_set() self.show_background_set()
self.show_reactions(soft) self.show_reactions(soft)
self.show_param_reactions(soft) # self.show_param_reactions(soft)
self.show_permanent_entities() self.show_permanent_entities()
self.show_meta_reactions() self.show_meta_reactions()
self.show_max_concentrations() self.show_max_concentrations()

38
rs_testing.py
View File

@@ -13,11 +13,42 @@ def run_tests():
# scalable_chain(print_system=True) # scalable_chain(print_system=True)
# example44() # example44()
# example44_param() # example44_param()
trivial_param() simple_param()
def trivial_param(): def trivial_param():
r = ReactionSystemWithConcentrationsParam()
r.add_bg_set_entity(("x", 3))
r.add_bg_set_entity(("c", 3))
r.add_bg_set_entity(("final", 1))
param_p1 = r.get_param("P1")
r.add_reaction(param_p1, [("c", 1)], [("final", 1)])
# r.add_reaction([("x", 1)], [("c", 1)], [("final", 1)])
c = ContextAutomatonWithConcentrations(r)
c.add_init_state("0")
c.add_state("1")
c.add_transition("0", [("x", 1)], "1")
c.add_transition("1", [], "1")
rc = ReactionSystemWithAutomaton(r, c)
rc.show()
smt_rsc = SmtCheckerRSCParam(rc)
f1 = Formula_rsLTL.f_F(
BagDescription.f_TRUE(),
BagDescription.f_entity("final") >= 1)
#
# WARNING: depth limit is set
#
smt_rsc.check_rsltl(formula=f1, max_level=10, print_witness=True)
def simple_param():
r = ReactionSystemWithConcentrationsParam() r = ReactionSystemWithConcentrationsParam()
r.add_bg_set_entity(("x", 3)) r.add_bg_set_entity(("x", 3))
r.add_bg_set_entity(("y", 3)) r.add_bg_set_entity(("y", 3))
@@ -29,7 +60,8 @@ def trivial_param():
r.add_reaction([("x", 1)], [("c", 1)], [("y", 2)]) r.add_reaction([("x", 1)], [("c", 1)], [("y", 2)])
r.add_reaction([("y", 1)], [("c", 1)], [("z", 1)]) r.add_reaction([("y", 1)], [("c", 1)], [("z", 1)])
r.add_reaction(param_p1, [("c", 1)], [("final", 1)]) # r.add_reaction([("y", 1)], [("c", 1)], r.get_param("P2"))
r.add_reaction(param_p1, [("x", 1)], [("final", 1)])
c = ContextAutomatonWithConcentrations(r) c = ContextAutomatonWithConcentrations(r)
c.add_init_state("0") c.add_init_state("0")

2
rssmt.py
View File

@@ -30,6 +30,7 @@ Author: Artur Męski <meski@ipipan.waw.pl> / <artur.meski@ncl.ac.uk>
################################################################## ##################################################################
def print_banner(): def print_banner():
print() print()
for line in rsmc_banner.split("\n"): for line in rsmc_banner.split("\n"):
@@ -38,6 +39,7 @@ def print_banner():
################################################################## ##################################################################
def main(): def main():
"""Main function""" """Main function"""

196
smt/smt_checker_rsc_param.py
View File

@@ -13,8 +13,8 @@ from logics import rsLTL_Encoder
from rs.reaction_system_with_concentrations_param import ParameterObj, is_param from rs.reaction_system_with_concentrations_param import ParameterObj, is_param
# def simplify(x): def simplify(x):
# return x return x
def z3_max(a, b): def z3_max(a, b):
return If(a > b, a, b) return If(a > b, a, b)
@@ -44,7 +44,8 @@ class SmtCheckerRSCParam(object):
self.v_improd = [] self.v_improd = []
self.v_improd_for_entities = [] self.v_improd_for_entities = []
self.v_param = [] # parameters:
self.v_param = dict()
self.next_level_to_encode = 0 self.next_level_to_encode = 0
@@ -58,11 +59,11 @@ class SmtCheckerRSCParam(object):
# improducible - entities that are never produces (there is no reaction that produces that entity) # improducible - entities that are never produces (there is no reaction that produces that entity)
self.loop_position = Int("loop_position") self.loop_position = Int("loop_position")
self.solver = Solver() #For("QF_FD") self.solver = Solver() #For("QF_FD")
self.verification_time = None self.verification_time = None
self.prepare_param_variables()
def reset(self): def reset(self):
"""Reinitialises the state of the checker""" """Reinitialises the state of the checker"""
@@ -74,7 +75,6 @@ class SmtCheckerRSCParam(object):
self.prepare_state_variables() self.prepare_state_variables()
self.prepare_context_variables() self.prepare_context_variables()
self.prepare_intermediate_product_variables() self.prepare_intermediate_product_variables()
self.prepare_param_variables()
self.next_level_to_encode += 1 self.next_level_to_encode += 1
def prepare_context_variables(self): def prepare_context_variables(self):
@@ -133,9 +133,24 @@ class SmtCheckerRSCParam(object):
reaction_id = self.rs.reactions.index(reaction) reaction_id = self.rs.reactions.index(reaction)
entities_dict = dict() entities_dict = dict()
if is_param(products):
for entity in self.rs.set_of_bgset_ids:
entity_name = self.rs.get_entity_name(entity)
new_var = Int("L"+ str(level) + "_ImProd" + "_r" +
str(reaction_id) + "_" + entity_name)
entities_dict[entity] = new_var
all_entities_dict.setdefault(entity, [])
all_entities_dict[entity].append(new_var)
else:
for entity, conc in products: for entity, conc in products:
new_var = Int("IP" + str(level) + "_R" + entity_name = self.rs.get_entity_name(entity)
str(reaction_id) + "_e" + str(entity)) new_var = Int("L"+ str(level) + "_ImProd" + "_r" +
str(reaction_id) + "_" + entity_name)
entities_dict[entity] = new_var entities_dict[entity] = new_var
all_entities_dict.setdefault(entity, []) all_entities_dict.setdefault(entity, [])
@@ -146,10 +161,6 @@ class SmtCheckerRSCParam(object):
self.v_improd.append(reactions_dict) self.v_improd.append(reactions_dict)
self.v_improd_for_entities.append(all_entities_dict) self.v_improd_for_entities.append(all_entities_dict)
# print(self.v_improd)
# print(self.v_improd_for_entities)
def prepare_param_variables(self): def prepare_param_variables(self):
""" """
Prepares variables for parameters Prepares variables for parameters
@@ -159,42 +170,55 @@ class SmtCheckerRSCParam(object):
background set. background set.
""" """
level = self.next_level_to_encode
param_vars_for_cur_level = dict()
for param_name in self.rs.parameters.keys(): for param_name in self.rs.parameters.keys():
# we start collecting bg-related vars for the given param # we start collecting bg-related vars for the given param
vars_for_param = [] vars_for_param = []
for entity in self.rs.background_set: for entity in self.rs.set_of_bgset_ids:
new_var = Int("L{:d}_Pm_{:s}".format(level, entity)) new_var = Int("Pm{:s}_{:s}".format(param_name, self.rs.get_entity_name(entity)))
vars_for_param.append(new_var) vars_for_param.append(new_var)
param_vars_for_cur_level[param_name] = vars_for_param self.v_param[param_name] = vars_for_param
self.v_param.append(param_vars_for_cur_level) def enc_param_concentration_levels_assertion(self):
"""
Assertions for the parameter variables
"""
enc_param_gz = True
enc_param_at_least_one = False
for param_vars in self.v_param.values():
for pvar in param_vars:
#
# TODO: fixed upper limit: 100
#
enc_param_gz = simplify(And(enc_param_gz, pvar >= 0, pvar < 100))
enc_param_at_least_one = simplify(Or(enc_param_at_least_one, pvar > 0))
print(self.v_param) return simplify(And(enc_param_gz, enc_param_at_least_one))
def enc_concentration_levels_assertion(self, level): def enc_concentration_levels_assertion(self, level):
""" """
Encodes assertions that (some) variables need to be >0 Encodes assertions that (some) variables need to be >=0
We do not need to actually control all the variables, We do not need to actually control all the variables,
only those that can possibly go below 0. only those that can possibly go below 0.
""" """
enc_nz = True enc_gz = True
for e_i in range(len(self.rs.background_set)): for e_i in self.rs.set_of_bgset_ids:
v = self.v[level][e_i] var = self.v[level][e_i]
v_ctx = self.v_ctx[level][e_i] var_ctx = self.v_ctx[level][e_i]
e_max = self.rs.get_max_concentration_level(e_i) e_max = self.rs.get_max_concentration_level(e_i)
enc_nz = simplify(And(enc_nz, v >= 0, v_ctx >= 0, v <= e_max, v_ctx <= e_max)) enc_gz = simplify(And(enc_gz, var >= 0, var_ctx >= 0, var <= e_max, var_ctx <= e_max))
return enc_nz vars_per_reaction = self.v_improd_for_entities[level + 1]
if e_i in vars_per_reaction:
for var_improd in vars_per_reaction[e_i]:
enc_gz = simplify(And(enc_gz, var_improd >= 0, var_improd <= e_max))
return enc_gz
def enc_init_state(self, level): def enc_init_state(self, level):
"""Encodes the initial state at the given level""" """Encodes the initial state at the given level"""
@@ -233,39 +257,45 @@ class SmtCheckerRSCParam(object):
# we need reaction_id to find the intermediate product variable # we need reaction_id to find the intermediate product variable
reaction_id = self.rs.reactions.index(reaction) reaction_id = self.rs.reactions.index(reaction)
# ** REACTANTS ** # ** REACTANTS *******************************************
enc_reactants = True enc_reactants = True
if is_param(reactants):
param_name = reactants.name
for entity in self.rs.set_of_bgset_ids:
enc_reactants = And(enc_reactants, Or(
self.v[level][entity] >= self.v_param[param_name][entity],
self.v_ctx[level][entity] >= self.v_param[param_name][entity]))
else:
for entity, conc in reactants: for entity, conc in reactants:
enc_reactants = And(enc_reactants, Or( enc_reactants = And(enc_reactants, Or(
self.v[level][entity] >= conc, self.v_ctx[level][entity] >= conc)) self.v[level][entity] >= conc,
self.v_ctx[level][entity] >= conc))
# ** INHIBITORS **
# ** INHIBITORS ******************************************
enc_inhibitors = True enc_inhibitors = True
if is_param(inhibitors):
param_name = inhibitors.name
for entity in self.rs.set_of_bgset_ids:
enc_inhibitors = And(enc_inhibitors, Or(
self.v[level][entity] < self.v_param[param_name][entity],
self.v_ctx[level][entity] < self.v_param[param_name][entity]))
else:
for entity, conc in inhibitors: for entity, conc in inhibitors:
enc_inhibitors = And(enc_inhibitors, And( enc_inhibitors = And(enc_inhibitors, And(
self.v[level][entity] < conc, self.v_ctx[level][entity] < conc)) self.v[level][entity] < conc,
self.v_ctx[level][entity] < conc))
# ** PRODUCTS ** # ** PRODUCTS *******************************************
# entities that were produced
produced = []
enc_products = True enc_products = True
if is_param(products):
param_name = products.name
for entity in self.rs.set_of_bgset_ids:
enc_products = And(enc_products,
self.v_improd[level + 1][reaction_id][entity] == self.v_param[param_name][entity])
else:
for entity, conc in products: for entity, conc in products:
enc_products = And(enc_products, self.v_improd[ enc_products = And(enc_products,
level + 1][reaction_id][entity] == conc) self.v_improd[level + 1][reaction_id][entity] == conc)
produced.append(entity)
print("----", str(produced))
# encoding of the zeros (for the products)
# we iterate through all the entities' indices
# for entity in range(len(self.rs.background_set)):
# if entity not in produced:
# enc_products = And(enc_products, self.v_improd[
# level + 1][reaction_id][entity] == 0)
# #
# (R and I) iff P # (R and I) iff P
@@ -275,20 +305,25 @@ class SmtCheckerRSCParam(object):
return enc_reaction return enc_reaction
def enc_single_param_reaction(self, level, p_reaction): def enc_general_reaction_enabledness(self, level):
"""
General enabledness condition for reactions
reactants, inhibitors, products = p_reaction The necessary condition for a reaction to be enabled is
that the state is not empty, i.e., at least one entity
is present in the current state.
if is_param(reactants): This condition must be used when there are parametric reactions
print("PARAM") because parameters could have all the entities set to zero and
that immediately allows for all the conditions on the reactants
to be fulfilled: (entity <= param) -> (0 <= 0)
"""
if is_param(inhibitors): enc_cond = False
print("PARAM") for entity in self.rs.set_of_bgset_ids:
enc_cond = simplify(Or(enc_cond, self.v[level][entity] > 0, self.v_ctx[level][entity] > 0))
if is_param(products): return enc_cond
print("PARAM")
return True
def enc_rs_trans(self, level): def enc_rs_trans(self, level):
"""Encodes the transition relation""" """Encodes the transition relation"""
@@ -313,11 +348,6 @@ class SmtCheckerRSCParam(object):
enc_reaction = self.enc_single_reaction(level, reaction) enc_reaction = self.enc_single_reaction(level, reaction)
enc_trans = simplify(And(enc_trans, enc_reaction)) enc_trans = simplify(And(enc_trans, enc_reaction))
for p_reaction in self.rs.parametric_reactions:
print(p_reaction)
enc_p_reaction = self.enc_single_param_reaction(level, p_reaction)
enc_trans = simplify(And(enc_trans, enc_p_reaction))
# Next we encode the MAX concentration values: # Next we encode the MAX concentration values:
# we collect those from the intermediate product variables # we collect those from the intermediate product variables
@@ -334,17 +364,16 @@ class SmtCheckerRSCParam(object):
# {products}[level+1] # {products}[level+1]
# #
current_v_improd_for_entities = self.v_improd_for_entities[level + 1] current_v_improd_for_entities = self.v_improd_for_entities[level + 1]
for entity, per_reaction_vars in current_v_improd_for_entities.items(): for entity, per_reaction_vars in current_v_improd_for_entities.items():
enc_max_prod = simplify( enc_max_prod = simplify(
And(enc_max_prod, self.v[level + 1][entity] == self.enc_max(per_reaction_vars))) And(enc_max_prod, self.v[level + 1][entity] == self.enc_max(per_reaction_vars)))
# for entity in self.improducible_entities: # make sure at least one entity is >0
# enc_max_prod = simplify( enc_general_cond = self.enc_general_reaction_enabledness(level)
# And(enc_max_prod, self.v[level + 1][entity] == 0))
enc_trans_with_max = simplify(And(enc_max_prod, enc_trans)) enc_trans_with_max = simplify(And(enc_general_cond, enc_max_prod, enc_trans))
print(enc_trans_with_max)
return enc_trans_with_max return enc_trans_with_max
@@ -483,9 +512,29 @@ class SmtCheckerRSCParam(object):
end="") end="")
print(" }") print(" }")
# Parameters
print("\n\n Parameters:\n")
for param_name in self.rs.parameters.keys():
print("{: >6}: ".format(param_name), end="")
print("{", end="")
params = self.v_param[param_name]
for entity in self.rs.set_of_bgset_ids:
var_rep = repr(m[params[entity]])
if not var_rep.isdigit():
raise RuntimeError(
"unexpected: representation is not a positive integer")
if int(var_rep) > 0:
print(
" " + str(self.rs.get_entity_name(entity)) + "=" + str(var_rep),
end="")
print(" }")
def check_rsltl( def check_rsltl(
self, formula, print_witness=True, print_time=True, print_mem=True, self, formula,
print_witness=True, print_time=True, print_mem=True,
max_level=None): max_level=None):
"""Bounded Model Checking for rsLTL properties""" """Bounded Model Checking for rsLTL properties"""
@@ -506,6 +555,7 @@ class SmtCheckerRSCParam(object):
self.prepare_all_variables() self.prepare_all_variables()
self.solver.add(self.enc_concentration_levels_assertion(0)) self.solver.add(self.enc_concentration_levels_assertion(0))
self.solver.add(self.enc_param_concentration_levels_assertion())
encoder = rsLTL_Encoder(self) encoder = rsLTL_Encoder(self)