initial commit with .py files

bftype.py

@@ -0,0 +1,9 @@
+#from enum import Enum
+
+class opertype:
+    """Boolean Formulae Types"""
+    unspec = 0
+    o_and = 1
+    o_or = 2
+    o_not = 3
+    o_var = 10

boolform.py

@@ -0,0 +1,69 @@
+#!/usr/bin/env python
+
+import bftype
+
+class BoolForm:
+
+    def __init__(self, a=None):
+        self.initAll()
+        if a != None and isinstance(a, int):
+            self.setupVarNode(a)
+
+    def initAll(self):
+        self.oper = bftype.opertype.unspec
+
+    def setupConjunction(self, lhs, rhs):
+        self.oper = bftype.opertype.o_and
+        self.sfLeft = lhs
+        self.sfRight = rhs
+
+    def setupDisjunction(self, lhs, rhs):
+        self.oper = bftype.opertype.o_or
+        self.sfLeft = lhs
+        self.sfRight = rhs
+
+    def setupNegation(self, form):
+        self.oper = bftype.opertype.o_not
+        self.subForm = form
+
+    def setupVarNode(self, var_id):
+        if not isinstance(var_id, int):
+            print("var_id must be an integer")
+            raise
+        self.oper = bftype.opertype.o_var
+        self.varid = var_id
+
+    def __add__(self, x):
+        r = BoolForm()
+        r.setupDisjunction(self,x)
+        return r
+
+    def __mul__(self, x):
+        r = BoolForm()
+        r.setupConjunction(self,x)
+        return r
+
+    def __neg__(self):
+        r = BoolForm()
+        r.setupNegation(self)
+        return r
+
+    def __str__(self):
+        if self.oper == bftype.opertype.o_or:
+            return "(" + str(self.sfLeft) + " or " + str(self.sfRight) + ")"
+        elif self.oper == bftype.opertype.o_and:
+            return "(" + str(self.sfLeft) + " and " + str(self.sfRight) + ")"
+        elif self.oper == bftype.opertype.o_not:
+            return "-" + str(self.subForm)
+        elif self.oper == bftype.opertype.o_var:
+            return str(self.varid)
+        else:
+            print("This should not happen. Unknown operator?")
+            raise
+
+if __name__ == "__main__":
+
+    y = (BoolForm(1) + -BoolForm(2) * BoolForm(3)) + BoolForm(2)
+    
+    print(y)
+

distrib_rctsys.py

@@ -0,0 +1,305 @@
+#!/usr/bin/env python
+
+"""
+--- Distributed Reaction Systems Manipulation
+"""
+
+from sys import exit
+
+class DistributedReactionSystem(object):
+
+    def __init__(self):
+
+        self.__reactions = []
+        self.__background_set = []
+        self.__reactions_by_prod = []
+        self.__states = []
+        self.__transitions = []
+        self.__init_state = None
+
+        self.__context_sets_size = None
+
+    @property
+    def background_set(self):
+        return self.__background_set
+    
+    @property
+    def set_of_background_ids(self):
+        return set(range(len(self.__background_set)))
+
+    @property
+    def components_count(self):
+        return len(self.__reactions)
+
+    def add_bg_set_entity(self, name):
+        if not self.is_in_background_set(name):
+            self.__background_set.append(name)
+        else:
+            print("The entity \"" + str(name) + "\" is already on the list")
+            exit(1)
+
+    def add_bg_set_entities(self, names):
+        for name in names:
+            self.add_bg_set_entity(name)
+
+    def is_in_background_set(self, entity):
+        """Checks if the given name is valid wrt the background set="""
+        if entity in self.__background_set:
+            return True
+        else:
+            return False
+
+    def get_entity_id(self, name):
+        try:
+            return self.__background_set.index(name)
+        except ValueError:
+            print("Undefined background set entity: " + repr(name))
+            exit(1)
+
+    def get_state_ids(self, state):
+        ids = []
+        for entity in state:
+            ids.append(self.get_entity_id(entity))
+
+        return ids
+
+    def get_entity_name(self, entity_id):
+        """Returns the string corresponding to the entity"""
+        return self.__background_set[entity_id]
+
+    def ensure_reactions(self, k):
+        while len(self.__reactions) <= k:
+            self.__reactions.append([])
+
+    def add_reaction(self, k, R, I, P):
+        """Adds a reaction"""
+        
+        if R == [] or P == []:
+            print("No reactants of products defined")
+            raise
+
+        reactants = []
+        for entity in R:
+            reactants.append(self.get_entity_id(entity))
+
+        inhibitors = []
+        for entity in I:
+            inhibitors.append(self.get_entity_id(entity))
+
+        products = []
+        for entity in P:
+            products.append(self.get_entity_id(entity))
+        
+        self.ensure_reactions(k)
+        self.__reactions[k].append((reactants, inhibitors, products))
+        
+    def entities_names_set_to_str(self, entities):
+        s = ""
+        for entity in entities:
+            s += entity + ", "
+        s = s[:-2]
+        return s
+
+    def entities_ids_set_to_str(self, entities):
+        s = ""
+        for entity in entities:
+            s += self.get_entity_name(entity) + ", "
+        s = s[:-2]
+        return s
+
+    def show_reactions(self):
+        print("[*] Reactions:")
+        for i in range(len(self.__reactions)):
+            local_reactions = self.__reactions[i]
+            print("    agent = " + str(i))
+            for rcts,inhib,prods in local_reactions:
+                print("\t - ( R={" + self.entities_ids_set_to_str(rcts) + "}, \tI={" + self.entities_ids_set_to_str(inhib) + "}, \tP={" + self.entities_ids_set_to_str(prods) + "} )")
+
+    def show_background_set(self):
+        print("[*] Background set: {" + self.entities_names_set_to_str(self.__background_set) + "}")
+
+    def show(self):
+        self.show_background_set()
+        self.show_reactions()
+        self.show_states()
+        self.show_transitions()
+        print()
+        
+    def reactions_by_product_are_cached(self, component_id):
+        if len(self.__reactions_by_prod) > component_id:
+            if self.__reactions_by_prod[component_id] != None:
+                return True
+        return False
+        
+    def get_reactions_by_product(self, component_id):
+        """Sorts reactions for the component given by component_id by their products and returns a dictionary of products"""
+        
+        if self.reactions_by_product_are_cached(component_id):
+            return self.__reactions_by_prod[component_id]
+
+        producible_entities = set()
+
+        for reactants,inhibitors,products in self.__reactions[component_id]:
+            producible_entities = producible_entities.union(set(products))
+
+        reactions_by_prod = dict()
+
+        for prod_entity in producible_entities:
+            reactions_by_prod[prod_entity] = []
+            for reactants,inhibitors,products in self.__reactions[component_id]:
+                if prod_entity in products:
+                    reactions_by_prod[prod_entity].append([reactants,inhibitors])
+
+        # save in cache
+        while len(self.__reactions_by_prod) <= component_id: # ensure that the size of cache is right
+            self.__reactions_by_prod.append(None)
+        self.__reactions_by_prod[component_id] = reactions_by_prod
+
+        return reactions_by_prod
+
+    @property
+    def states(self):
+        return self.__states
+    
+    @property
+    def transitions(self):
+        return self.__transitions
+    
+    def add_state(self, name):
+        if name not in self.__states:
+            self.__states.append(name)
+        else:
+            print("\'%s\' already added. skipping..." % (name,))
+            
+    def add_init_state(self, name):
+        self.add_state(name)
+        self.__init_state = self.__states.index(name)
+
+    def get_init_state_name(self):
+        if self.__init_state == None:
+            return None
+        return self.__states[self.__init_state]
+
+    def is_state(self, name):
+        if name in self.__states:
+            return True
+        else:
+            return False
+
+    def get_state_id(self, name):
+        try:
+            return self.__states.index(name)
+        except ValueError:
+            print("Undefined context automaton state: " + repr(name))
+            exit(1)
+
+    def get_init_state_id(self):
+        return self.__init_state
+    
+    def print_states(self):
+         for state in self.__states:
+            print(state)
+    
+    def is_valid_context_sets(self, context_sets):
+        for c in context_sets:
+            if not self.is_valid_context(c):
+                return False
+        return True
+    
+    def is_valid_context(self, context):
+        if set(context).issubset(self.__background_set):
+            return True
+        else:
+            return False
+        
+    def add_transition(self, src, label, dst):
+        recipients = label[0]
+        context_sets = label[1]
+        if not type(context_sets) is list:
+            print("Context sets must be of type list")
+            exit(1)
+        for s in context_sets:
+            if not type(s) is set and not type(s) is list:
+                print("Each context set must be of type set or list")
+                exit(1)
+                
+        if self.__context_sets_size == None:
+            self.__context_sets_size = len(context_sets)
+        else:
+            if len(context_sets) != self.__context_sets_size:
+                print("Inconsistent size of the context sets: " + str(len(context_sets)) + " != " + str(self.__context_sets_size))
+                exit(1)
+                
+        if not self.is_valid_context_sets(context_sets):
+            raise RuntimeError("one of the entities in the context set is unknown (undefined)!")
+            
+        if not self.is_state(src):
+            raise RuntimeError("\"" + src + "\" is an unknown (undefined) state")
+
+        if not self.is_state(dst):
+            raise RuntimeError("\"" + dst + "\" is an unknown (undefined) state")
+
+        new_context_sets = []
+        for c in context_sets:
+            new_c_set = set()
+            for e in set(c):
+                new_c_set.add(self.get_entity_id(e))
+            new_context_sets.append(new_c_set)
+
+        recipients = set(recipients)
+
+        self.__transitions.append((self.get_state_id(src),(recipients,new_context_sets),self.get_state_id(dst)))
+
+    def context2str(self, context_sets):
+        """Converts the set of entities ids into the string with their names"""
+        s = "["
+        for ctx in context_sets:
+            if len(ctx) == 0:
+                s += " 0"
+            else:
+                s += " {"
+                for c in ctx:
+                    s += " " + self.get_entity_name(c) 
+                s += " }"
+        s += " ]"
+        return s
+        
+    def show_transitions(self):
+        print("[*] Context automaton transitions:")
+        for transition in self.__transitions:
+            str_transition = str(transition[0]) + " --( "
+            str_transition += str(list(transition[1][0])) + "<=" + self.context2str(transition[1][1])
+            str_transition += " )--> " + str(transition[2])
+            print("\t- " + str_transition)
+    
+    def show_states(self):
+        init_state_name = self.get_init_state_name()
+        print("[*] Context automaton states:")
+        for state in self.__states:
+            print("\t- " + state, end="")
+            if state == init_state_name:
+                print(" [init]")
+            else:
+                print()
+
+    def sanity_check(self):
+        """Performs a sanity check of the defined distributed reaction system"""
+
+        print("[i] Performing sanity check of the DRS")
+
+        if self.__reactions == []:
+            print("No reactions defined")
+            exit(1)
+
+        if self.__background_set == []:
+            print("Empty background set")
+            exit(1)
+
+        if self.__init_state == None:
+            print("Initial state not specified")
+            exit(1)
+
+        if self.__context_sets_size != len(self.__reactions):
+            print("Inconsistent sizes of the context sets with respect to the number of components/agents/processes: " + str(self.__context_sets_size) + " != " + str(len(self.__reactions)))
+            exit(1)
+

