Meta-reactions

rctsys.py

@@ -37,6 +37,10 @@ class ContextAutomaton(object):
         else:
             print("\'%s\' already added. skipping..." % (name,))
     
+    def add_states(self, states_set):
+        for st in states_set:
+            self.add_state(st)
+    
     def add_init_state(self, name):
         self.add_state(name)
         self._init_state = self._states.index(name)
@@ -59,6 +63,9 @@ class ContextAutomaton(object):
             print("Undefined context automaton state: " + repr(name))
             exit(1)
 
+    def get_state_name(self, state_id):
+        return self._states[state_id]
+
     def get_init_state_id(self):
         return self._init_state
     
@@ -166,6 +173,26 @@ class ContextAutomatonWithConcentrations(ContextAutomaton):
         
         self._transitions.append((self.get_state_id(src),new_context_set,self.get_state_id(dst)))
 
+    def get_automaton_with_flat_contexts(self, ordinary_reaction_system):
+        
+        ca = ContextAutomaton(ordinary_reaction_system)
+        ca._states = self._states
+        ca._init_state = self._init_state
+        
+        for src,ctx,dst in self._transitions:
+            
+            new_ctx = set()
+            
+            for ent,conc in ctx:
+                for i in range(1,conc+1):
+                    n = self._reaction_system.get_entity_name(ent) + "_" + str(i)
+                    ca._reaction_system.ensure_bg_set_entity(n)
+                    new_ctx.add(n)
+            
+            ca.add_transition(ca.get_state_name(src),new_ctx,ca.get_state_name(dst))
+            
+        return ca 
+
 class ReactionSystem(object):
 
     def __init__(self):
@@ -219,13 +246,11 @@ class ReactionSystem(object):
         """Returns the string corresponding to the entity"""
         return self.background_set[entity_id]
 
-
     def add_reaction(self, R, I, P):
         """Adds a reaction"""
         
         if R == [] or P == []:
-            print("No reactants of products defined")
+            raise RuntimeError("No reactants of products defined")
-            raise
 
         reactants = []
         for entity in R:
@@ -348,6 +373,7 @@ class ReactionSystemWithConcentrations(ReactionSystem):
     def __init__(self):
 
         self.reactions = []
+        self.meta_reactions = dict()
         self.background_set = []
 
         self.context_entities = []
@@ -380,12 +406,11 @@ class ReactionSystemWithConcentrations(ReactionSystem):
         if elem[1] < 1:
             raise RuntimeError("Unexpected concentration level in state: " + str(elem))
 
-    def add_reaction(self, R, I, P):
+    def process_rip(self, R, I, P):
-        """Adds a reaction"""
+        """Chcecks concentration levels and converts entities names into their ids"""
 
-        if R == [] or P == []:
+        if R == []:
-            print("No reactants of products defined")
+            raise RuntimeError("No reactants defined")
-            raise
         
         reactants = []
         for r in R:
@@ -395,7 +420,6 @@ class ReactionSystemWithConcentrations(ReactionSystem):
             reactants.append((self.get_entity_id(entity),level))
             if self.max_concentration < level:
                 self.max_concentration = level
-            
         inhibitors = []
         for i in I:
             self.is_valid_entity_with_concentration(i)
@@ -404,7 +428,6 @@ class ReactionSystemWithConcentrations(ReactionSystem):
             inhibitors.append((self.get_entity_id(entity),level))
             if self.max_concentration < level:
                 self.max_concentration = level
-
         products = []
         for p in P:
             self.is_valid_entity_with_concentration(p)
@@ -412,7 +435,23 @@ class ReactionSystemWithConcentrations(ReactionSystem):
             entity,level = p
             products.append((self.get_entity_id(entity),level))
         
-        self.reactions.append((reactants, inhibitors, products))
+        return reactants,inhibitors,products
+
+    def add_reaction(self, R, I, P):
+        """Adds a reaction"""
+        
+        if P == []:
+            raise RuntimeError("No products defined")
+        reaction = self.process_rip(R,I,P)
+        self.reactions.append(reaction)
+
+    def add_reaction_inc(self, incr_entity, R, I):
+        """Adds a macro/meta reaction for increasing the value of incr_entity"""
+
+        reactants,inhibitors,products = self.process_rip(R,I,[])
+        incr_entity_id = self.get_entity_id(incr_entity)
+        self.meta_reactions.setdefault(incr_entity_id,[])
+        self.meta_reactions[incr_entity_id].append(("inc", reactants, inhibitors))
             
     def set_context_entities(self, entities):
         raise NotImplementedError
@@ -441,9 +480,16 @@ class ReactionSystemWithConcentrations(ReactionSystem):
     def show_background_set(self):
         print("[*] Background set: {" + self.entities_names_set_to_str(self.background_set) + "}")
 
