diff --git a/paynt/cli.py b/paynt/cli.py
index 0b401075..a9c3a5e8 100644
--- a/paynt/cli.py
+++ b/paynt/cli.py
@@ -1,12 +1,13 @@
 from . import version
 
-import paynt.utils.profiler
+import paynt.utils.timer
 import paynt.parser.sketch
 
 import paynt.quotient
 import paynt.quotient.pomdp
 import paynt.quotient.decpomdp
 import paynt.quotient.storm_pomdp_control
+import paynt.quotient.mdp
 
 import paynt.synthesizer.all_in_one
 import paynt.synthesizer.synthesizer
@@ -133,6 +134,8 @@ def setup_logger(log_path = None):
     help="decision tree synthesis: tree depth")
 @click.option("--tree-enumeration", is_flag=True, default=False,
     help="decision tree synthesis: if set, all trees of size at most tree_depth will be enumerated")
+@click.option("--add-dont-care-action", is_flag=True, default=False,
+    help="decision tree synthesis: if set, an explicit action simulating a random action selection will be added to each state")
 
 @click.option(
     "--constraint-bound", type=click.FLOAT, help="bound for creating constrained POMDP for Cassandra models",
@@ -156,7 +159,7 @@ def paynt_run(
     export_fsc_storm, export_fsc_paynt, export_evaluation,
     all_in_one, all_in_one_maxmem,
     mdp_split_wrt_mdp, mdp_discard_unreachable_choices, mdp_use_randomized_abstraction,
-    tree_depth, tree_enumeration,
+    tree_depth, tree_enumeration, add_dont_care_action,
     constraint_bound,
     ce_generator,
     profiling
@@ -166,7 +169,7 @@ def paynt_run(
     if profiling:
         profiler = cProfile.Profile()
         profiler.enable()
-    paynt.utils.profiler.GlobalTimeLimit.start(timeout)
+    paynt.utils.timer.GlobalTimer.start(timeout)
 
     logger.info("This is Paynt version {}.".format(version()))
 
@@ -183,6 +186,7 @@ def paynt_run(
 
     paynt.synthesizer.decision_tree.SynthesizerDecisionTree.tree_depth = tree_depth
     paynt.synthesizer.decision_tree.SynthesizerDecisionTree.tree_enumeration = tree_enumeration
+    paynt.quotient.mdp.MdpQuotient.add_dont_care_action = add_dont_care_action
 
     storm_control = None
     if storm_pomdp:
diff --git a/paynt/quotient/mdp.py b/paynt/quotient/mdp.py
index 4bf2c022..38da49b4 100644
--- a/paynt/quotient/mdp.py
+++ b/paynt/quotient/mdp.py
@@ -70,10 +70,11 @@ class DecisionTreeNode:
     def __init__(self, parent):
         self.parent = parent
         self.variable_index = None
-        self.child_false = None
         self.child_true = None
-        self.hole = None
+        self.child_false = None
         self.identifier = None
+        self.holes = None
+        self.hole_assignment = None
 
     @property
     def is_terminal(self):
@@ -87,13 +88,10 @@ def child_nodes(self):
     def is_true_child(self):
         return self is self.parent.child_true
 
-    def add_decision(self):
-        '''
-        Associate (an index of) a variable with the node.
-        '''
-        if self.is_terminal:
-            self.child_false = DecisionTreeNode(self)
-            self.child_true = DecisionTreeNode(self)
+    def add_children(self):
+        assert self.is_terminal
+        self.child_true = DecisionTreeNode(self)
+        self.child_false = DecisionTreeNode(self)
 
     def assign_identifiers(self, identifier=0):
         self.identifier = identifier
@@ -103,6 +101,29 @@ def assign_identifiers(self, identifier=0):
         identifier = self.child_false.assign_identifiers(identifier+1)
         return identifier
 
+    def associate_holes(self, node_hole_info):
+        self.holes = [hole for hole,_,_ in node_hole_info[self.identifier]]
+        if self.is_terminal:
+            return
+        self.child_true.associate_holes(node_hole_info)
+        self.child_false.associate_holes(node_hole_info)
+
+    def associate_assignment(self, assignment):
+        self.hole_assignment = [assignment.hole_options(hole)[0] for hole in self.holes]
+        if self.is_terminal:
+            return
+        self.child_true.associate_assignment(assignment)
+        self.child_false.associate_assignment(assignment)
+
+    def apply_hint(self, subfamily, tree_hint):
+        if self.is_terminal or tree_hint.is_terminal:
+            return
+        for hole_index,option in enumerate(tree_hint.hole_assignment):
+            hole = self.holes[hole_index]
+            subfamily.hole_set_options(hole,[option])
+        self.child_true.apply_hint(subfamily,tree_hint.child_true)
+        self.child_false.apply_hint(subfamily,tree_hint.child_false)
+
 
 class DecisionTree:
 
@@ -118,7 +139,8 @@ def set_depth(self, depth:int):
         self.reset()
         for level in range(depth):
             for node in self.collect_terminals():
-                node.add_decision()
+                node.add_children()
+        self.root.assign_identifiers()
 
     def collect_nodes(self, node_condition=None):
         if node_condition is None:
@@ -139,7 +161,6 @@ def collect_nonterminals(self):
         return self.collect_nodes(lambda node : not node.is_terminal)
 
     def to_list(self):
-        self.root.assign_identifiers()
         num_nodes = len(self.collect_nodes())
         node_info = [ None for node in range(num_nodes) ]
         for node in self.collect_nodes():
@@ -149,22 +170,18 @@ def to_list(self):
             node_info[node.identifier] = (parent,child_true,child_false)
         return node_info
 
-    def collect_variables(self):
-        node_to_variable = None
-        nodes = self.collect_nodes()
-        node_to_variable = [len(self.variables) for node in range(len(nodes))]
-        for node in nodes:
-            if node.variable_index is not None:
-                node_to_variable[node.identifier] = node.variable_index
-        return node_to_variable
 
 class MdpQuotient(paynt.quotient.quotient.Quotient):
 
+    # if set, an explicit action simulating a random action selection will be added to each state
+    add_dont_care_action = False
+
     def __init__(self, mdp, specification):
         super().__init__(specification=specification)
         updated = payntbind.synthesis.restoreActionsInAbsorbingStates(mdp)
         if updated is not None: mdp = updated
-        mdp = payntbind.synthesis.addDontCareAction(mdp)
+        if MdpQuotient.add_dont_care_action:
+            mdp = payntbind.synthesis.addDontCareAction(mdp)
 
         self.quotient_mdp = mdp
         self.choice_destinations = payntbind.synthesis.computeChoiceDestinations(mdp)
@@ -181,7 +198,10 @@ def decide(self, node, var_name):
     '''
     Build the design space and coloring corresponding to the current decision tree.
     '''
-    def build_coloring(self):
+    def set_depth(self, depth):
+        logger.debug(f"synthesizing tree of depth {depth}")
+        self.decision_tree.set_depth(depth)
+
         # logger.debug("building coloring...")
         variables = self.decision_tree.variables
         variable_name = [v.name for v in variables]
@@ -197,24 +217,23 @@ def build_coloring(self):
         self.is_variable_hole = [False for hole in hole_info]
         domain_max = max([len(domain) for domain in variable_domain])
         bound_domain = list(range(domain_max))
+        node_hole_info = [[] for node in self.decision_tree.collect_nodes()]
         for hole,info in enumerate(hole_info):
-            hole_name,hole_type = info
+            node,hole_name,hole_type = info
+            node_hole_info[node].append( (hole,hole_name,hole_type) )
             if hole_type == "__action__":
                 self.is_action_hole[hole] = True
                 option_labels = self.action_labels
             elif hole_type == "__decision__":
                 self.is_decision_hole[hole] = True
                 option_labels = variable_name
-            elif hole_type == "__bound__":
-                self.is_bound_hole[hole] = True
-                option_labels = bound_domain
             else:
                 self.is_variable_hole[hole] = True
                 variable = variable_name.index(hole_type)
                 option_labels = variables[variable].hole_domain
             self.family.add_hole(hole_name, option_labels)
         self.splitter_count = [0] * self.family.num_holes
-
+        self.decision_tree.root.associate_holes(node_hole_info)
 
 
     def build_unsat_result(self):
diff --git a/paynt/quotient/storm_pomdp_control.py b/paynt/quotient/storm_pomdp_control.py
index 64e1c4e7..1fe9f8b5 100644
--- a/paynt/quotient/storm_pomdp_control.py
+++ b/paynt/quotient/storm_pomdp_control.py
@@ -2,7 +2,7 @@
 import stormpy.pomdp
 import payntbind
 
-import paynt.utils.profiler
+import paynt.utils.timer
 
 from os import makedirs
 
@@ -129,7 +129,7 @@ def run_storm_analysis(self):
         belmc = stormpy.pomdp.BeliefExplorationModelCheckerDouble(self.pomdp, options)
 
         logger.info("starting Storm POMDP analysis")
-        storm_timer = paynt.utils.profiler.Timer()
+        storm_timer = paynt.utils.timer.Timer()
         storm_timer.start()
         result = belmc.check(self.spec_formulas[0], self.paynt_export)   # calls Storm
         storm_timer.stop()
diff --git a/paynt/synthesizer/all_in_one.py b/paynt/synthesizer/all_in_one.py
index 0d21b2ce..a55a6a52 100644
--- a/paynt/synthesizer/all_in_one.py
+++ b/paynt/synthesizer/all_in_one.py
@@ -1,7 +1,7 @@
 import stormpy
 import payntbind
 
-import paynt.utils.profiler
+import paynt.utils.timer
 
 import logging
 logger = logging.getLogger(__name__)
@@ -17,7 +17,7 @@ def __init__(self, all_in_one_program, specification, approach, memory_limit_mb,
         self.approach = approach
         self.family = family
 
-        build_timer = paynt.utils.profiler.Timer()
+        build_timer = paynt.utils.timer.Timer()
         logger.info(f"building {self.approach} all in one MDP")
 
         build_timer.start()
@@ -41,7 +41,7 @@ def __init__(self, all_in_one_program, specification, approach, memory_limit_mb,
 
     def synthesize(self, optimum_threshold=None):
 
-        all_in_one_timer = paynt.utils.profiler.Timer()
+        all_in_one_timer = paynt.utils.timer.Timer()
         logger.info(f"starting all in one analysis")
 
         all_in_one_timer.start()
diff --git a/paynt/synthesizer/decision_tree.py b/paynt/synthesizer/decision_tree.py
index b0ed0422..b6b6a6ee 100644
--- a/paynt/synthesizer/decision_tree.py
+++ b/paynt/synthesizer/decision_tree.py
@@ -1,5 +1,6 @@
 import paynt.synthesizer.synthesizer_ar
 import paynt.quotient.mdp
+import paynt.utils.timer
 
 import stormpy
 import payntbind
@@ -13,9 +14,13 @@ class SynthesizerDecisionTree(paynt.synthesizer.synthesizer_ar.SynthesizerAR):
     tree_depth = 0
     # if set, all trees of size at most tree_depth will be enumerated
     tree_enumeration = False
+    # TODO
+    use_tree_hint = True
 
     def __init__(self, *args):
         super().__init__(*args)
+        self.tree_hint = None
+        self.tree_hint_size = None
 
     @property
     def method_name(self):
@@ -34,7 +39,7 @@ def verify_hole_selection(self, family, hole_selection):
             return
         assignment_value = res.optimality_result.value
         if spec.optimality.improves_optimum(assignment_value):
-            logger.info(f"harmonization achieved value {res.optimality_result.value}")
+            # logger.info(f"harmonization achieved value {res.optimality_result.value}")
             self.num_harmonization_succeeded += 1
             family.analysis_result.improving_assignment = assignment
             family.analysis_result.improving_value = assignment_value
@@ -43,7 +48,7 @@ def verify_hole_selection(self, family, hole_selection):
 
             family_value = family.analysis_result.optimality_result.primary.value
             if(abs(family_value-assignment_value) < 1e-3):
-                logger.info(f"harmonization leads to family skip")
+                # logger.info(f"harmonization leads to family skip")
                 self.num_harmonization_skip += 1
                 family.analysis_result.can_improve = False
 
@@ -89,16 +94,107 @@ def verify_family(self, family):
         self.harmonize_inconsistent_scheduler(family)
 
 
+    def presplit(self, subfamily):
+        family = self.quotient.family
+
+        # fix subfamily first
+        for hole in range(family.num_holes):
+            all_options = family.hole_options(hole)
+            if subfamily.hole_num_options(hole) == family.hole_num_options(hole):
+                continue
+            if not self.quotient.is_variable_hole[hole]:
+                continue
+            options = [option for option in all_options if option<=subfamily.hole_options(hole)[0]]
+            subfamily.hole_set_options(hole,options)
+
+        subfamilies = [subfamily]
+        prototype_family = subfamily.copy()
+        for hole in range(family.num_holes):
+            if prototype_family.hole_num_options(hole) == family.hole_num_options(hole):
+                continue;
+            prototype_options = prototype_family.hole_options(hole)
+            suboptions_list = []
+            all_options = family.hole_options(hole)
+            if not self.quotient.is_variable_hole[hole]:
+                complement = [option for option in all_options if option not in prototype_options]
+                suboptions_list.append(complement)
+            else:
+                complement = [option for option in all_options if option>prototype_options[-1]]
+                if len(complement)>0:
+                    suboptions_list.append(complement)
+
+            for suboptions in suboptions_list:
+                new_subfamily = prototype_family.copy()
+                new_subfamily.hole_set_options(hole,suboptions)
+                subfamilies.append(new_subfamily)
+            prototype_family.hole_set_options(hole,all_options)
+        return subfamilies
+
     def synthesize_tree(self, depth:int):
-        logger.debug(f"synthesizing tree of depth {depth}")
-        self.quotient.decision_tree.set_depth(depth)
-        self.quotient.build_coloring()
+        self.quotient.set_depth(depth)
         self.synthesize(keep_optimum=True)
 
+    def synthesize_tree_sequence(self, opt_result_value):
+        self.tree_hint = None
+
+        for depth in range(SynthesizerDecisionTree.tree_depth+1):
+            self.quotient.set_depth(depth)
+            best_assignment_old = self.best_assignment
+
+            family = self.quotient.family
+            self.explored = 0
+            self.stat = paynt.synthesizer.statistic.Statistic(self)
+            print()
+            self.stat.start(family)
+            families = [family]
+
+            if self.tree_hint is not None and SynthesizerDecisionTree.use_tree_hint:
+                subfamily = family.copy()
+                # find the correct application point
+                application_node = self.quotient.decision_tree.root
+                # for _ in range(self.tree_hint_size,depth):
+                #     application_node = application_node.child_true
+                application_node.apply_hint(subfamily,self.tree_hint)
+
+                # presplit
+                # families = self.presplit(subfamily)
+                # assert family.size == sum([f.size for f in families])
+
+                # hint,
+                families = [subfamily,family]
+
+            for family in families:
+                self.synthesize_one(family)
+            self.stat.finished_synthesis()
+            self.stat.print()
+
+            new_assignment_synthesized = self.best_assignment != best_assignment_old
+            if new_assignment_synthesized:
+                logger.info("printing synthesized assignment below:")
+                logger.info(self.best_assignment)
+
+                if self.best_assignment is not None and self.best_assignment.size == 1:
+                    dtmc = self.quotient.build_assignment(self.best_assignment)
+                    result = dtmc.check_specification(self.quotient.specification)
+                    logger.info(f"double-checking specification satisfiability: {result}")
+
+                if abs( (self.best_assignment_value-opt_result_value)/opt_result_value ) < 1e-3:
+                    break
+
+                self.tree_hint = self.quotient.decision_tree.root
+                self.tree_hint.associate_assignment(self.best_assignment)
+                self.tree_hint_size = depth
+
+            if self.resource_limit_reached():
+                break
+
+
+
     def run(self, optimum_threshold=None, export_evaluation=None):
         paynt_mdp = paynt.models.models.Mdp(self.quotient.quotient_mdp)
         mc_result = paynt_mdp.model_check_property(self.quotient.get_property())
-        logger.info(f"the optimal scheduler has value: {mc_result}")
+        opt_result_value = mc_result.value
+        logger.info(f"the optimal scheduler has value: {opt_result_value}")
 
         self.num_families_considered = 0
         self.num_families_skipped = 0
@@ -110,10 +206,10 @@ def run(self, optimum_threshold=None, export_evaluation=None):
 
         if self.quotient.specification.has_optimality:
             epsilon = 1e-1
-            mc_result_positive = mc_result.value > 0
+            mc_result_positive = opt_result_value > 0
             if self.quotient.specification.optimality.maximizing == mc_result_positive:
                 epsilon *= -1
-            # optimum_threshold = mc_result.value * (1 + epsilon)
+            # optimum_threshold = opt_result_value * (1 + epsilon)
 
         self.set_optimality_threshold(optimum_threshold)
         self.best_assignment = None
@@ -122,21 +218,17 @@ def run(self, optimum_threshold=None, export_evaluation=None):
         if not SynthesizerDecisionTree.tree_enumeration:
             self.synthesize_tree(SynthesizerDecisionTree.tree_depth)
         else:
-            for depth in range(SynthesizerDecisionTree.tree_depth+1):
-                self.synthesize_tree(depth)
-                if self.global_resource_limit_reached():
-                    break
-                # if self.best_assignment is not None:
-                if self.best_assignment_value is not None and abs( (self.best_assignment_value-mc_result.value)/mc_result.value ) <1e-2:
-                    break
+            self.synthesize_tree_sequence(opt_result_value)
 
-        logger.info(f"the optimal scheduler has value: {mc_result}")
+        logger.info(f"the optimal scheduler has value: {opt_result_value}")
         if self.best_assignment is not None:
             logger.info(f"admissible assignment found: {self.best_assignment}")
             if self.quotient.specification.has_optimality:
                 logger.info(f"best assignment has value {self.quotient.specification.optimality.optimum}")
         else:
             logger.info("no admissible assignment found")
+        time_total = paynt.utils.timer.GlobalTimer.read()
+        logger.info(f"synthesis time: {round(time_total, 2)} s")
 
         print()
         logger.info(f"families considered: {self.num_families_considered}")
@@ -148,7 +240,6 @@ def run(self, optimum_threshold=None, export_evaluation=None):
         logger.info(f"harmonizations lead to family skip: {self.num_harmonization_skip}")
 
         print()
-        time_total = self.stat.synthesis_timer_total.read()
         for name,time in self.quotient.coloring.getProfilingInfo():
             time_percent = round(time/time_total*100,1)
             print(f"{name} = {time} s ({time_percent} %)")
diff --git a/paynt/synthesizer/policy_tree.py b/paynt/synthesizer/policy_tree.py
index a536f172..362d4ece 100644
--- a/paynt/synthesizer/policy_tree.py
+++ b/paynt/synthesizer/policy_tree.py
@@ -7,7 +7,7 @@
 import paynt.quotient.quotient
 import paynt.verification.property_result
 from paynt.verification.property import Property
-import paynt.utils.profiler
+import paynt.utils.timer
 
 import paynt.family.smt
 import paynt.synthesizer.conflict_generator.dtmc
@@ -408,7 +408,7 @@ def merge_compatible_policies(self, policy_indices):
     
     def postprocess(self, quotient, prop):
 
-        postprocessing_timer = paynt.utils.profiler.Timer()
+        postprocessing_timer = paynt.utils.timer.Timer()
         postprocessing_timer.start()
         logger.info("post-processing the policy tree...")
 
diff --git a/paynt/synthesizer/statistic.py b/paynt/synthesizer/statistic.py
index fe518d9a..66db8617 100644
--- a/paynt/synthesizer/statistic.py
+++ b/paynt/synthesizer/statistic.py
@@ -1,6 +1,6 @@
 import stormpy.storage
 
-import paynt.utils.profiler
+import paynt.utils.timer
 import paynt.synthesizer.synthesizer
 import paynt.models.models
 
@@ -28,7 +28,7 @@ class Statistic:
 
     # parameters
     status_period_seconds = 3
-    synthesis_timer_total = paynt.utils.profiler.Timer()
+    synthesis_timer_total = paynt.utils.timer.Timer()
     
     def __init__(self, synthesizer):
         
@@ -59,11 +59,12 @@ def __init__(self, synthesizer):
         self.num_policies_merged = None
 
         self.family_size = None
-        self.synthesis_timer = paynt.utils.profiler.Timer()
+        self.synthesis_timer = paynt.utils.timer.Timer()
         self.status_horizon = Statistic.status_period_seconds
 
 
     def start(self, family):
+        logger.info("synthesis initiated, design space: {}".format(family.size_or_order))
         self.family_size = family.size
         self.synthesis_timer.start()
         if not self.synthesis_timer_total.running:
diff --git a/paynt/synthesizer/synthesizer.py b/paynt/synthesizer/synthesizer.py
index 38bea32c..b2b9dec2 100644
--- a/paynt/synthesizer/synthesizer.py
+++ b/paynt/synthesizer/synthesizer.py
@@ -1,4 +1,5 @@
 import paynt.synthesizer.statistic
+import paynt.utils.timer
 
 import logging
 logger = logging.getLogger(__name__)
@@ -68,6 +69,7 @@ def choose_synthesizer(quotient, method, fsc_synthesis=False, storm_control=None
     def __init__(self, quotient):
         self.quotient = quotient
         self.stat = None
+        self.synthesis_timer = None
         self.explored = None
         self.best_assignment = None
         self.best_assignment_value = None
@@ -77,20 +79,21 @@ def method_name(self):
         ''' to be overridden '''
         pass
 
-    def global_time_limit_reached(self):
-        if paynt.utils.profiler.GlobalTimeLimit.limit_reached():
+    def time_limit_reached(self):
+        if (self.synthesis_timer is not None and self.synthesis_timer.time_limit_reached()) or \
+            paynt.utils.timer.GlobalTimer.time_limit_reached():
             logger.info("time limit reached, aborting...")
             return True
         return False
 
-    def global_memory_limit_reached(self):
-        if paynt.utils.profiler.GlobalMemoryLimit.limit_reached():
+    def memory_limit_reached(self):
+        if paynt.utils.timer.GlobalMemoryLimit.limit_reached():
             logger.info("memory limit reached, aborting...")
             return True
         return False
 
-    def global_resource_limit_reached(self):
-        return self.global_time_limit_reached() or self.global_memory_limit_reached()
+    def resource_limit_reached(self):
+        return self.time_limit_reached() or self.memory_limit_reached()
 
     def set_optimality_threshold(self, optimum_threshold):
         if self.quotient.specification.has_optimality and optimum_threshold is not None:
@@ -145,14 +148,18 @@ def synthesize_one(self, family):
         ''' to be overridden '''
         pass
 
-    def synthesize(self, family=None, optimum_threshold=None, keep_optimum=False, return_all=False, print_stats=True):
+    def synthesize(
+        self, family=None, optimum_threshold=None, keep_optimum=False, return_all=False, print_stats=True, timeout=None
+    ):
         '''
         :param family family of assignment to search in
+        :param families alternatively, a list of families can be given
         :param optimum_threshold known bound on the optimum value
         :param keep_optimum if True, the optimality specification will not be reset upon finish
         :param return_all if True and the synthesis returns a family, all assignments will be returned instead of an
             arbitrary one
         :param print_stats if True, synthesis stats will be printed upon completion
+        :param timeout synthesis time limit, seconds
         '''
         if family is None:
             family = self.quotient.family
@@ -160,9 +167,10 @@ def synthesize(self, family=None, optimum_threshold=None, keep_optimum=False, re
             family.constraint_indices = list(range(len(self.quotient.specification.constraints)))
         
         self.set_optimality_threshold(optimum_threshold)
+        self.synthesis_timer = paynt.utils.timer.Timer(timeout)
+        self.synthesis_timer.start()
         self.stat = paynt.synthesizer.statistic.Statistic(self)
         self.explored = 0
-        logger.info("synthesis initiated, design space: {}".format(family.size_or_order))
         self.stat.start(family)
         self.synthesize_one(family)
         if self.best_assignment is not None and self.best_assignment.size > 1 and not return_all:
diff --git a/paynt/synthesizer/synthesizer_ar.py b/paynt/synthesizer/synthesizer_ar.py
index cc942cc6..9fa0d8aa 100644
--- a/paynt/synthesizer/synthesizer_ar.py
+++ b/paynt/synthesizer/synthesizer_ar.py
@@ -100,7 +100,7 @@ def update_optimum(self, family):
         self.quotient.specification.optimality.update_optimum(iv)
         self.best_assignment = ia
         self.best_assignment_value = iv
-        # logger.info(f"new optimum achieved: {iv}")
+        logger.info(f"value {round(iv,4)} achieved after {round(self.synthesis_timer.read(),2)} seconds")
         if isinstance(self.quotient, paynt.quotient.pomdp.PomdpQuotient):
             self.stat.new_fsc_found(family.analysis_result.improving_value, ia, self.quotient.policy_size(ia))
 
@@ -108,7 +108,7 @@ def update_optimum(self, family):
     def synthesize_one(self, family):
         families = [family]
         while families:
-            if self.global_resource_limit_reached():
+            if self.resource_limit_reached():
                 break
             family = families.pop(-1)
             self.verify_family(family)
diff --git a/paynt/synthesizer/synthesizer_decpomdp.py b/paynt/synthesizer/synthesizer_decpomdp.py
index 84ffe523..6fcd52e0 100644
--- a/paynt/synthesizer/synthesizer_decpomdp.py
+++ b/paynt/synthesizer/synthesizer_decpomdp.py
@@ -8,8 +8,6 @@
 
 import paynt.quotient.quotient
 import paynt.quotient.decpomdp
-from ..utils.profiler import Timer
-
 import paynt.verification.property
 
 import logging 
diff --git a/paynt/synthesizer/synthesizer_hybrid.py b/paynt/synthesizer/synthesizer_hybrid.py
index 8936a107..da048586 100644
--- a/paynt/synthesizer/synthesizer_hybrid.py
+++ b/paynt/synthesizer/synthesizer_hybrid.py
@@ -3,8 +3,7 @@
 import paynt.synthesizer.synthesizer_cegis
 
 import paynt.family.smt
-
-from ..utils.profiler import Timer
+import paynt.utils.timer
 
 import logging
 logger = logging.getLogger(__name__)
@@ -25,8 +24,8 @@ class StageControl:
 
     def __init__(self, family_size):
         # timings
-        self.timer_ar = Timer()
-        self.timer_cegis = Timer()
+        self.timer_ar = paynt.utils.timer.Timer()
+        self.timer_cegis = paynt.utils.timer.Timer()
 
         self.family_size = family_size
         self.pruned_ar = 0
diff --git a/paynt/synthesizer/synthesizer_pomdp.py b/paynt/synthesizer/synthesizer_pomdp.py
index 45fc93f6..f469711f 100644
--- a/paynt/synthesizer/synthesizer_pomdp.py
+++ b/paynt/synthesizer/synthesizer_pomdp.py
@@ -8,7 +8,7 @@
 
 import paynt.quotient.quotient
 import paynt.quotient.pomdp
-from ..utils.profiler import Timer
+import paynt.utils.timer
 
 import paynt.verification.property
 
@@ -49,7 +49,7 @@ def __init__(self, quotient, method, storm_control):
             self.synthesis_terminate = False
             self.synthesizer = SynthesizerARStorm       # SAYNT only works with abstraction refinement
             if self.storm_control.iteration_timeout is not None:
-                self.saynt_timer = Timer()
+                self.saynt_timer = paynt.utils.timer.Timer()
                 self.synthesizer.saynt_timer = self.saynt_timer
                 self.storm_control.saynt_timer = self.saynt_timer
 
diff --git a/paynt/utils/profiler.py b/paynt/utils/profiler.py
deleted file mode 100644
index 3ea52431..00000000
--- a/paynt/utils/profiler.py
+++ /dev/null
@@ -1,129 +0,0 @@
-import time
-import psutil
-import os
-
-
-class Timer:
-    def __init__(self):
-        self.running = False
-        self.timer = None   # last timestamp
-        self.time = 0       # total time measured
-
-    @staticmethod
-    def timestamp():
-        return time.perf_counter()
-
-    def reset(self):
-        self.__init__()
-
-    def start(self):
-        if self.running:
-            return
-        self.timer = self.timestamp()
-        self.running = True
-
-    def stop(self):
-        if not self.running:
-            return
-        self.time += self.timestamp() - self.timer
-        self.timer = None
-        self.running = False
-
-    def read(self):
-        if not self.running:
-            return self.time
-        else:
-            return self.time + (self.timestamp() - self.timer)
-
-class Profiler:
-
-    percentage_filter = 5
-
-    @staticmethod
-    def initialize():
-        Profiler.timers = {}        # dictionary of all timers
-        Profiler.running = None     # currently running timer
-        Profiler.paused = []        # stack of paused timers
-        
-        Profiler.timer_total = Timer()
-        Profiler.timer_total.start()
-
-
-    @staticmethod
-    def is_running():
-        return Profiler.running is not None
-
-    @staticmethod
-    def stop():
-        if not Profiler.is_running():
-            return
-        Profiler.running.stop()
-        Profiler.running = None
-
-    @staticmethod
-    def pause():
-        if not Profiler.is_running():
-            return
-        Profiler.running.stop()
-        Profiler.paused += [Profiler.running]
-        Profiler.running = None
-
-    @staticmethod
-    def resume():
-        Profiler.stop()
-        if Profiler.paused:
-            Profiler.running = Profiler.paused.pop(-1)
-            Profiler.running.start()
-    
-    @staticmethod
-    def start(timer_name = "-"):
-        if Profiler.is_running():
-            Profiler.pause()
-        Profiler.timers[timer_name] = Profiler.timers.get(timer_name, Timer())
-        Profiler.timers[timer_name].start()
-        Profiler.running = Profiler.timers[timer_name]
-
-    @staticmethod
-    def print_all():
-        time_total = Profiler.timer_total.read()
-        covered = 0
-        print("profiling report:")
-        for timer_name, timer in Profiler.timers.items():
-            t = timer.read()
-            covered += t
-            percentage = round(t / time_total * 100, 1)
-            if percentage > Profiler.percentage_filter:
-                print(f'> {timer_name} : {percentage}%')
-        print(f"> covered {round(covered / time_total * 100, 0)}% of {round(time_total, 1)} sec")
-
-    @staticmethod
-    def print():
-        Profiler.stop()
-        Profiler.timer_total.stop()
-        Profiler.print_all()
-
-
-class GlobalTimeLimit:
-
-    time_limit_s = None
-
-    @classmethod
-    def start(cls, time_limit_s=None):
-        cls.time_limit_s = time_limit_s
-        cls.global_timer = Timer()
-        cls.global_timer.start()
-
-    @classmethod
-    def limit_reached(cls):
-        return cls.time_limit_s is not None and cls.global_timer.read() > cls.time_limit_s
-
-
-class GlobalMemoryLimit:
-
-    memory_limit_mb = None
-
-    @classmethod
-    def limit_reached(cls):
-        process = psutil.Process(os.getpid())
-        allocated_mb = process.memory_info().rss / (1024 * 1024)
-        return cls.memory_limit_mb is not None and allocated_mb > cls.memory_limit_mb
diff --git a/paynt/utils/timer.py b/paynt/utils/timer.py
new file mode 100644
index 00000000..a6478b39
--- /dev/null
+++ b/paynt/utils/timer.py
@@ -0,0 +1,70 @@
+import time
+import psutil
+import os
+
+
+class Timer:
+
+    def __init__(self, time_limit_seconds=None):
+        self.running = False
+        self.last_timestamp = None
+        self.time = 0 # total time measured
+        self.time_limit_seconds = time_limit_seconds
+
+    @staticmethod
+    def timestamp():
+        return time.perf_counter()
+
+    def reset(self):
+        self.__init__()
+
+    def start(self):
+        if self.running:
+            return
+        self.last_timestamp = self.timestamp()
+        self.running = True
+
+    def stop(self):
+        if not self.running:
+            return
+        self.time += self.timestamp() - self.last_timestamp
+        self.last_timestamp = None
+        self.running = False
+
+    def read(self):
+        if not self.running:
+            return self.time
+        else:
+            return self.time + (self.timestamp() - self.last_timestamp)
+
+    def time_limit_reached(self):
+        return self.time_limit_seconds is not None and self.read() > self.time_limit_seconds
+
+
+class GlobalTimer:
+
+    global_timer = None
+
+    @classmethod
+    def start(cls, time_limit_seconds=None):
+        cls.global_timer = Timer(time_limit_seconds)
+        cls.global_timer.start()
+
+    @classmethod
+    def read(cls):
+        return cls.global_timer.read()
+
+    @classmethod
+    def time_limit_reached(cls):
+        return cls.global_timer is not None and cls.global_timer.time_limit_reached()
+
+
+class GlobalMemoryLimit:
+
+    memory_limit_mb = None
+
+    @classmethod
+    def limit_reached(cls):
+        process = psutil.Process(os.getpid())
+        allocated_mb = process.memory_info().rss / (1024 * 1024)
+        return cls.memory_limit_mb is not None and allocated_mb > cls.memory_limit_mb
diff --git a/payntbind/src/synthesis/quotient/ColoringSmt.cpp b/payntbind/src/synthesis/quotient/ColoringSmt.cpp
index ee0c4e92..ae3ad7f9 100644
--- a/payntbind/src/synthesis/quotient/ColoringSmt.cpp
+++ b/payntbind/src/synthesis/quotient/ColoringSmt.cpp
@@ -118,15 +118,13 @@ ColoringSmt<ValueType>::ColoringSmt(
     }
 
     // collect all path expressions
+    std::vector<z3::expr_vector> path_step_expression;
     z3::expr_vector path_expression(ctx);
     for(std::vector<bool> const& path: getRoot()->paths) {
-        z3::expr_vector expression(ctx);
-        getRoot()->loadPathExpression(path,expression);
-        z3::expr_vector steps(ctx);
-        for(uint64_t step = 0; step < expression.size()-1; ++step) {
-            steps.push_back(expression[step]);
-        }
-        path_expression.push_back(not z3::mk_and(steps) or expression.back());
+        z3::expr_vector step_expression(ctx);
+        getRoot()->loadPathExpression(path,step_expression);
+        path_step_expression.push_back(step_expression);
+        path_expression.push_back(z3::mk_or(step_expression));
         const TreeNode *node = getRoot()->getNodeOfPath(path,path.size()-1);
         const TerminalNode * terminal = dynamic_cast<const TerminalNode *>(node);
         path_action_hole.push_back(terminal->action_hole.hole);
@@ -155,45 +153,48 @@ ColoringSmt<ValueType>::ColoringSmt(
         return;
     }
 
-    timers["ColoringSmt:: create harmonizing variants"].start();
+    timers["ColoringSmt:: create harmonizing variants (1)"].start();
     // create harmonizing variants
     std::vector<const TreeNode::Hole *> all_holes(family.numHoles());
     getRoot()->loadAllHoles(all_holes);
-    std::vector<std::vector<uint64_t>> path_holes(numPaths());
-    for(uint64_t path = 0; path < numPaths(); ++path) {
-        getRoot()->loadPathHoles(getRoot()->paths[path],path_holes[path]);
-    }
     std::vector<z3::expr_vector> hole_what;
     std::vector<z3::expr_vector> hole_with;
     for(const TreeNode::Hole *hole: all_holes) {
         z3::expr_vector what(ctx); what.push_back(hole->solver_variable); hole_what.push_back(what);
         z3::expr_vector with(ctx); with.push_back(hole->solver_variable_harm); hole_with.push_back(with);
     }
+    
+    std::vector<std::vector<std::vector<uint64_t>>> path_step_holes(numPaths());
+    for(uint64_t path = 0; path < numPaths(); ++path) {
+        getRoot()->loadPathStepHoles(getRoot()->paths[path],path_step_holes[path]);
+    }
 
     z3::expr_vector path_expression_harmonizing(ctx);
     for(uint64_t path = 0; path < numPaths(); ++path) {
-        z3::expr e = path_expression[path];
         z3::expr_vector variants(ctx);
-        variants.push_back(e);
-        for(uint64_t hole: path_holes[path]) {
-            variants.push_back(
-                (harmonizing_variable == (int)hole) and e.substitute(hole_what[hole],hole_with[hole])
-            );
+        variants.push_back(path_expression[path]);
+        for(uint64_t step = 0; step < path_step_expression[path].size(); ++step) {
+            for(uint64_t hole: path_step_holes[path][step]) {
+                variants.push_back(
+                    (harmonizing_variable == (int)hole) and path_step_expression[path][step].substitute(hole_what[hole],hole_with[hole])
+                );
+            }
         }
         path_expression_harmonizing.push_back(z3::mk_or(variants));
     }
+    timers["ColoringSmt:: create harmonizing variants (1)"].stop();
 
     for(uint64_t choice = 0; choice < numChoices(); ++choice) {
-        z3::expr_vector path_variants(ctx);
-        for(uint64_t path = 0; path < numPaths(); ++path) {
-            path_variants.push_back(path_expression_harmonizing[path].substitute(substitution_variables,choice_substitution_expr[choice]));
-            // path_variants.push_back(path_expression_harmonizing[path].substitute(substitution_variables,choice_substitution_expr[choice]).simplify());
+        choice_path_expresssion_harm.push_back(z3::expr_vector(ctx));
+    }
+    timers["ColoringSmt:: create harmonizing variants (2)"].start();
+    for(uint64_t path = 0; path < numPaths(); ++path) {
+        for(uint64_t choice = 0; choice < numChoices(); ++choice) {
+            choice_path_expresssion_harm[choice].push_back(path_expression_harmonizing[path].substitute(substitution_variables,choice_substitution_expr[choice]));
         }
-        choice_path_expresssion_harm.push_back(path_variants);
     }
-    timers["ColoringSmt:: create harmonizing variants"].stop();
+    timers["ColoringSmt:: create harmonizing variants (2)"].stop();
 
-    // std::cout << __FUNCTION__ << " end" << std::endl;
     timers[__FUNCTION__].stop();
 }
 
@@ -235,8 +236,8 @@ bool ColoringSmt<ValueType>::check() {
     return sat;
 }
 
-template<typename ValueType>std::vector<std::pair<std::string,std::string>> ColoringSmt<ValueType>::getFamilyInfo() {
-    std::vector<std::pair<std::string,std::string>> hole_info(family.numHoles());
+template<typename ValueType>std::vector<std::tuple<uint64_t,std::string,std::string>> ColoringSmt<ValueType>::getFamilyInfo() {
+    std::vector<std::tuple<uint64_t,std::string,std::string>> hole_info(family.numHoles());
     getRoot()->loadHoleInfo(hole_info);
     return hole_info;
 }
@@ -320,8 +321,8 @@ BitVector ColoringSmt<ValueType>::selectCompatibleChoices(Family const& subfamil
         if(not any_choice_enabled) {
             if(subfamily.isAssignment()) {
                 // STORM_LOG_WARN("Hole assignment does not induce a DTMC, enabling first action...");
-                uint64_t choice = row_groups[state]; // pick the first choice
-                // uint64_t choice = row_groups[state+1]-1; // pick the last choice executing the random choice
+                // uint64_t choice = row_groups[state]; // pick the first choice
+                uint64_t choice = row_groups[state+1]-1; // pick the last choice executing the random choice
                 selection.set(choice,true);
                 for(uint64_t dst: choice_destinations[choice]) {
                     if(not state_reached[dst]) {
@@ -370,8 +371,7 @@ void ColoringSmt<ValueType>::loadUnsatCore(z3::expr_vector const& unsat_core_exp
     for(z3::expr expr: unsat_core_expr) {
         std::istringstream iss(expr.decl().name().str());
         char prefix; iss >> prefix;
-        STORM_LOG_THROW(prefix == 'h' or prefix == 'p', storm::exceptions::UnexpectedException, "unexpected expression label");
-        if(prefix == 'h') {
+        if(prefix == 'h' or prefix == 'z') {
             // uint64_t hole; iss >> hole;
             // this->unsat_core.emplace_back(hole,numPaths());
             continue;
diff --git a/payntbind/src/synthesis/quotient/ColoringSmt.h b/payntbind/src/synthesis/quotient/ColoringSmt.h
index aaf0832f..668dc5e0 100644
--- a/payntbind/src/synthesis/quotient/ColoringSmt.h
+++ b/payntbind/src/synthesis/quotient/ColoringSmt.h
@@ -30,7 +30,7 @@ class ColoringSmt {
     );
 
     /** For each hole, get a list of name-type pairs.  */
-    std::vector<std::pair<std::string,std::string>> getFamilyInfo();
+    std::vector<std::tuple<uint64_t,std::string,std::string>> getFamilyInfo();
 
     /** Whether a check for an admissible family member is done before choice selection. */
     const bool CHECK_FAMILY_CONSISTENCE = true;
diff --git a/payntbind/src/synthesis/quotient/TreeNode.cpp b/payntbind/src/synthesis/quotient/TreeNode.cpp
index 101e765d..db3611c6 100644
--- a/payntbind/src/synthesis/quotient/TreeNode.cpp
+++ b/payntbind/src/synthesis/quotient/TreeNode.cpp
@@ -147,8 +147,8 @@ void TerminalNode::createHoles(Family& family) {
     action_hole.createSolverVariable();
 }
 
-void TerminalNode::loadHoleInfo(std::vector<std::pair<std::string,std::string>> & hole_info) const {
-    hole_info[action_hole.hole] = std::make_pair(action_hole.name,"__action__");
+void TerminalNode::loadHoleInfo(std::vector<std::tuple<uint64_t,std::string,std::string>> & hole_info) const {
+    hole_info[action_hole.hole] = std::make_tuple(identifier,action_hole.name,"__action__");
 }
 
 void TerminalNode::createPaths(z3::expr_vector const& substitution_variables) {
@@ -171,10 +171,11 @@ void TerminalNode::loadAllHoles(std::vector<const Hole *> & holes) const {
     holes[action_hole.hole] = &action_hole;
 }
 
-void TerminalNode::loadPathHoles(std::vector<bool> const& path, std::vector<uint64_t> & holes) const {
-    holes.push_back(action_hole.hole);
+void TerminalNode::loadPathStepHoles(std::vector<bool> const& path, std::vector<std::vector<uint64_t>> & step_holes) const {
+    step_holes.push_back({action_hole.hole});
 }
 
+
 void TerminalNode::addFamilyEncoding(Family const& subfamily, z3::solver& solver) const {
     action_hole.addDomainEncoding(subfamily,solver);
 }
@@ -265,11 +266,11 @@ void InnerNode::createHoles(Family& family) {
     child_false->createHoles(family);
 }
 
-void InnerNode::loadHoleInfo(std::vector<std::pair<std::string,std::string>> & hole_info) const {
-    hole_info[decision_hole.hole] = std::make_pair(decision_hole.name,"__decision__");
+void InnerNode::loadHoleInfo(std::vector<std::tuple<uint64_t,std::string,std::string>> & hole_info) const {
+    hole_info[decision_hole.hole] = std::make_tuple(identifier,decision_hole.name,"__decision__");
     for(uint64_t variable = 0; variable < numVariables(); ++variable) {
         Hole const& hole = variable_hole[variable];
-        hole_info[hole.hole] = std::make_pair(hole.name,variable_name[variable]);
+        hole_info[hole.hole] = std::make_tuple(identifier,hole.name,variable_name[variable]);
     }
     child_true->loadHoleInfo(hole_info);
     child_false->loadHoleInfo(hole_info);
@@ -289,8 +290,9 @@ void InnerNode::createPaths(z3::expr_vector const& substitution_variables) {
         step_true_options.push_back( dv == (int)variable and substitution_variables[variable] <= vv);
         step_false_options.push_back(dv == (int)variable and substitution_variables[variable]  > vv);
     }
-    step_true = z3::mk_or(step_true_options);
-    step_false = z3::mk_or(step_false_options);
+    // mind the negation below
+    step_true = not z3::mk_or(step_true_options);
+    step_false = not z3::mk_or(step_false_options);
 
     // create paths
     for(bool condition: {true,false}) {
@@ -347,14 +349,17 @@ void InnerNode::loadAllHoles(std::vector<const Hole *> & holes) const {
     child_false->loadAllHoles(holes);
 }
 
-void InnerNode::loadPathHoles(std::vector<bool> const& path, std::vector<uint64_t> & holes) const {
+void InnerNode::loadPathStepHoles(std::vector<bool> const& path, std::vector<std::vector<uint64_t>> & step_holes) const {
+    std::vector<uint64_t> holes;
     holes.push_back(decision_hole.hole);
     for(Hole const& hole: variable_hole) {
         holes.push_back(hole.hole);
     }
-    getChild(path[depth])->loadPathHoles(path,holes);
+    step_holes.push_back(holes);
+    getChild(path[depth])->loadPathStepHoles(path,step_holes);
 }
 
+
 void InnerNode::addFamilyEncoding(Family const& subfamily, z3::solver& solver) const {
     decision_hole.addDomainEncoding(subfamily,solver);
     for(Hole const& hole: variable_hole) {
diff --git a/payntbind/src/synthesis/quotient/TreeNode.h b/payntbind/src/synthesis/quotient/TreeNode.h
index f3fab9a3..4937ec02 100644
--- a/payntbind/src/synthesis/quotient/TreeNode.h
+++ b/payntbind/src/synthesis/quotient/TreeNode.h
@@ -110,7 +110,7 @@ class TreeNode {
     /** Create all holes and solver variables associated with this node. */
     virtual void createHoles(Family& family) {}
     /** Collect name and type (action,decision, or variable) of each hole. */
-    virtual void loadHoleInfo(std::vector<std::pair<std::string,std::string>> & hole_info) const {}
+    virtual void loadHoleInfo(std::vector<std::tuple<uint64_t,std::string,std::string>> & hole_info) const {}
 
     /** Create a list of paths from this node. */
     virtual void createPaths(z3::expr_vector const& substitution_variables) {}
@@ -120,7 +120,7 @@ class TreeNode {
     virtual void loadPathExpression(std::vector<bool> const& path, z3::expr_vector & expression) const {}
     /** TODO */
     virtual void loadAllHoles(std::vector<const Hole *> & holes) const {};
-    virtual void loadPathHoles(std::vector<bool> const& path, std::vector<uint64_t> & holes) const {};
+    virtual void loadPathStepHoles(std::vector<bool> const& path, std::vector<std::vector<uint64_t>> & step_holes) const {};
 
     /** Add encoding of hole option in the given family. */
     virtual void addFamilyEncoding(Family const& subfamily, z3::solver & solver) const {}
@@ -165,12 +165,12 @@ class TerminalNode: public TreeNode {
     );
 
     void createHoles(Family& family) override;
-    void loadHoleInfo(std::vector<std::pair<std::string,std::string>> & hole_info) const override;
+    void loadHoleInfo(std::vector<std::tuple<uint64_t,std::string,std::string>> & hole_info) const override;
     void createPaths(z3::expr_vector const& substitution_variables) override;
     void createPathsHarmonizing(z3::expr_vector const& substitution_variables, z3::expr const& harmonizing_variable) override;
     void loadPathExpression(std::vector<bool> const& path, z3::expr_vector & expression) const override;
     void loadAllHoles(std::vector<const Hole *> & holes) const override;
-    void loadPathHoles(std::vector<bool> const& path, std::vector<uint64_t> & holes) const override;
+    void loadPathStepHoles(std::vector<bool> const& path, std::vector<std::vector<uint64_t>> & step_holes) const override;
 
     void addFamilyEncoding(Family const& subfamily, z3::solver & solver) const override;
     bool isPathEnabledInState(
@@ -215,12 +215,12 @@ class InnerNode: public TreeNode {
     );
 
     void createHoles(Family& family) override;
-    void loadHoleInfo(std::vector<std::pair<std::string,std::string>> & hole_info) const override;
+    void loadHoleInfo(std::vector<std::tuple<uint64_t,std::string,std::string>> & hole_info) const override;
     void createPaths(z3::expr_vector const& substitution_variables) override;
     void createPathsHarmonizing(z3::expr_vector const& substitution_variables, z3::expr const& harmonizing_variable) override;
     void loadPathExpression(std::vector<bool> const& path, z3::expr_vector & expression) const override;
     void loadAllHoles(std::vector<const Hole *> & holes) const override;
-    void loadPathHoles(std::vector<bool> const& path, std::vector<uint64_t> & holes) const override;
+    void loadPathStepHoles(std::vector<bool> const& path, std::vector<std::vector<uint64_t>> & step_holes) const override;
 
     void addFamilyEncoding(Family const& subfamily, z3::solver & solver) const override;
     bool isPathEnabledInState(