example.py

@@ -0,0 +1,39 @@
+from z3 import *
+
+p0 = Bool('p0')
+p1 = Bool('p1')
+p2 = Bool('p2') 
+p3 = Bool('p3')
+p4 = Bool('p4')
+
+p0e = Bool('p0e')
+p1e = Bool('p1e')
+p2e = Bool('p2e')
+p3e = Bool('p3e')
+p4e = Bool('p4e')
+
+p0p = Bool('p0p')
+p1p = Bool('p1p')
+p2p = Bool('p2p')
+p3p = Bool('p3p')
+p4p = Bool('p4p')
+
+Cd_a1 = And( Or(p1,p1e), Or(p4,p4e), Not( Or(p2,p2e) ) )
+Cd_a2 = And( Or(p2,p2e),             Not( Or(p3,p3e) ) )
+Cd_a3 = And( Or(p1,p1e), Or(p3,p3e), Not( Or(p2,p2e) ) )
+Cd_a4 = And( Or(p3,p3e),             Not( Or(p2,p2e) ) )
+
+En1 = Or(Cd_a1,Cd_a2,Cd_a3,Cd_a4)
+En2 = Or(Cd_a1,Cd_a3)
+En3 = Cd_a2
+En4 = Cd_a2
+
+Pr1 = Or( And(En1,p1p) , And( Not(En1), Not(p1p)) )
+Pr2 = Or( And(En2,p2p) , And( Not(En2), Not(p2p)) )
+Pr3 = Or( And(En3,p3p) , And( Not(En3), Not(p3p)) )
+Pr4 = Or( And(En4,p4p) , And( Not(En4), Not(p4p)) )
+
+PrConjunction = And( Pr1, Pr2, Pr3, Pr4 )
+
+print PrConjunction
+

rctsys.py