+    def show_meta_reactions(self):
+        print("[*] Meta reactions:")
+        for incr_ent,reactions in self.meta_reactions.items():
+            for r_type,reactants,inhibitors in reactions:
+                print("\t - [ Type=" + repr(r_type) + " Parameter=( " + self.get_entity_name(incr_ent) + " ) ] -- ( R={" + self.state_to_str(reactants) + "}, \tI={" + self.state_to_str(inhibitors) + "} )")
+
     def show(self, soft=False):
         self.show_background_set()
         self.show_reactions(soft)
+        self.show_meta_reactions()
         
     def get_reactions_by_product(self):
         """Sorts reactions by their products and returns a dictionary of products"""
@@ -525,6 +571,17 @@ class ReactionSystemWithAutomaton(object):
     def sanity_check(self):
         pass
         
+    def get_ordinary_reaction_system_with_automaton(self):
+        
+        if not self.is_with_concentrations():
+            raise RuntimeError("Not RS/CA with concentrations")
+        
+        ors = self.rs.get_reaction_system()
+        oca = self.ca.get_automaton_with_flat_contexts(ors)
+        
+        return ReactionSystemWithAutomaton(ors, oca)
+        
+
 # class ReactionSystemWithConcentrationWithAutomaton(ReactionSystemWithAutomaton):
 #
 #     def __init__(self, reaction_system, context_automaton):

rssmt.py

@@ -65,31 +65,37 @@ def main():
     # smt = SmtCheckerDistribRS(drs,debug_level=3)
     # smt.check_reachability(rs_examples.drs_toy_ex1_property1(),max_level=10,print_time=True)
 
+    
     r = ReactionSystemWithConcentrations()
-    r.add_bg_set_entities(["a","b","c","d"])
 
-    r.add_reaction([("c",1)],[("b",2)],[("c",1),("b",1)])
+    r.add_bg_set_entity("e")
-    r.add_reaction([("b",2)],[("c",1)],[("c",5),("b",1)])
+    r.add_bg_set_entity("inc")    
-    r.show()
+    r.add_reaction_inc("e",[("e",1),("inc",1)],[("e",N)])
-    
+    # for i in range(1,N):
-    print("================================================")
+        # r.add_reaction([("e",i),("inc",1)],[("e",N)],[("e",i+1)])
-    
+    # r.show()
-    ordinary_rs = r.get_reaction_system()
-    ordinary_rs.show()
-    
-    # print(r.get_reactions_by_product())
     
     c = ContextAutomatonWithConcentrations(r)
-    c.add_init_state("1")
+    c.add_init_state("init")
-    c.add_transition("1", [("c",1)], "1")
+    c.add_state("working")
-    c.add_transition("1", [("d",2)], "1")
+    c.add_transition("init", [("e",1),("inc",1)], "working")
+    c.add_transition("working", [("inc",1)], "working")
      # c.show()
 
-    # rc = ReactionSystemWithAutomaton(r,c)
+    rc = ReactionSystemWithAutomaton(r,c)
-    #
+    
-    # smt = SmtCheckerRSC(rc)
+    rc.show()
-    #
+    
-    # smt.check_reachability([('c',1)],print_time=True,max_level=20)
+    smt_rsc = SmtCheckerRSC(rc)    
+    smt_rsc.check_reachability([('e',N)],print_time=True,max_level=N)
+
+    # orc = rc.get_ordinary_reaction_system_with_automaton()
+    # orc.show()
+    # smt_tr_rs = SmtCheckerPGRS(orc)
+    # smt_tr_rs.check_reachability(['e_' + str(N)],print_time=True)
+
+    print("Reaction System with Concentrations:", smt_rsc.get_verification_time())
+    # print("Reaction System from translating RSC:", smt_tr_rs.get_verification_time())
 
 if __name__ == "__main__":
     try:

smtcheckerpgrs.py

@@ -22,6 +22,8 @@ class SmtCheckerPGRS(object):
 
         self.solver = Solver()
         
+        self.verification_time = None
+
     def prepare_all_variables(self):
         """Encodes all the variables"""
 
@@ -200,7 +202,7 @@ class SmtCheckerPGRS(object):
                         print(" " + self.rs.get_entity_name(var_id), end="")
                 print(" }")
 
-    def check_reachability(self, state, print_witness=True, print_time=False):
+    def check_reachability(self, state, exact_state=False, print_witness=True, print_time=False):
         """Main testing function"""
 
         if print_time:
@@ -211,14 +213,18 @@ class SmtCheckerPGRS(object):
         current_level = 0
 
         while True:
-            print("\r[i] Level: " + str(current_level), end="")
+            print("-----[ Working at level=" + str(current_level) + " ]-----")
             stdout.flush()
 
             self.prepare_all_variables()
 
             # reachability test:
+            print("[i] Adding the reachability test...")    
             self.solver.push()
+            if exact_state:
                 self.solver.add(self.enc_state(current_level,state))
+            else:
+                self.solver.add(self.enc_non_exclusive_state(current_level,state))
 
             result = self.solver.check()
             if result == sat:
@@ -229,12 +235,17 @@ class SmtCheckerPGRS(object):
             else:
                 self.solver.pop()
            
+            print("[i] Unrolling the transition relation")
             self.solver.add(self.enc_transition_relation(current_level))
 
+            print("-----[ level=" + str(current_level) + " done ]")
             current_level += 1
 
         if print_time:
             stop = time()
+            self.verification_time = stop-start
             print()
-            print("[i] Time: " + repr(stop-start))
+            print("[i] Time: " + repr(self.verification_time))
             
+    def get_verification_time(self):
+        return self.verification_time

smtcheckerrsc.py

@@ -5,6 +5,7 @@ SMT-based Model Checking Module for RS with Concentrations and Context Automaton
 from z3 import *
 from time import time
 from sys import stdout
+from itertools import chain
 
 class SmtCheckerRSC(object):
 
@@ -25,6 +26,8 @@ class SmtCheckerRSC(object):
 
         self.solver = Solver()
         
+        self.verification_time = None
+
     def prepare_all_variables(self):
         """Encodes all the variables"""
 
@@ -72,10 +75,20 @@ class SmtCheckerRSC(object):
     def enc_produced_concentration(self, level, prod_entity):
         """Encodes the produced concentrations for the given level and entity"""
 
+        rcts_for_prod_entity = []
+        if prod_entity in self.rs.get_reactions_by_product():
             rcts_for_prod_entity = self.rs.get_reactions_by_product()[prod_entity]
 
-        if rcts_for_prod_entity == []:
+        meta_reactions = []
-            return simplify(self.v[level+1][prod_entity] == 0)
+        if prod_entity in self.rs.meta_reactions:
+            meta_reactions = self.rs.meta_reactions[prod_entity]
+
+        if rcts_for_prod_entity == [] and meta_reactions == []:
+            return simplify(self.v[level+1][prod_entity] == 0) # this should never happen
+            
+        enc_enabledness = False
+
+        # ----------- ordinary reactions --------------------------------------------
 
         enc_rct_prod = False
         for reactants,inhibitors,products in rcts_for_prod_entity:
@@ -92,8 +105,33 @@ class SmtCheckerRSC(object):
 
             enc_products = self.v[level+1][products[0][0]] == products[0][1]
             
+            enc_enabledness = simplify(Or(enc_enabledness, And(enc_reactants, enc_inhibitors)))
             enc_rct_prod = simplify(Or(enc_rct_prod, And(enc_reactants, enc_inhibitors, enc_products)))
             
+        # -------- meta reactions ---------------------------------------------------
+        
+        for r_type,reactants,inhibitors in meta_reactions:
+            enc_reactants   = True
+            enc_inhibitors  = True
+            
+            for reactant,concentration in reactants:
+                enc_reactants = simplify(And(enc_reactants, 
+                                            Or(self.v[level][reactant] >= concentration, self.v_ctx[level][reactant] >= concentration)))
+            for inhibitor,concentration in inhibitors:
+                enc_inhibitors = simplify(And(enc_inhibitors, 
+                                             And(self.v[level][inhibitor] < concentration, self.v_ctx[level][inhibitor] < concentration)))
+
+            enc_products = self.v[level+1][prod_entity] == self.v[level][prod_entity]+1
+
+            enc_enabledness = simplify(Or(enc_enabledness, And(enc_reactants, enc_inhibitors)))
+
+            enc_rct_prod = simplify(Or(enc_rct_prod, And(enc_reactants, enc_inhibitors, enc_products)))
+            
+        # -----------------------------------------------------------------------------
+            
+        enc_when_to_produce_zero_conc = simplify(And(Not(enc_enabledness), self.v[level+1][prod_entity] == 0))
+
+        enc_rct_prod = Or(enc_rct_prod, enc_when_to_produce_zero_conc)
         return enc_rct_prod
 
     # def enc_entity_production(self, level, prod_entity):
@@ -106,6 +144,7 @@ class SmtCheckerRSC(object):
     #
     #     return simplify(enc_ent_prod)
     
+        
     def enc_transition_relation(self, level):
         return simplify(And(self.enc_rs_trans(level), self.enc_automaton_trans(level)))
 
@@ -117,8 +156,9 @@ class SmtCheckerRSC(object):
         enc_trans = True
 
         reactions = self.rs.get_reactions_by_product()
+        meta_reactions = self.rs.meta_reactions
 
-        for prod_entity in reactions:
+        for prod_entity in chain(reactions, meta_reactions):
             unused_entities.remove(prod_entity)
     
             enc_trans = simplify(And(enc_trans, self.enc_produced_concentration(level, prod_entity)))
@@ -232,7 +272,6 @@ class SmtCheckerRSC(object):
         current_level = 0
 
         self.prepare_all_variables()
-        print(And(self.enc_transition_relation(current_level), self.enc_init_state(0)))
         
         while True:
             self.prepare_all_variables()
@@ -270,6 +309,10 @@ class SmtCheckerRSC(object):
 
         if print_time:
             stop = time()
+            self.verification_time = stop-start
             print()
-            print("[i] Time: " + repr(stop-start))
+            print("[i] Time: " + repr(self.verification_time))
+            
+    def get_verification_time(self):
+        return self.verification_time