@@ -0,0 +1,319 @@
+#!/usr/bin/env python
+
+"""
+--- Reaction Systems Manipulation
+"""
+
+#
+# TODO:
+#
+#  *  zamiast zbiorów - listy
+#  *  index pobieramy za pomoca .index()
+#  *  nie dopuszczamy modyfikacji listy stanow/background_set, itp.
+#  *  uzywamy setterow i getterow
+#
+
+from sys import exit
+
+class ContextAutomaton(object):
+    
+    def __init__(self, reaction_system):
+        self.__states = []
+        self.__transitions = []
+        self.__init_state = None
+        self.__reaction_system = reaction_system
+    
+    @property
+    def states(self):
+        return self.__states
+    
+    @property
+    def transitions(self):
+        return self.__transitions
+    
+    def add_state(self, name):
+        if name not in self.__states:
+            self.__states.append(name)
+        else:
+            print("\'%s\' already added. skipping..." % (name,))
+            
+    def add_init_state(self, name):
+        self.add_state(name)
+        self.__init_state = self.__states.index(name)
+
+    def get_init_state_name(self):
+        if self.__init_state == None:
+            return None
+        return self.__states[self.__init_state]
+
+    def is_state(self, name):
+        if name in self.__states:
+            return True
+        else:
+            return False
+
+    def get_state_id(self, name):
+        try:
+            return self.__states.index(name)
+        except ValueError:
+            print("Undefined context automaton state: " + repr(name))
+            exit(1)
+
+    def get_init_state_id(self):
+        return self.__init_state
+    
+    def print_states(self):
+         for state in self.__states:
+            print(state)
+            
+    def is_valid_context(self, context):
+        if set(context).issubset(self.__reaction_system.background_set):
+            return True
+        else:
+            return False
+        
+    def add_transition(self, src, context_set, dst):
+        if not type(context_set) is set and not type(context_set) is list:
+            print("Contexts set must be of type set or list")
+            
+        if not self.is_valid_context(context_set):
+            raise RuntimeError("one of the entities in the context set is unknown (undefined)!")
+            
+        if not self.is_state(src):
+            raise RuntimeError("\"" + src + "\" is an unknown (undefined) state")
+
+        if not self.is_state(dst):
+            raise RuntimeError("\"" + dst + "\" is an unknown (undefined) state")
+            
+        new_context_set = set()
+        for e in set(context_set):
+            new_context_set.add(self.__reaction_system.get_entity_id(e))
+
+        self.__transitions.append((self.get_state_id(src),new_context_set,self.get_state_id(dst)))
+
+    def context2str(self, ctx):
+        """Converts the set of entities ids into the string with their names"""
+        if len(ctx) == 0:
+            return "0"
+        s = "{"
+        for c in ctx:
+            s += " " + self.__reaction_system.get_entity_name(c) 
+        s += " }"
+        return s
+        
+    def show_transitions(self):
+        print("[*] Context automaton transitions:")
+        for transition in self.__transitions:
+            str_transition = str(transition[0]) + " --( "
+            str_transition += self.context2str(transition[1])
+            str_transition += " )--> " + str(transition[2])
+            print("\t- " + str_transition)
+    
+    def show_states(self):
+        init_state_name = self.get_init_state_name()
+        print("[*] Context automaton states:")
+        for state in self.__states:
+            print("\t- " + state, end="")
+            if state == init_state_name:
+                print(" [init]")
+            else:
+                print()
+
+    def show(self):
+        self.show_states()
+        self.show_transitions()
+
+
+class ReactionSystem(object):
+
+    def __init__(self):
+
+        self.reactions = []
+        self.background_set = []
+        self.init_contexts = []
+        self.context_entities = []
+
+        self.reactions_by_agents = [] # each element is 'reactions_by_prod'
+
+        self.reactions_by_prod = None
+
+    def add_bg_set_entity(self, name):
+        if not self.is_in_background_set(name):
+            self.background_set.append(name)
+        else:
+            print("The entity", name, "is already on the list")
+            raise
+
+    def add_bg_set_entities(self, names):
+        for name in names:
+            self.add_bg_set_entity(name)
+
+    def is_in_background_set(self, entity):
+        """Checks if the given name is valid wrt the background set="""
+        if entity in self.background_set:
+            return True
+        else:
+            return False
+
+    def get_entity_id(self, name):
+        try:
+            return self.background_set.index(name)
+        except ValueError:
+            print("Undefined background set entity: " + repr(name))
+            exit(1)
+
+    def get_state_ids(self, state):
+        ids = []
+        for entity in state:
+            ids.append(self.get_entity_id(entity))
+
+        return ids
+
+    def get_entity_name(self, entity_id):
+        """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
+
+        reactants = []
+        for entity in R:
+            reactants.append(self.get_entity_id(entity))
+
+        inhibitors = []
+        for entity in I:
+            inhibitors.append(self.get_entity_id(entity))
+
+        products = []
+        for entity in P:
+            products.append(self.get_entity_id(entity))
+
+        self.reactions.append((reactants, inhibitors, products))
+
+    def add_initial_context_set(self, context_set):
+        if context_set == []:
+            print("Empty context set is not allowed")
+            raise
+
+        integers = []
+        for entity in context_set:
+            if not entity in self.background_set:
+                print("The entity", entity, "is not in the background set")
+                raise
+            else:
+                integers.append(self.get_entity_id(entity))
+
+        self.init_contexts.append(integers)
+
+    def set_context_entities(self, entities):
+
+        for entity in entities:
+            entity_id = self.get_entity_id(entity)
+            self.context_entities.append(entity_id)
+
+    def entities_names_set_to_str(self, entities):
+        s = ""
+        for entity in entities:
+            s += entity + ", "
+        s = s[:-2]
+        return s
+
+    def entities_ids_set_to_str(self, entities):
+        s = ""
+        for entity in entities:
+            s += self.get_entity_name(entity) + ", "
+        s = s[:-2]
+        return s
+
+    def show_reactions(self):
+        print("[*] Reactions:")
+        for reaction in self.reactions:
+            print("\t - ( R={" + self.entities_ids_set_to_str(reaction[0]) + "}, \tI={" + self.entities_ids_set_to_str(reaction[1]) + "}, \tP={" + self.entities_ids_set_to_str(reaction[2]) + "} )")
+
+    def show_background_set(self):
+        print("[*] Background set: {" + self.entities_names_set_to_str(self.background_set) + "}")
+
+    def show_initial_contexts(self):
+        if len(self.init_contexts) > 0:
+            print("[*] Initial context sets:")
+            for ctx in self.init_contexts:
+                  print("\t - {" + self.entities_ids_set_to_str(ctx) + "}")
+
+    def show_context_entities(self):
+        if len(self.context_entities) > 0:
+            print("[*] Context entities: " + self.entities_ids_set_to_str(self.context_entities))
+
+    def show(self):
+
+        self.show_background_set()
+        self.show_initial_contexts()
+        self.show_reactions()
+        self.show_context_entities()
+        
+    def get_reactions_by_product(self):
+        """Sorts reactions by their products and returns a dictionary of products"""
+
+        if self.reactions_by_prod != None:
+            return self.reactions_by_prod
+
+        producible_entities = set()
+
+        for reaction in self.reactions:
+            producible_entities = producible_entities.union(set(reaction[2]))
+
+        reactions_by_prod = {} 
+
+        for prod_entity in producible_entities:
+            reactions_by_prod[prod_entity] = []
+            for reaction in self.reactions:
+                if prod_entity in reaction[2]:
+                    reactions_by_prod[prod_entity].append([reaction[0],reaction[1]])
+
+        # save in cache
+        self.reactions_by_prod = reactions_by_prod
+
+        return reactions_by_prod
+
+    def sanity_check(self):
+        """Performs a sanity check on the defined reaction system"""
+
+        if self.reactions == []:
+            print("No reactions defined")
+            exit(1)
+
+        if self.background_set == []:
+            print("Empty background set")
+            exit(1)
+
+        if self.init_contexts == []:
+            print("No initial context defined")
+            exit(1)
+
+        if self.context_entities == []:
+            print("WARNING: no context entities!")
+
+
+class ReactionSystemWithAutomaton(object):
+    
+    def __init__(self, reaction_system, context_automaton):
+        self.rs = reaction_system
+        self.ca = context_automaton
+
+    def show(self):
+        self.rs.show()
+        self.ca.show()
+        
+    def sanity_check(self):
+        if self.rs.context_entities != []:
+            raise RuntimeError("Context entities should not be defined within RS!")
+
+        if self.rs.init_contexts != []:
+            raise RuntimeError("Initial context should not be defined!")
+
+        
+        # rs should have no initial context, no context entities, etc.

rs_examples.py

@@ -0,0 +1,202 @@
+#!/usr/bin/env python
+
+from rctsys import ReactionSystem
+from rctsys import ReactionSystemWithAutomaton,ContextAutomaton
+from distrib_rctsys import DistributedReactionSystem
+
+def toy_ex1():
+
+    rs = ReactionSystem()
+
+    #rs.add_bg_set_entities(["a","b","x","c"])
+
+    #rs.add_reaction(["a","b"],["x"],["c"])
+    #rs.add_reaction(["c"],["a"],["a"])
+    #rs.add_reaction(["a"],["x"],["c"])
+    #rs.set_context_entities(["b"])
+
+    #rs.add_initial_context_set(["a","b"])
+
+    #rs.add_bg_set_entities(["a", "b", "c", "d", "e", "f", "x", "y","1","2","3","4","5","6"])
+    rs.add_bg_set_entities(list("qwertyuiopasdfghjklzxcvbnm"))
+
+    rs.add_reaction(["a"], ["x"], ["b"])
+    rs.add_reaction(["b"], ["x"], ["c"])
+    rs.add_reaction(["c","f"], ["x"], ["d"])
+    rs.add_reaction(["d"], ["x"], ["e"])
+
+    rs.set_context_entities(["f","g","h"])
+    rs.add_initial_context_set(["a"])
+
+    return rs
+
+def toy_ex2():
+
+    rs = ReactionSystem()
+
+    rs.add_bg_set_entities(["a","b","c","d","d_comp","e","f","x"])
+    rs.add_reaction(["a"],["x"],["b"])
+    rs.add_reaction(["b"],["a"],["c"])
+    rs.add_reaction(["c"],["x"],["d"])
+    rs.add_reaction(["d","d_comp"],["e"],["f"])
+    rs.add_reaction(["f"],["a"],["e"])
+    rs.set_context_entities(["d_comp", "e", "x"])
+    rs.add_initial_context_set(["a"])
+
+    return rs
+
+def toy_ex3():
+
+    rs = ReactionSystem()
+
+    rs.add_bg_set_entities(["1","2","3","4"])
+    rs.add_reaction(["1","4"],["2"],["1","2"])
+    rs.add_reaction(["2"],["3"],["1","3","4"])
+    rs.add_reaction(["1","3"],["2"],["1","2"])
+    rs.add_reaction(["3"],["2"],["1"])
+    rs.add_initial_context_set(["1","4"])
+    rs.set_context_entities(["4"])
+
+    return rs
+
+def bitctr(bits):
+
+    rs = ReactionSystem()
+
+    n = bits
+
+    for i in range(0,bits):
+        rs.addBgSetEntity("p" + str(i))
+
+    rs.addBgSetEntity("dec")
+    rs.addBgSetEntity("inc")
+
+    # (1) no dec, no inc
+    for j in range(0,bits):
+        rs.add_reaction(["p"+str(j)], ["dec","inc"], ["p"+str(j)])
+
+    # (2) increment
+    rs.add_reaction(["inc"],["dec","p0"],["p0"])
+    for j in range(1,bits):
+        R = ["inc"]
+        for k in range(0,j):
+            R.append("p"+str(k))
+        I = ["dec","p"+str(j)]
+        P = ["p"+str(j)]
+        rs.add_reaction(R, I, P)
+
+    for j in range(0,bits):
+        for k in range(j+1,bits):
+            rs.add_reaction(["inc","p"+str(k)], ["dec","p"+str(j)], ["p"+str(k)])
+
+    # (3) decrement
+    for j in range(0,bits):
+        R=["dec"]
+        I=["inc"]
+        for k in range(0,j+1):
+            I.append("p"+str(k))
+        P=["p"+str(j)]
+        rs.add_reaction(R, I, P)
+
+    for j in range(0,bits):
+        for k in range(j+1,bits):
+            rs.add_reaction(["dec","p"+str(j),"p"+str(k)], ["inc"], ["p"+str(k)])
+
+    rs.set_context_entities(["dec","inc"])
+    rs.add_initial_context_set(["inc"])
+    
+    return rs
+    
+
+def ca_toy_ex1():
+
+    rs = ReactionSystem()
+
+    rs.add_bg_set_entities(list("axbcdfegh"))
+
+    rs.add_reaction(["a"], ["x"], ["b"])
+    rs.add_reaction(["b"], ["x"], ["c"])
+    rs.add_reaction(["c","f"], ["x"], ["d"])
+    rs.add_reaction(["d"], ["x"], ["e"])
+
+    #rs.set_context_entities(["f","g","h"])
+    #rs.add_initial_context_set(["a"])
+
+    ca = ContextAutomaton(rs)
+    ca.add_init_state("1")
+    ca.add_state("2")
+    ca.add_transition("1", ["a"], "2")
+    ca.add_transition("2", [], "2")
+    ca.add_transition("2", ["f"], "2")
+    ca.add_transition("2", ["x"], "2")
+
+    rsca = ReactionSystemWithAutomaton(rs,ca)
+
+    return rsca
+
+def ca_toy_ex1_property1():
+    return ["e"]
+
+
+def drs_toy_ex1():
+    
+    drs = DistributedReactionSystem()
+    drs.add_bg_set_entities(list("axbcd"))
+    drs.add_reaction(0, ["a"], ["x"], ["b"])
+    drs.add_reaction(0, ["b"], ["c"], ["d"])
+    drs.add_reaction(1, ["a"], [], ["b"])
+    
+    drs.add_init_state("init")
+    drs.add_state("working")
+    drs.add_transition("init", ([0],[ ["a"], ["a"] ]), "working")
+    drs.add_transition("working", ([1],[ ["b"], ["a"] ]), "working")
+    
+    return drs
+
+def drs_toy_ex1_property1():
+    return [["b"],["b"]]
+    
+def drs_mutex(k):
+    
+    drs = DistributedReactionSystem()
+
+    ctrl = k
+    # agents = 0 ... (k-1) 
+
+    drs.add_bg_set_entities(["lock", "busy", "release", "req", "out", "in"])
+    
+    drs.add_reaction(ctrl, ["lock"], ["release"], ["lock"])
+    drs.add_reaction(ctrl, ["req"], [], ["lock"])
+
+    for i in range(k):
+        #drs.add_reaction(i, ["out"], [], ["req"])
+        drs.add_reaction(i, ["out","busy"], [], ["out"])
+        drs.add_reaction(i, ["out"], [], ["req"])
+        drs.add_reaction(i, ["req"], ["lock"], ["in"])
+        drs.add_reaction(i, ["in"], ["busy"], ["out","release"])
+        drs.add_reaction(i, ["in","busy"], [], ["in"])
+        drs.add_reaction(i, ["req", "lock"], [], ["req"])
+
+    drs.add_init_state("0")
+    drs.add_state("1")
+    
+    all_out = [["out"] for i in range(k)]
+    all_out.append([]) # empty context for the controller
+
+    all_empty = [[] for i in range(k+1)]
+    
+    drs.add_transition("0", ([i for i in range(k+1)], all_out), "1")
+    for i in range(k):
+        drs.add_transition("1", ([i,ctrl], all_empty), "1")
+        one_busy = [[] for x in range(k+1)]
+        one_busy[i] = ["busy"]
+        drs.add_transition("1", ([i,ctrl], one_busy), "1")
+    
+    return drs
+    
+def drs_mutex_property1(k):
+    state = [["in"]]
+    state.extend([[] for i in range(k)])
+    
+    return state
+    

rssmt.py

@@ -0,0 +1,80 @@
+#!/usr/bin/env python
+# -*- coding: utf-8 -*-
+
+"""
+    Reaction Systems SMT-Based Model Checking Module
+
+"""
+
+from rctsys import ReactionSystem
+from smtchecker import SmtChecker
+from smtcheckerpgrs import SmtCheckerPGRS
+from smtcheckerdistribrs import SmtCheckerDistribRS
+import sys
+import rs_examples
+
+import resource
+
+profiling = False
+
+version = "0.002"
+rsmc_banner = """
+ *** Reaction Systems SMT-Based Model Checking
+
+ *** Version: """ + version + """
+ *** Author:  Artur Męski <meski@ipipan.waw.pl> / <artur.meski@ncl.ac.uk>
+"""
+
+
+def main():
+    
+    if len(sys.argv) < 1+1:
+        print("provide N")
+        exit(1)
+    
+    N=int(sys.argv[1])
+        
+    print(rsmc_banner)
+
+    #rs = rs_examples.toy_ex3()
+    #rs = rs_examples.bitctr(16)
+    
+    #smt = SmtChecker(rs)
+    #smt.check_reachability(["p0","p1","p2","p3","p4","p5","p6","p7","p8","p9","p10","p11","p12","p13","p14","p15"], print_time=True)
+    #smt.check_reachability(["1","3","4"], print_time=True)
+
+    # PGRS:
+    # rsca = rs_examples.ca_toy_ex1()
+    # rsca.show()
+    # smt = SmtCheckerPGRS(rsca)
+    # smt.check_reachability(rs_examples.ca_toy_ex1_property1(), print_time=True)
+
+    # Distributed RS:
+
+    drs = rs_examples.drs_mutex(N)
+    #drs.show()
+
+    smt = SmtCheckerDistribRS(drs,debug_level=0)
+    smt.check_reachability(rs_examples.drs_mutex_property1(N),
+        exclusive_state=False,
+        max_level=10,
+        print_time=True,
+        print_mem=True)
+
+    # drs = rs_examples.drs_toy_ex1()
+    # drs.show()
+    #
+    # smt = SmtCheckerDistribRS(drs,debug_level=3)
+    # smt.check_reachability(rs_examples.drs_toy_ex1_property1(),max_level=10,print_time=True)
+
+if __name__ == "__main__":
+    try:
+        if profiling:
+            import profile
+            profile.run('main()')
+        else:
+            main()
+    except KeyboardInterrupt:
+        print("\nQuitting...")
+        sys.exit(99)
+

run.sh

@@ -0,0 +1,7 @@
+#!/bin/sh
+
+for i in 2 3 4 5 6 7 8 9 10 15 20 25 30 35 40 45 50 55 60 65 70 75 80 85 90 95 100;
+do
+  ./rssmt.py $i | tee "out/result_$i.out" 
+  echo "$i: $(tail -1 out/result_$i.out)" >> out/all.out
+done

smtchecker.py

@@ -0,0 +1,530 @@
+"""
+SMT-based Model Checking Module for RS
+"""
+
+from z3 import *
+from time import time
+from sys import stdout
+
+class SmtChecker(object):
+
+    def __init__(self, rs):
+
+        ##############################################################
+        # Encoded RS
+        ##############################################################
+        rs.sanity_check()
+        self.reaction_system = rs
+
+        ##############################################################
+        # SMT variables
+        ##############################################################
+        self.v = []
+        self.v_init = []
+
+        #self.vSucc = []
+        #self.vSuccInit = []
+
+        self.v_ctx = []
+
+        self.next_level_to_encode = 0
+
+        ##############################################################
+        # SMT solver instance
+        ##############################################################
+        self.solver = Solver()
+
+    #def smtVar(self, level, entityID, primed=False):
+    #    return "?"
+
+    def prepare_context_variables(self):
+        """Encodes all the context variables"""
+
+        level = self.next_level_to_encode
+
+        variables = []
+        for entity in self.reaction_system.background_set:
+            variables.append(Bool("C"+str(level)+"_"+entity))
+
+        self.v_ctx.append(variables)
+
+    def prepare_state_variables(self):
+        """Encodes all the state variables (including successors)"""
+
+        level = self.next_level_to_encode
+
+        variables = []
+        #variablesSucc = []
+        for entity in self.reaction_system.background_set:
+            variables.append(Bool("L"+str(level)+"_"+entity))
+            #variablesSucc.append(Bool("R"+str(level)+"_"+entity))
+
+        self.v.append(variables)
+        #self.vSucc.append(variablesSucc)
+
+        self.v_init.append(Bool("L"+str(level)+"_Init"))
+        #self.vSuccInit.append(Bool("R"+str(level)+"_Init"))
+
+    def prepare_all_variables(self):
+        """Encodes all the variables"""
+
+        self.prepare_state_variables()
+        self.prepare_context_variables()
+        self.next_level_to_encode += 1
+
+    def enc_init_state(self, level):
+        """Encodes the initial state at the given level"""
+
+        init_state_enc = self.v_init[level] 
+
+        for v in self.v[level]:
+            init_state_enc = simplify(And(init_state_enc, Not(v)))
+
+        return init_state_enc
+
+    def enc_init_contexts(self, level):
+        """Encodes the initial contexts set at the given level"""
+
+        init_contexts_set_enc = False # Or
+
+        for ctx in self.reaction_system.init_contexts:
+            single_ctx_enc = True # And
+            
+            not_ctx_entities = list(range(0, len(self.reaction_system.background_set)))
+            for entity in ctx:
+                single_ctx_enc = simplify(And(single_ctx_enc, self.v_ctx[level][entity]))
+                not_ctx_entities.remove(entity)
+
+            for entity in not_ctx_entities:
+                single_ctx_enc = simplify(And(single_ctx_enc, Not(self.v_ctx[level][entity])))
+
+            init_contexts_set_enc = simplify(Or(init_contexts_set_enc, single_ctx_enc))
+
+        #print("initContextSetEnc: " + repr(initContextsSetEnc))
+        return init_contexts_set_enc
+
+    def enc_not_allowed_contexts(self, level):
+        """Encodes all the context entities that are not allowed in the context sets"""
+
+        bg_set = set(range(0,len(self.reaction_system.background_set)))
+        ctx_ent_set = set(self.reaction_system.context_entities)
+
+        not_ctx_ent_set = bg_set.difference(ctx_ent_set)
+
+        enc = True
+        for entity in not_ctx_ent_set:
+            enc = simplify(And(enc, Not(self.v_ctx[level][entity])))
+
+        return enc
+
+    def enc_enabledness(self, level, prod_entity):
+        """Encodes the enabledness condition for a given level and a given entity"""
+
+        rcts_for_prod_entity = self.reaction_system.get_reactions_by_product()[prod_entity]
+
+        if rcts_for_prod_entity == []:
+            return False
+
+        #encInitEnab = simplify(Or(And(self.vInit[level], self.encInitContexts(level)), 
+        #        And(Not(self.vInit[level]), self.encNotAllowedContexts(level))))
+
+        enc_rct_prod = False
+        for ri_pair in rcts_for_prod_entity: # reactants-inhibitors pair
+            enc_reactants = True
+            enc_inhibitors = True
+            for reactant in ri_pair[0]:
+                enc_reactants = simplify(And(enc_reactants, 
+                                            Or(self.v[level][reactant], self.v_ctx[level][reactant])))
+            for inhibitor in ri_pair[1]:
+                enc_inhibitors = simplify(And(enc_inhibitors, 
+                                             Not(Or(self.v[level][inhibitor], self.v_ctx[level][inhibitor]))))
+
+            enc_rct_prod = simplify(Or(enc_rct_prod, And(enc_reactants, enc_inhibitors)))
+
+        #print("encEnabledness(" + repr(prodEntity) + "): " + repr(simplify(And(encInitEnab, encRctProd))))
+        #return simplify(And(encInitEnab, encRctProd))
+        return enc_rct_prod
+
+    def enc_entity_production(self, level, prod_entity):
+        """Encodes the production of a given entity from a given level at level+1"""
+
+        enc_enab_cond = self.enc_enabledness(level, prod_entity)
+
+        enc_ent_prod = Or(And(enc_enab_cond, self.v[level+1][prod_entity]),
+                And(Not(enc_enab_cond), Not(self.v[level+1][prod_entity])))
+
+        return enc_ent_prod
+
+    def enc_transition_relation(self, level):
+        """Encodes the transition relation"""
+
+        unused_entities = list(range(0,len(self.reaction_system.background_set)))
+
+        enc_trans = True
+
+        for prod_entity in self.reaction_system.get_reactions_by_product():
+            unused_entities.remove(prod_entity)
+    
+            enc_trans = simplify(And(enc_trans, self.enc_entity_production(level, prod_entity)))
+
+        enc_trans = simplify(And(enc_trans, Not(self.v_init[level+1])))
+
+        for prod_entity in unused_entities:
+            enc_trans = simplify(And(enc_trans, Not(self.v[level+1][prod_entity])))
+
+        if level == 0:
+            enc_init_enab = self.enc_init_contexts(level)
+        else: # level > 0:
+            enc_init_enab = self.enc_not_allowed_contexts(level)
+
+        #encInitEnab = simplify(Or(And(self.vInit[level], self.encInitContexts(level)), 
+        #        And(Not(self.vInit[level]), self.encNotAllowedContexts(level))))
+
+        enc_trans = simplify(And(enc_trans, enc_init_enab))
+
+        return enc_trans
+
+    def enc_state(self, level, state):
+        """Encodes the state at the given level"""
+
+        enc = Not(self.v_init[level])
+
+        state_ids = self.reaction_system.get_state_ids(state)
+
+        for entity in state_ids:
+            enc = And(enc, self.v[level][entity])
+
+        not_in_state = set(range(0, len(self.reaction_system.background_set)))
+        not_in_state = not_in_state.difference(set(state_ids))
+
+        for entity in not_in_state:
+            enc = And(enc, Not(self.v[level][entity]))
+
+        return enc
+
+    def decode_witness(self, max_level):
+
+        m = self.solver.model()
+
+        for level in range(0,max_level+1):
+
+            print("\n[Level=" + repr(level) + "]")
+
+            if repr(m[self.v_init[level]]) == "True":
+                print("** Initial state")
+
+            print("State:\n{"),
+            for var_id in range(0, len(self.v[level])):
+                if repr(m[self.v[level][var_id]]) == "True":
+                    print("\t" + self.reaction_system.get_entity_name(var_id)),
+            print("}")
+
+            if level != max_level:
+                print("Context set:"),
+                print("{"),
+                for var_id in range(0, len(self.v[level])):
+                    if repr(m[self.v_ctx[level][var_id]]) == "True":
+                        print("\t" + self.reaction_system.get_entity_name(var_id)),
+                print("}")
+
+
+    def check_reachability(self, state, print_witness=True, print_time=False):
+        """Main testing function"""
+
+        if print_time:
+            start = time()
+
+        self.prepare_all_variables()
+        self.solver.add(self.enc_init_state(0))
+        current_level = 0
+
+        while True:
+            #print("Level: " + str(current_level))
+            print("\rLevel: " + str(current_level)),
+            stdout.flush()
+
+            self.prepare_all_variables()
+
+            # reachability test:
+            self.solver.push()
+            self.solver.add(self.enc_state(current_level,state))
+
+            result = self.solver.check()
+            print(result)
+            if result == sat:
+                print("\nSAT")
+                if print_witness:
+                    self.decode_witness(current_level)
+                break
+            else:
+                self.solver.pop()
+           
+            self.solver.add(self.enc_transition_relation(current_level))
+
+            current_level += 1
+
+        if print_time:
+            stop = time()
+            print("Time: " + repr(stop-start))
+
+
+class SmtCheckerPGRS(object):
+
+    def __init__(self, rs):
+
+        ##############################################################
+        # Encoded RS
+        ##############################################################
+        rs.sanity_check()
+        self.reaction_system = rs
+
+        ##############################################################
+        # SMT variables
+        ##############################################################
+        self.v = []
+        self.v_init = []
+
+        #self.vSucc = []
+        #self.vSuccInit = []
+
+        self.v_ctx = []
+
+        self.next_level_to_encode = 0
+
+        ##############################################################
+        # SMT solver instance
+        ##############################################################
+        self.solver = Solver()
+
+    #def smtVar(self, level, entityID, primed=False):
+    #    return "?"
+
+    def prepare_context_variables(self):
+        """Encodes all the context variables"""
+
+        level = self.next_level_to_encode
+
+        variables = []
+        for entity in self.reaction_system.background_set:
+            variables.append(Bool("C"+str(level)+"_"+entity))
+
+        self.v_ctx.append(variables)
+
+    def prepare_state_variables(self):
+        """Encodes all the state variables (including successors)"""
+
+        level = self.next_level_to_encode
+
+        variables = []
+        #variablesSucc = []
+        for entity in self.reaction_system.background_set:
+            variables.append(Bool("L"+str(level)+"_"+entity))
+            #variablesSucc.append(Bool("R"+str(level)+"_"+entity))
+
+        self.v.append(variables)
+        #self.vSucc.append(variablesSucc)
+
+        self.v_init.append(Bool("L"+str(level)+"_Init"))
+        #self.vSuccInit.append(Bool("R"+str(level)+"_Init"))
+
+    def prepare_all_variables(self):
+        """Encodes all the variables"""
+
+        self.prepare_state_variables()
+        self.prepare_context_variables()
+        self.next_level_to_encode += 1
+
+    def enc_init_state(self, level):
+        """Encodes the initial state at the given level"""
+
+        init_state_enc = self.v_init[level] 
+
+        for v in self.v[level]:
+            init_state_enc = simplify(And(init_state_enc, Not(v)))
+
+        return init_state_enc
+
+    def enc_init_contexts(self, level):
+        """Encodes the initial contexts set at the given level"""
+
+        init_contexts_set_enc = False # Or
+
+        for ctx in self.reaction_system.init_contexts:
+            single_ctx_enc = True # And
+            
+            not_ctx_entities = list(range(0, len(self.reaction_system.background_set)))
+            for entity in ctx:
+                single_ctx_enc = simplify(And(single_ctx_enc, self.v_ctx[level][entity]))
+                not_ctx_entities.remove(entity)
+
+            for entity in not_ctx_entities:
+                single_ctx_enc = simplify(And(single_ctx_enc, Not(self.v_ctx[level][entity])))
+
+            init_contexts_set_enc = simplify(Or(init_contexts_set_enc, single_ctx_enc))
+
+        #print("initContextSetEnc: " + repr(initContextsSetEnc))
+        return init_contexts_set_enc
+
+    def enc_not_allowed_contexts(self, level):
+        """Encodes all the context entities that are not allowed in the context sets"""
+
+        bg_set = set(range(0,len(self.reaction_system.background_set)))
+        ctx_ent_set = set(self.reaction_system.context_entities)
+
+        not_ctx_ent_set = bg_set.difference(ctx_ent_set)
+
+        enc = True
+        for entity in not_ctx_ent_set:
+            enc = simplify(And(enc, Not(self.v_ctx[level][entity])))
+
+        return enc
+
+    def enc_enabledness(self, level, prod_entity):
+        """Encodes the enabledness condition for a given level and a given entity"""
+
+        rcts_for_prod_entity = self.reaction_system.get_reactions_by_product()[prod_entity]
+
+        if rcts_for_prod_entity == []:
+            return False
+
+        #encInitEnab = simplify(Or(And(self.vInit[level], self.encInitContexts(level)), 
+        #        And(Not(self.vInit[level]), self.encNotAllowedContexts(level))))
+
+        enc_rct_prod = False
+        for ri_pair in rcts_for_prod_entity: # reactants-inhibitors pair
+            enc_reactants = True
+            enc_inhibitors = True
+            for reactant in ri_pair[0]:
+                enc_reactants = simplify(And(enc_reactants, 
+                                            Or(self.v[level][reactant], self.v_ctx[level][reactant])))
+            for inhibitor in ri_pair[1]:
+                enc_inhibitors = simplify(And(enc_inhibitors, 
+                                             Not(Or(self.v[level][inhibitor], self.v_ctx[level][inhibitor]))))
+
+            enc_rct_prod = simplify(Or(enc_rct_prod, And(enc_reactants, enc_inhibitors)))
+
+        #print("encEnabledness(" + repr(prodEntity) + "): " + repr(simplify(And(encInitEnab, encRctProd))))
+        #return simplify(And(encInitEnab, encRctProd))
+        return enc_rct_prod
+
+    def enc_entity_production(self, level, prod_entity):
+        """Encodes the production of a given entity from a given level at level+1"""
+
+        enc_enab_cond = self.enc_enabledness(level, prod_entity)
+
+        enc_ent_prod = Or(And(enc_enab_cond, self.v[level+1][prod_entity]),
+                And(Not(enc_enab_cond), Not(self.v[level+1][prod_entity])))
+
+        return enc_ent_prod
+
+    def enc_transition_relation(self, level):
+        """Encodes the transition relation"""
+
+        unused_entities = list(range(0,len(self.reaction_system.background_set)))
+
+        enc_trans = True
+
+        for prod_entity in self.reaction_system.get_reactions_by_product():
+            unused_entities.remove(prod_entity)
+    
+            enc_trans = simplify(And(enc_trans, self.enc_entity_production(level, prod_entity)))
+
+        enc_trans = simplify(And(enc_trans, Not(self.v_init[level+1])))
+
+        for prod_entity in unused_entities:
+            enc_trans = simplify(And(enc_trans, Not(self.v[level+1][prod_entity])))
+
+        if level == 0:
+            enc_init_enab = self.enc_init_contexts(level)
+        else: # level > 0:
+            enc_init_enab = self.enc_not_allowed_contexts(level)
+
+        #encInitEnab = simplify(Or(And(self.vInit[level], self.encInitContexts(level)), 
+        #        And(Not(self.vInit[level]), self.encNotAllowedContexts(level))))
+
+        enc_trans = simplify(And(enc_trans, enc_init_enab))
+
+        return enc_trans
+
+    def enc_state(self, level, state):
+        """Encodes the state at the given level"""
+
+        enc = Not(self.v_init[level])
+
+        state_ids = self.reaction_system.get_state_ids(state)
+
+        for entity in state_ids:
+            enc = And(enc, self.v[level][entity])
+
+        not_in_state = set(range(0, len(self.reaction_system.background_set)))
+        not_in_state = not_in_state.difference(set(state_ids))
+
+        for entity in not_in_state:
+            enc = And(enc, Not(self.v[level][entity]))
+
+        return enc
+
+    def decode_witness(self, max_level):
+
+        m = self.solver.model()
+
+        for level in range(0,max_level+1):
+
+            print("\n[Level=" + repr(level) + "]")
+
+            if repr(m[self.v_init[level]]) == "True":
+                print("** Initial state")
+
+            print("State:\n{"),
+            for var_id in range(0, len(self.v[level])):
+                if repr(m[self.v[level][var_id]]) == "True":
+                    print("\t" + self.reaction_system.get_entity_name(var_id)),
+            print("}")
+
+            if level != max_level:
+                print("Context set:"),
+                print("{"),
+                for var_id in range(0, len(self.v[level])):
+                    if repr(m[self.v_ctx[level][var_id]]) == "True":
+                        print("\t" + self.reaction_system.get_entity_name(var_id)),
+                print("}")
+
+
+    def check_reachability(self, state, print_witness=True, print_time=False):
+        """Main testing function"""
+
+        if print_time:
+            start = time()
+
+        self.prepare_all_variables()
+        self.solver.add(self.enc_init_state(0))
+        current_level = 0
+
+        while True:
+            #print("Level: " + str(current_level))
+            print("\rLevel: " + str(current_level)),
+            stdout.flush()
+
+            self.prepare_all_variables()
+
+            # reachability test:
+            self.solver.push()
+            self.solver.add(self.enc_state(current_level,state))
+
+            result = self.solver.check()
+            print(result)
+            if result == sat:
+                print("\nSAT")
+                if print_witness:
+                    self.decode_witness(current_level)
+                break
+            else:
+                self.solver.pop()
+           
+            self.solver.add(self.enc_transition_relation(current_level))
+
+            current_level += 1
+
+        if print_time:
+            stop = time()
+            print("Time: " + repr(stop-start))

smtcheckerdistribrs.py

@@ -0,0 +1,432 @@
+"""
+SMT-based Model Checking Module for RS with Context Automaton
+"""
+
+from z3 import *
+from time import time,sleep
+from sys import stdout
+import resource
+
+# def simplify(x):
+#     return x
+
+class SmtCheckerDistribRS(object):
+
+    def __init__(self, drs, debug_level=1):
+
+        print("[i] Initialising the SMT module")
+
+        drs.sanity_check()
+        
+        self.solver = Solver()
+        
+        self.debug_level = debug_level
+        self.drs = drs
+        self.n_components = self.drs.components_count
+
+        # encoding variables
+        self.v = [] # level -> component -> variables
+        self.v_ctx = []
+        self.v_act = [] # indicators of which component is active
+        self.ca_state = []
+
+        self.level_to_encode = 0
+        
+    def prepare_all_variables(self):
+        """Encodes the required variables"""
+
+        self.prepare_state_variables()
+        self.prepare_context_variables()
+        self.prepare_activity_variables()
+        
+        self.level_to_encode += 1 # prepare for the next invocation
+        
+    def prepare_context_variables(self):
+        """Encodes the context variables"""
+
+        level = self.level_to_encode
+
+        if self.debug_level > 1:
+            print("[ii] Preparing context variables for level=" + str(level))
+            
+        level_variables = []
+        
+        for i in range(self.n_components):
+
+            comp_variables = []
+            
+            for entity in self.drs.background_set:
+                comp_variables.append(Bool("Ctx"+str(level)+"V"+str(i)+"_"+entity))
+            
+            level_variables.append(comp_variables)
+
+        self.v_ctx.append(level_variables)
+        
+    def prepare_activity_variables(self):
+        """Encodes the activity variables"""
+
+        level = self.level_to_encode
+        
+        if self.debug_level > 1:
+            print("[ii] Preparing activity variables for level=" + str(level))
+
+        level_variables = []
+        
+        for i in range(self.n_components):
+            # L - level, A - activity indicator
+            level_variables.append(Bool("L"+str(level)+"A"+str(i)))
+            
+        self.v_act.append(level_variables)
+
+    def prepare_state_variables(self):
+        """Encodes all the state variables"""
+
+        level = self.level_to_encode
+
+        if self.debug_level > 1:
+            print("[ii] Preparing state variables for level=" + str(level))
+
+        level_variables = [] # level vars
+        
+        for i in range(self.n_components):
+
+            comp_variables = []
+            
+            for entity in self.drs.background_set:
+                # L - level, V - component
+                comp_variables.append(Bool("L"+str(level)+"V"+str(i)+"_"+entity)) 
+            
+            level_variables.append(comp_variables)
+
+        self.v.append(level_variables)
+
+        # single state variable for CA
+        self.ca_state.append(Int("CA"+str(level)+"_state"))
+
+    def enc_init_state(self, level):
+        """Encodes the initial state at the given level"""
+
+        if self.debug_level > 1:
+            print("[ii] Encoding the initial state for level=" + str(level))
+
+        rs_init_state_enc = True
+
+        for i in range(self.n_components):
+            for v in self.v[level][i]:
+                rs_init_state_enc = simplify(And(rs_init_state_enc, Not(v))) # the initial state is empty
+
+        ca_init_state_enc = self.ca_state[level] == self.drs.get_init_state_id()
+        init_state_enc = simplify(And(rs_init_state_enc, ca_init_state_enc))
+
+        # print ("init_state_enc:\n", init_state_enc)
+
+        return init_state_enc
+
+    def enc_enabledness(self, level, prod_entity, component_id):
+        """Encodes the enabledness condition for a given level and a given entity"""
+
+        rcts_for_prod_entity = self.drs.get_reactions_by_product(component_id)[prod_entity]
+
+        if rcts_for_prod_entity == []:
+            return False
+
+        enc_rct_prod = False
+        for reactants,inhibitors in rcts_for_prod_entity:
+
+    
+            enc_reactants = True
+            for reactant in reactants:
+                
+                enc_active_reactants = False
+                for i in range(self.n_components):
+                    enc_active_reactants = simplify(Or(enc_active_reactants, 
+                            And(  Or(self.v[level][i][reactant], self.v_ctx[level][i][reactant]),  self.v_act[level][i]  )))
+                
+                enc_reactants = And(enc_reactants, enc_active_reactants)
+            
+    
+            enc_inhibitors = True
+            for inhibitor in inhibitors:
+            
+                enc_active_inhibitors = True
+                for i in range(self.n_components):
+                    enc_active_inhibitors = simplify(And(enc_active_inhibitors, 
+                            And(
+                                Or( And(Not(self.v[level][i][inhibitor]), Not(self.v_ctx[level][i][inhibitor])), Not(self.v_act[level][i]) )
+                            )))
+
+                enc_inhibitors = simplify(And(enc_inhibitors, enc_active_inhibitors))
+
+            # print("--> enc_inhibitors\n", enc_inhibitors)
+            enc_rct_prod = Or(enc_rct_prod, And(enc_reactants, enc_inhibitors))
+
+        # print("enc_rct_prod:\n", enc_rct_prod)
+
+        enc_rct_prod = simplify(enc_rct_prod)
+
+        return enc_rct_prod
+
+    def enc_entity_production(self, level, prod_entity, component_id):
+        """Encodes the production of a given entity at level+1 from a given level"""
+
+        enc_enab_cond = self.enc_enabledness(level, prod_entity, component_id)
+
+        enc_base_ent_prod = simplify(Or(And(enc_enab_cond, self.v[level+1][component_id][prod_entity]),
+                And(Not(enc_enab_cond), Not(self.v[level+1][component_id][prod_entity]))))
+
+        enc_active_ent_prod = simplify(And(self.v_act[level][component_id], enc_base_ent_prod))
+        enc_inactive_ent_prod = simplify(And(
+                Not(self.v_act[level][component_id]), 
+                self.v[level][component_id][prod_entity] == self.v[level+1][component_id][prod_entity]))
+        
+        enc_ent_prod = Or(enc_active_ent_prod, enc_inactive_ent_prod)
+
+        # print("enc_ent_prod:\n", enc_ent_prod)
+
+        return simplify(enc_ent_prod)
+        
+    def enc_rs_trans(self, level):
+        """Encodes the transition relation"""
+
+        enc_trans = True
+
+        for component_id in range(self.n_components):
+
+            print("\rEncoding for reactions: %d/%d" % (component_id,self.n_components-1), flush=True, end="")
+                
+            unused_entities = list(range(len(self.drs.background_set)))
+
+            for prod_entity in self.drs.get_reactions_by_product(component_id):
+                unused_entities.remove(prod_entity)
+    
+                enc_trans = simplify(And(enc_trans, self.enc_entity_production(level, prod_entity, component_id)))
+
+            for prod_entity in unused_entities:
+                enc_trans = simplify(And(enc_trans, Not(self.v[level+1][component_id][prod_entity])))
+
+        print()
+        # print("enc_rs_trans:\n", enc_trans)
+        
+        enc_trans = simplify(enc_trans)
+        
+        return enc_trans
+    
+    def enc_automaton_trans(self, level):
+        """Encodes the transition relation for the context automaton"""
+        
+        enc_trans = False
+    
+        i = 0
+        for src,(components,ctx_set),dst in self.drs.transitions:
+                src_enc = self.ca_state[level] == src
+                dst_enc = self.ca_state[level+1] == dst
+                
+                print("\rEncoding for context automaton: %d/%d" % (i,len(self.drs.transitions)-1), flush=True, end="")
+                
+                i = i + 1
+                
+                # contexts {
+                ctx_set_enc = True
+                for comp_id in range(self.n_components):
+
+                    all_ent = set(range(len(self.drs.background_set)))
+                    incl_ctx = ctx_set[comp_id]
+                    excl_ctx = all_ent - incl_ctx
+                
+                    ctx_enc = True
+                
+                    for c in incl_ctx:
+                        ctx_enc = And(ctx_enc, self.v_ctx[level][comp_id][c])
+                    for c in excl_ctx:
+                        ctx_enc = And(ctx_enc, Not(self.v_ctx[level][comp_id][c]))
+                        
+                    ctx_set_enc = And(ctx_set_enc, ctx_enc)
+                # } contexts
+                
+                # active components {
+                all_active = set(range(self.n_components))
+                incl_comp = components
+                excl_comp = all_active - incl_comp
+                       
+                active_components_enc = True
+                for comp_id in incl_comp:
+                    active_components_enc = And(active_components_enc, self.v_act[level][comp_id])
+                for comp_id in excl_comp:
+                    active_components_enc = And(active_components_enc, Not(self.v_act[level][comp_id]))
+                # } active components
+                
+                cur_trans = And(src_enc, ctx_set_enc, active_components_enc, dst_enc)
+                
+                enc_trans = Or(enc_trans, cur_trans)
+        
+        print()
+        
+        enc_trans = simplify(enc_trans)
+        
+        # print("enc_automaton_trans:\n", enc_trans)
+        return enc_trans
+
+    def enc_transition_relation(self, level):
+        
+        rs_enc = self.enc_rs_trans(level)
+        aut_enc = self.enc_automaton_trans(level)
+        
+        print("Conjunction...", flush=True, end="")
+        
+        c = simplify(And(rs_enc, aut_enc))
+        
+        print("done.")
+        
+        return c
+        
+    def enc_state(self, level, global_state):
+        """Encodes the state at the given level"""
+
+        if len(global_state) != self.n_components:
+            print("EEE: Wrong size of the global state! " + "(is " + str(len(global_state)) + ", should be " + str(self.n_components) + ")")
+            exit(1)
+
+        enc = True
+
+        if self.debug_level > 2:
+            print("[iii] Encoding exclusive/exact global state " + str(global_state) + " for level=" + str(level))
+
+        for i in range(self.n_components):
+            local_state = global_state[i]
+            
+            local_state_ids = self.drs.get_state_ids(local_state)
+
+            for entity in local_state_ids:
+                enc = And(enc, self.v[level][i][entity])
+
+            not_in_state = self.drs.set_of_background_ids - set(local_state_ids)
+
+            for e_id in not_in_state:
+                enc = simplify(And(enc, Not(self.v[level][i][e_id])))
+        
+        simplify(enc)
+        # print("state:\n", enc)
+        return enc
+
+    def enc_inclusive_state(self, level, global_state):
+        """Encodes the state at the given level"""
+
+        if len(global_state) != self.n_components:
+            print("EEE: Wrong size of the global state! " + "(is " + str(len(global_state)) + ", should be " + str(self.n_components) + ")")
+            exit(1)
+
+        enc = True
+
+        if self.debug_level > 2:
+            print("[iii] Encoding inclusive/general global state " + str(global_state) + " for level=" + str(level))
+
+        for i in range(self.n_components):
+            local_state = global_state[i]
+            
+            local_state_ids = self.drs.get_state_ids(local_state)
+
+            for entity in local_state_ids:
+                enc = And(enc, self.v[level][i][entity])
+        
+        simplify(enc)
+        return enc
+
+    def decode_witness(self, max_level, print_model=False):
+
+        m = self.solver.model()
+
+        print("\nWitness:")
+
+        if print_model:
+            print(m)
+
+        for level in range(max_level+1):
+
+            print("\n[Level=" + repr(level) + "]")
+
+            #print(m)
+            #print(self.v[level][0][2])
+            #print(m[self.v[level][0][2]])
+
+            print("  State: {", end=""),
+            for c in range(self.n_components):
+                print(" <", end="")
+                for var_id in range(len(self.v[level][c])):
+                    if repr(m[self.v[level][c][var_id]]) == "True":
+                        print(" " + self.drs.get_entity_name(var_id), end="")
+                print(" >", end="")
+            print(" }")
+
+            if level != max_level:
+                print("  Context set: ", end="")
+                print("{", end="")
+                for c in range(self.n_components):
+                    print(" <", end="")
+                    for var_id in range(len(self.v[level][c])):
+                        if repr(m[self.v_ctx[level][c][var_id]]) == "True":
+                            print(" " + self.drs.get_entity_name(var_id), end="")
+                    print(" >", end="")
+                print(" }")
+            
+                print("  Active components:", end="")
+                for c in range(self.n_components):
+                    if repr(m[self.v_act[level][c]]) == "True":
+                        print(" " + str(c), end="")
+                print()
+            
+
+
+    def check_reachability(self, state, exclusive_state=False, print_witness=True, print_time=False, print_mem=False, max_level=100):
+        """Reachability checking"""
+
+        print("[i] Checking reachability...")
+
+        if print_time:
+            start = time()
+
+        self.prepare_all_variables()
+        self.solver.add(self.enc_init_state(0))
+        current_level = 0
+
+        while True:
+
+            self.prepare_all_variables()
+                   
+            print("-----[ Working at level=" + str(current_level) + " ]-----")
+            stdout.flush()
+       
+             # reachability test:
+            self.solver.push()
+            print("[i] Adding the reachability test...")
+            if exclusive_state == False:
+                self.solver.add(self.enc_inclusive_state(current_level,state))
+            else:
+                self.solver.add(self.enc_state(current_level,state))
+        
+            result = self.solver.check()
+            if result == sat:
+                print("\n[+] SAT at level=" + str(current_level))
+                if print_witness:
+                    self.decode_witness(current_level)
+                break
+            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 current_level > max_level:
+                print("Stopping at level=" + str(max_level))
+                break
+
+        if print_time:
+            stop = time()
+            print()
+            print("==== Time: " + repr(stop-start))
+
+        if print_mem:
+            usage=resource.getrusage(resource.RUSAGE_SELF)
+            print("MEM: usertime=%s systime=%s mem=%sMB" % (usage[0],usage[1], (resource.getrusage(resource.RUSAGE_SELF).ru_maxrss / 1024**2)))

smtcheckerpgrs.py

@@ -0,0 +1,240 @@
+"""
+SMT-based Model Checking Module for RS with Context Automaton
+"""
+
+from z3 import *
+from time import time
+from sys import stdout
+
+class SmtCheckerPGRS(object):
+
+    def __init__(self, rsca):
+
+        rsca.sanity_check()
+        
+        self.rs = rsca.rs
+        self.ca = rsca.ca        
+
+        self.v = []
+        self.v_ctx = []
+        self.ca_state = []
+        self.next_level_to_encode = 0
+
+        self.solver = Solver()
+
+    def prepare_all_variables(self):
+        """Encodes all the variables"""
+
+        self.prepare_state_variables()
+        self.prepare_context_variables()
+        self.next_level_to_encode += 1
+        
+    def prepare_context_variables(self):
+        """Encodes all the context variables"""
+
+        level = self.next_level_to_encode
+
+        variables = []
+        for entity in self.rs.background_set:
+            variables.append(Bool("C"+str(level)+"_"+entity))
+
+        self.v_ctx.append(variables)
+
+    def prepare_state_variables(self):
+        """Encodes all the state variables"""
+
+        level = self.next_level_to_encode
+
+        variables = []
+        for entity in self.rs.background_set:
+            variables.append(Bool("L"+str(level)+"_"+entity))
+        self.v.append(variables)
+        
+        self.ca_state.append(Int("CA"+str(level)+"_state"))
+
+    def enc_init_state(self, level):
+        """Encodes the initial state at the given level"""
+
+        rs_init_state_enc = True
+
+        for v in self.v[level]:
+            rs_init_state_enc = simplify(And(rs_init_state_enc, Not(v))) # the initial state is empty
+
+        ca_init_state_enc = self.ca_state[level] == self.ca.get_init_state_id()
+        
+        init_state_enc = simplify(And(rs_init_state_enc, ca_init_state_enc))
+
+        return init_state_enc
+
+    def enc_enabledness(self, level, prod_entity):
+        """Encodes the enabledness condition for a given level and a given entity"""
+
+        rcts_for_prod_entity = self.rs.get_reactions_by_product()[prod_entity]
+
+        if rcts_for_prod_entity == []:
+            return False
+
+        enc_rct_prod = False
+        for reactants,inhibitors in rcts_for_prod_entity:
+            enc_reactants = True
+            enc_inhibitors = True
+            for reactant in reactants:
+                enc_reactants = simplify(And(enc_reactants, 
+                                            Or(self.v[level][reactant], self.v_ctx[level][reactant])))
+            for inhibitor in inhibitors:
+                enc_inhibitors = simplify(And(enc_inhibitors, 
+                                             Not(Or(self.v[level][inhibitor], self.v_ctx[level][inhibitor]))))
+
+            enc_rct_prod = simplify(Or(enc_rct_prod, And(enc_reactants, enc_inhibitors)))
+
+        return enc_rct_prod
+
+    def enc_entity_production(self, level, prod_entity):
+        """Encodes the production of a given entity from a given level at level+1"""
+
+        enc_enab_cond = self.enc_enabledness(level, prod_entity)
+
+        enc_ent_prod = Or(And(enc_enab_cond, self.v[level+1][prod_entity]),
+                And(Not(enc_enab_cond), Not(self.v[level+1][prod_entity])))
+
+        return simplify(enc_ent_prod)
+        
+    def enc_transition_relation(self, level):
+        return simplify(And(self.enc_rs_trans(level), self.enc_automaton_trans(level)))
+
+    def enc_rs_trans(self, level):
+        """Encodes the transition relation"""
+
+        unused_entities = list(range(len(self.rs.background_set)))
+
+        enc_trans = True
+
+        for prod_entity in self.rs.get_reactions_by_product():
+            unused_entities.remove(prod_entity)
+    
+            enc_trans = simplify(And(enc_trans, self.enc_entity_production(level, prod_entity)))
+
+        for prod_entity in unused_entities:
+            enc_trans = simplify(And(enc_trans, Not(self.v[level+1][prod_entity])))
+
+        return enc_trans
+    
+    def enc_automaton_trans(self, level):
+        """Encodes the transition relation for the context automaton"""
+        
+        enc_trans = False
+    
+        for src,ctx,dst in self.ca.transitions:
+                src_enc = self.ca_state[level] == src
+                dst_enc = self.ca_state[level+1] == dst
+                
+                all_ent = set(range(len(self.rs.background_set)))
+                incl_ctx = ctx
+                excl_ctx = all_ent - incl_ctx
+                
+                ctx_enc = True
+                
+                for c in incl_ctx:
+                    ctx_enc = simplify(And(ctx_enc, self.v_ctx[level][c]))
+                for c in excl_ctx:
+                    ctx_enc = simplify(And(ctx_enc, Not(self.v_ctx[level][c])))
+                
+                cur_trans = simplify(And(src_enc, ctx_enc, dst_enc))
+                
+                enc_trans = simplify(Or(enc_trans, cur_trans))
+        
+        return enc_trans
+
+    def enc_state(self, level, state):
+        """Encodes the state at the given level"""
+
+        enc = True
+
+        state_ids = self.rs.get_state_ids(state)
+
+        for entity in state_ids:
+            enc = And(enc, self.v[level][entity])
+
+        not_in_state = set(range(len(self.rs.background_set)))
+        not_in_state = not_in_state.difference(set(state_ids))
+
+        for entity in not_in_state:
+            enc = And(enc, Not(self.v[level][entity]))
+
+        return enc
+
+    def enc_non_exclusive_state(self, level, state):
+        """Encodes the state at the given level"""
+
+        enc = True
+
+        state_ids = self.rs.get_state_ids(state)
+
+        for entity in state_ids:
+            enc = And(enc, self.v[level][entity])
+
+        return enc
+
+    def decode_witness(self, max_level, print_model=False):
+
+        m = self.solver.model()
+
+        if print_model:
+            print(m)
+
+        for level in range(max_level+1):
+
+            print("\n[Level=" + repr(level) + "]")
+
+            print("  State: {", end=""),
+            for var_id in range(len(self.v[level])):
+                if repr(m[self.v[level][var_id]]) == "True":
+                    print(" " + self.rs.get_entity_name(var_id), end="")
+            print(" }")
+
+            if level != max_level:
+                print("  Context set: ", end="")
+                print("{", end="")
+                for var_id in range(len(self.v[level])):
+                    if repr(m[self.v_ctx[level][var_id]]) == "True":
+                        print(" " + self.rs.get_entity_name(var_id), end="")
+                print(" }")
+
+    def check_reachability(self, state, print_witness=True, print_time=False):
+        """Main testing function"""
+
+        if print_time:
+            start = time()
+
+        self.prepare_all_variables()
+        self.solver.add(self.enc_init_state(0))
+        current_level = 0
+
+        while True:
+            print("\r[i] Level: " + str(current_level), end="")
+            stdout.flush()
+
+            self.prepare_all_variables()
+
+            # reachability test:
+            self.solver.push()
+            self.solver.add(self.enc_state(current_level,state))
+
+            result = self.solver.check()
+            if result == sat:
+                print("\n[+] SAT at level=" + str(current_level))
+                if print_witness:
+                    self.decode_witness(current_level)
+                break
+            else:
+                self.solver.pop()
+           
+            self.solver.add(self.enc_transition_relation(current_level))
+
+            current_level += 1
+
+        if print_time:
+            stop = time()
+            print()
+            print("[i] Time: " + repr(stop-start))
+            

testca.py

@@ -0,0 +1,16 @@
+#!/usr/bin/env python
+
+from rctsys import *
+import rs_examples
+
+rs = rs_examples.toy_ex3()
+
+ca = ContextAutomaton(rs)
+
+ca.add_state("st1")
+ca.add_init_state("in")
+
+ca.add_transition("st1", {"1","2"}, "in")
+
+ca.show()
+