class RxnConSystem:  # pylint: disable=too-many-instance-attributes
    """RxnConSystem holds all reactions and contingencies pertaining to a rxncon system."""
    def __init__(self, reactions: List[Reaction], contingencies: List[Contingency]) -> None:
        self.reactions = reactions
        self.contingencies = contingencies

        self._components = []  # type: List[Spec]
        self._states = []  # type: List[State]
        self._produced_states = []  # type: List[State]
        self._consumed_states = []  # type: List[State]
        self._synthesised_states = []  # type: List[State]
        self._global_states = []  # type: List[State]

        # Components get synthesised in a 'fully neutral' state. What this state is depends on all reactions
        # in the system, so this 'fully neutral' state has to be expanded into the combination of all neutral
        # states pertaining to the component.
        self._expand_fully_neutral_states()
        # These states are lazily generated.
        self._calculate_produced_states()
        self._calculate_consumed_states()
        self._calculate_synthesised_states()
        self._calculate_global_states()
        self._calculate_states()

        # Contingencies containing states at a non-elemental resolution are expanded to an OR of elemental
        # contingencies.
        self._expand_non_elemental_states()
        self._structure_contingencies()

        self.validate()

    def reaction_number(self, reaction: Reaction) -> int:
        return self.reactions.index(reaction)

    def validate(self) -> None:
        if not self.reactions:
            raise AssertionError('No reactions, boring!')

        missing_states = self._missing_states()
        if missing_states:
            raise AssertionError('State(s) {0} appear(s) in contingencies, but is not produced or consumed'
                                 .format(', '.join(str(state) for state in missing_states)))

        missing_reactions = self._missing_reactions()
        if missing_reactions:
            raise AssertionError('Reactions(s) {0} appear(s) in contingencies, but are not defined in reaction list'
                                 .format(', '.join(str(reaction) for reaction in missing_reactions)))

        unsatisfiable_contingencies = self._unsatisfiable_contingencies()
        if unsatisfiable_contingencies:
            reason_str = '\n'.join('{} : {}'.format(rxn, reason) for rxn, reason in unsatisfiable_contingencies)
            raise AssertionError('Unsatisfiable reaction contingencies:\n{}'.format(reason_str))

    def components(self) -> List[Spec]:
        if not self._components:
            self._calculate_components()
        return self._components

    def contingencies_for_reaction(self, reaction: Reaction) -> List[Contingency]:
        assert reaction in self.reactions
        return [x for x in self.contingencies if x.reaction == reaction]

    def q_contingencies_for_reaction(self, reaction: Reaction) -> List[Contingency]:
        assert reaction in self.reactions
        return [x for x in self.contingencies if x.reaction == reaction and x.contingency_type
                in [ContingencyType.positive, ContingencyType.negative]]

    def s_contingencies_for_reaction(self, reaction: Reaction) -> List[Contingency]:
        assert reaction in self.reactions
        return [x for x in self.contingencies if x.reaction == reaction and x.contingency_type
                in [ContingencyType.requirement, ContingencyType.inhibition]]

    @property
    def states(self) -> List[State]:
        if not self._states:
            self._calculate_states()
        return self._states

    @property
    def produced_states(self) -> List[State]:
        if not self._produced_states:
            self._calculate_produced_states()
        return self._produced_states

    @property
    def consumed_states(self) -> List[State]:
        if not self._consumed_states:
            self._calculate_consumed_states()
        return self._consumed_states

    @property
    def synthesised_states(self) -> List[State]:
        if not self._synthesised_states:
            self._calculate_synthesised_states()
        return self._synthesised_states

    @property
    def global_states(self) -> List[State]:
        if not self._global_states:
            self._calculate_global_states()
        return self._global_states

    def states_for_component(self, component: Spec) -> List[State]:
        """Returns all the States that live on a certain Component."""
        assert component.is_component_spec
        return [x for x in self.states if component.to_non_struct_spec() in x.components]

    def states_for_component_grouped(self, component: Spec) -> Dict[Spec, List[State]]:
        """Returns all the States that live on a certain component, grouped by their Spec.
        Within a group, each State is mutually exclusive with each other State."""
        states = self.states_for_component(component)  # type: List[State]
        grouped = defaultdict(list)  # type: Dict[Spec, List[State]]

        while states:
            state = states.pop()

            for spec in state.specs:
                if spec.to_component_spec() != component:
                    continue

                grouped[spec].append(state)

        return grouped

    def complement_states(self, state: State) -> List[State]:
        """Returns all States mutually exclusive with the State given. For multi-component States,
        all components are included."""
        states = []  # type: List[State]
        for component in state.components:
            states += self.complement_states_for_component(component, state)

        return states

    def complement_states_for_component(self, component: Spec, state: State) -> List[State]:
        """Returns all States mutually exclusive with the State given that live on the Component given."""
        if not state.is_structured:
            for group in self.states_for_component_grouped(component).values():
                if state in group:
                    return [x for x in group if x != state]
            raise AssertionError
        else:
            # The structure information is first thrown away to determine the mutually exclusive
            # states, and then applied again (as far as possible).
            complements = self.complement_states_for_component(component.to_non_struct_spec(),
                                                               state.to_non_structured())
            return [x.to_structured_from_state(state) for x in complements]

    def _calculate_components(self) -> None:
        """Determines all components in the system, and stores it in self._components."""
        components = []  # type: List[Spec]
        for reaction in self.reactions:
            components += [spec.to_non_struct_spec() for spec in reaction.components_lhs] + \
                          [spec.to_non_struct_spec() for spec in reaction.components_rhs]

        self._components = list(set(components))

    def _calculate_states(self) -> None:
        self._states = list(
            set(self.produced_states + self.consumed_states + self.synthesised_states + self.global_states))

    def _calculate_produced_states(self) -> None:
        states = []  # type: List[State]
        for reaction in self.reactions:
            states += [state.to_non_structured() for state in reaction.produced_states]

        self._produced_states = list(set(states))

    def _calculate_consumed_states(self) -> None:
        states = []  # type: List[State]
        for reaction in self.reactions:
            states += [state.to_non_structured() for state in reaction.consumed_states]

        self._consumed_states = list(set(states))

    def _calculate_synthesised_states(self) -> None:
        states = []  # type: List[State]
        for reaction in self.reactions:
            states += [state.to_non_structured() for state in reaction.synthesised_states]

        self._synthesised_states = list(set(states))

    def _calculate_global_states(self) -> None:
        states = []  # type: List[State]
        for contingency in self.contingencies:
            states += [state for state in contingency.effector.states if state.is_global]

        self._global_states = list(set(states))

    def _expand_fully_neutral_states(self) -> None:
        for reaction in self.reactions:
            self._expand_reaction_terms(reaction.terms_lhs)
            self._expand_reaction_terms(reaction.terms_rhs)
            reaction.invalidate_state_cache()

    def _expand_reaction_terms(self, terms: List[ReactionTerm]) -> None:
        """Expands the reaction rules that contain a FullyNeutralState, which appear in synthesis reactions.
        This state is expanded in a combination of all neutral states for the component."""
        for term in terms:
            if FullyNeutralState() in term.states:
                existing_states = [state for state in term.states if state != FullyNeutralState()]
                new_states = [state for component in term.specs for state in self.states_for_component(component)
                              if state.is_neutral and state != FullyNeutralState() and state not in existing_states and
                              not any(state.is_mutually_exclusive_with(existing) for existing in existing_states)]

                term.states = existing_states + new_states

    def _expand_non_elemental_states(self) -> None:
        """Expands the non-elemental States appearing in contingencies, these get expanded into a Union of
        the elemental States of which the non-elemental State is a superset."""
        def expanded_effector(effector: Effector) -> Effector:
            if isinstance(effector, StateEffector):
                if effector.expr.is_elemental:
                    return effector
                else:
                    elemental_states = [state.to_structured_from_state(effector.expr)
                                        for state in self.states if state.is_subset_of(effector.expr)]
                    assert elemental_states, 'Could not find elemental states which are subset of the non-elemental ' \
                                             'state {}'.format(effector.expr)
                    assert all(state.is_elemental for state in elemental_states)

                    LOGGER.info('expanded_effector: {} -> {}'
                                .format(str(effector.expr), ' | '.join(str(x) for x in elemental_states)))
                    return OrEffector(*(StateEffector(x) for x in elemental_states), name=str(effector.expr))
            elif isinstance(effector, AndEffector):
                return AndEffector(*(expanded_effector(x) for x in effector.exprs), name=effector.name,
                                   equivs=effector.equivs)
            elif isinstance(effector, OrEffector):
                return OrEffector(*(expanded_effector(x) for x in effector.exprs), name=effector.name,
                                  equivs=effector.equivs)
            elif isinstance(effector, NotEffector):
                return NotEffector(expanded_effector(effector.expr), name=effector.name)
            else:
                raise AssertionError

        for contingency in self.contingencies:
            contingency.effector = expanded_effector(contingency.effector)

    def _structure_contingencies(self) -> None:
        """Structures (i.e. augment with topological data) all contingencies."""
        structured_conts = []

        for rxn in self.reactions:
            counter_start = 2
            for con in self.contingencies_for_reaction(rxn):
                structured_conts.append(con.to_structured(counter_start))
                counter_start += 100

        self.contingencies = structured_conts

    def _missing_states(self) -> List[State]:
        """Returns the States that appear in contingencies, but which are not produced, consumed, synthesised
        or degraded by the Reactions."""
        required_states = []  # type: List[State]
        for contingency in self.contingencies:
            required_states += [state.to_non_structured() for state in contingency.effector.states if
                                not state.is_global]

        return [state for state in required_states if state not in self.states]

    def _missing_reactions(self) -> List[Reaction]:
        """Returns the Reactions that appear in the contingency list, but which are not in the Reaction list."""
        required_reactions = []
        for contingency in self.contingencies:
            required_reactions.append(contingency.reaction)

        return [reaction for reaction in required_reactions if reaction not in self.reactions]

    def _unsatisfiable_contingencies(self) -> List[Tuple[Reaction, str]]:
        """Determines the contingencies that are not satisfiable, returns a list of (Reaction, str) where
        the str contains the human-readable reason for the contingency not being satisfiable."""
        unsatisfiable = []

        for reaction in self.reactions:
            contingencies = self.s_contingencies_for_reaction(reaction)

            # Make sure the contingency does not contain the states produced / consumed by the reaction.
            states = (state for contingency in contingencies
                      for state in contingency.effector.states)

            for state in states:
                # We need to be talking about the states mentioning the reactants.
                if not all(spec.struct_index in (0, 1) for spec in state.specs):
                    continue

                if any(count > 1 for count in Counter(reaction.components_rhs).values()):
                    continue

                # States appear non-structured in the '.produced_states' etc. properties.
                state = state.to_non_structured()

                if state in reaction.produced_states and state not in reaction.synthesised_states:
                    unsatisfiable.append((reaction, 'Produced state {} appears in contingencies'.format(str(state))))
                if state in reaction.consumed_states and state not in reaction.degraded_states:
                    unsatisfiable.append((reaction, 'Consumed state {} appears in contingencies'.format(str(state))))

            # Make sure at least one solution is there (this might still contain mutually exclusive states)
            total_set = UniversalSet()  # type: Set[State]
            for contingency in contingencies:
                total_set = Intersection(total_set,
                                         contingency.to_venn_set())  # pylint: disable=redefined-variable-type

            solutions = total_set.calc_solutions()
            if len(solutions) == 0:
                unsatisfiable.append((reaction, 'Zero consistent solutions found.'))

            # Make sure that at least one solution doesn't contain mutually exclusive states.
            local_unsatisfiables = []
            at_least_one_consistent_soln = False

            for solution in solutions:
                trues = [state for state, val in solution.items() if val]
                if any(state.is_mutually_exclusive_with(other) for state, other in product(trues, trues)):
                    state, other = next((state, other) for state, other in product(trues, trues) if
                                        state.is_mutually_exclusive_with(other))
                    local_unsatisfiables.append(
                        (reaction, 'State {} mutually exclusive with {}.'.format(str(state), str(other))))
                else:
                    at_least_one_consistent_soln = True

            if not at_least_one_consistent_soln:
                unsatisfiable += local_unsatisfiables

        return unsatisfiable class BooleanModel:
    """Holds all data describing a Boolean model: a list of targets, a list of update rules and
    a list of initial conditions."""
    def __init__(self, targets: List['Target'], update_rules: List['UpdateRule'],
                 initial_conditions: 'BooleanModelState') -> None:
        self.update_rules = sorted(update_rules)
        self.initial_conditions = initial_conditions
        self._state_targets = {str(x): x for x in targets if isinstance(x, StateTarget)}
        self._reaction_targets = {str(x): x for x in targets if isinstance(x, ReactionTarget)}
        self._knockout_targets = {str(x): x for x in targets if isinstance(x, KnockoutTarget)}
        self._overexpression_targets = {str(x): x for x in targets if isinstance(x, OverexpressionTarget)}
        self._validate_update_rules()
        self._validate_initial_conditions()

        self.current_state = None  # type: Optional[BooleanModelState]

    def set_initial_condition(self, target: 'Target', value: bool) -> None:
        self.initial_conditions.set_target(target, value)

    def update_rule_by_target(self, target: 'Target') -> 'UpdateRule':
        for rule in self.update_rules:
            if rule.target == target:
                return rule

        raise KeyError

    def state_target_by_name(self, name: str) -> 'StateTarget':
        return self._state_targets[name]

    def reaction_target_by_name(self, name: str) -> 'ReactionTarget':
        return self._reaction_targets[name]

    def knockout_target_by_name(self, name: str) -> 'KnockoutTarget':
        return self._knockout_targets[name]

    def overexpression_target_by_name(self, name: str) -> 'OverexpressionTarget':
        return self._overexpression_targets[name]

    def step(self) -> None:
        """Takes one timestep in the Boolean model. This is rather inefficient, but not meant for
        actual simulations, this is only used in the unit tests that test all different motifs
        and their desired steady states."""
        if not self.current_state:
            self.current_state = deepcopy(self.initial_conditions)
        else:
            new_state = dict()
            for rule in self.update_rules:
                new_state[rule.target] = rule.factor.eval_boolean_func(self.current_state.target_to_value)

            self.current_state = BooleanModelState(new_state)

    def calc_steady_state(self) -> 'BooleanModelState':
        """Calculates the steady state by taking max MAX_STEADY_STATE_ITERS steps. If no steady state
        found, raises."""
        iters = 0

        while iters < MAX_STEADY_STATE_ITERS:
            prev = deepcopy(self.current_state)
            self.step()
            if prev == self.current_state:
                assert isinstance(prev, BooleanModelState)
                return prev

            iters += 1

        raise AssertionError('Could not find steady state.')

    def _validate_update_rules(self) -> None:
        """Assert that all targets appearing on the RHS in an update rule have their own LHS."""
        all_lhs_targets = []  # type: List[Target]
        all_rhs_targets = []  # type: List[Target]
        for rule in self.update_rules:
            all_lhs_targets.append(rule.target)
            all_rhs_targets += rule.factor_targets

        assert all(x in all_lhs_targets for x in all_rhs_targets)

    def _validate_initial_conditions(self) -> None:
        self.initial_conditions.validate_by_model(self) class Target(metaclass=ABCMeta):
    """Abstract base class for the different targets."""
    def __hash__(self) -> int:
        return hash(str(self))

    def __repr__(self) -> str:
        return str(self) class ReactionTarget(Target):
    """Reaction target of the boolean model. For all non-degrading reactions the relation between
    rxncon reactions and Boolean targets is 1:1. The relation for degrading reactions is more difficult
    since (1) the contingencies determine what the reaction degrades (which obviously becomes problematic
    in the case of a logical disjunction), and (2) the degradation of bonds should produce empty binding
    partners. We refer to our paper."""
    def __init__(self, reaction_parent: Reaction, contingency_variant: Optional[int]=None,
                 interaction_variant: Optional[int] = None, contingency_factor: VennSet['StateTarget']=None) -> None:
        self.reaction_parent = reaction_parent  # type: Reaction
        self.produced_targets = [StateTarget(x) for x in reaction_parent.produced_states]  # type: List[StateTarget]
        self.consumed_targets = [StateTarget(x) for x in reaction_parent.consumed_states]  # type: List[StateTarget]
        self.synthesised_targets = [StateTarget(x) for x in
                                    reaction_parent.synthesised_states]  # type: List[StateTarget]
        self.degraded_targets = [StateTarget(x) for x in reaction_parent.degraded_states]  # type: List[StateTarget]

        self.contingency_variant_index = contingency_variant
        self.interaction_variant_index = interaction_variant

        if contingency_factor is None:
            self.contingency_factor = UniversalSet()  # type: VennSet[StateTarget]
        else:
            self.contingency_factor = contingency_factor  # type: VennSet[StateTarget]

    def __hash__(self) -> int:
        return hash(str(self))

    def __eq__(self, other: object) -> bool:
        if not isinstance(other, Target):
            return NotImplemented
        return isinstance(other, ReactionTarget) and self.reaction_parent == other.reaction_parent and \
            self.contingency_variant_index == other.contingency_variant_index and \
            self.interaction_variant_index == other.interaction_variant_index

    def __str__(self) -> str:
        suffix = ''
        if self.interaction_variant_index is not None and self.contingency_variant_index is not None:
            suffix = '#c{}/i{}'.format(self.contingency_variant_index, self.interaction_variant_index)
        elif self.contingency_variant_index is not None and self.interaction_variant_index is None:
            suffix = '#c{}'.format(self.contingency_variant_index)
        elif self.interaction_variant_index is not None and self.contingency_variant_index is None:
            suffix = '#i{}'.format(self.interaction_variant_index)

        return str(self.reaction_parent) + suffix

    def __repr__(self) -> str:
        return str(self)

    def produces(self, state_target: 'StateTarget') -> bool:
        return state_target in self.produced_targets

    def consumes(self, state_target: 'StateTarget') -> bool:
        return state_target in self.consumed_targets

    def synthesises(self, state_target: 'StateTarget') -> bool:
        return state_target in self.synthesised_targets

    def degrades(self, state_target: 'StateTarget') -> bool:
        return state_target in self.degraded_targets

    @property
    def components_lhs(self) -> List[Spec]:
        return self.reaction_parent.components_lhs

    @property
    def components_rhs(self) -> List[Spec]:
        return self.reaction_parent.components_rhs

    @property
    def degraded_components(self) -> List[Spec]:
        return [component for component in self.components_lhs if component not in self.components_rhs]

    @property
    def synthesised_components(self) -> List[Spec]:
        return [component for component in self.components_rhs if component not in self.components_lhs]

    def degrades_component(self, spec: Spec) -> bool:
        assert spec.is_component_spec
        return spec in self.degraded_components

    def synthesises_component(self, spec: Spec) -> bool:
        assert spec.is_component_spec
        return spec in self.synthesised_components

    def is_output(self) -> bool:
        return isinstance(self.reaction_parent, OutputReaction) class KnockoutTarget(ComponentStateTarget):
    """When enabled, KnockoutTargets are ANDed into the update rule for each target, one KnockoutTarget
    per component. The KnockoutTarget itself has a trivial update rule. By changing the initial conditions
    for the KnockoutTarget we can knock out all states for a particular component, making them always FALSE."""
    def __eq__(self, other: object) -> bool:
        if not isinstance(other, Target):
            return NotImplemented
        return isinstance(other, KnockoutTarget) and self.component == other.component

    def __str__(self) -> str:
        return 'Knockout<{}>'.format(str(self.component))

    def __repr__(self) -> str:
        return str(self)

    def __hash__(self) -> int:
        return hash(str(self)) class OverexpressionTarget(ComponentStateTarget):
    """Similar to KnockoutTarget, but now (a set of) states per component can be made always TRUE."""
    def __eq__(self, other: object) -> bool:
        if not isinstance(other, Target):
            return NotImplemented
        return isinstance(other, OverexpressionTarget) and self.component == other.component

    def __str__(self) -> str:
        return 'Overexpression<{}>'.format(str(self.component))

    def __repr__(self) -> str:
        return str(self)

    def __hash__(self) -> int:
        return hash(str(self)) class StateTarget(Target):
    """State target of the Boolean model. The relation between rxncon states and Boolean state targets
    is generally 1:1, but not quite: the components that carry no internal state are assigned so-called
    ComponentStateTargets (see next class.)"""
    def __init__(self, state_parent: State) -> None:
        self.state_parent = state_parent

    def __hash__(self) -> int:
        return hash(str(self))

    def __str__(self) -> str:
        return str(self.state_parent)

    def __eq__(self, other: object) -> bool:
        if not isinstance(other, Target):
            return NotImplemented
        return isinstance(other, StateTarget) and self.state_parent == other.state_parent

    def is_produced_by(self, reaction_target: ReactionTarget) -> bool:
        return reaction_target.produces(self)

    def is_consumed_by(self, reaction_target: ReactionTarget) -> bool:
        return reaction_target.consumes(self)

    def is_synthesised_by(self, reaction_target: ReactionTarget) -> bool:
        return reaction_target.synthesises(self)

    def is_degraded_by(self, reaction_target: ReactionTarget) -> bool:
        return reaction_target.degrades(self)

    def is_input(self) -> bool:
        return self.state_parent.is_global

    @property
    def components(self) -> List[Spec]:
        return self.state_parent.components

    def shares_component_with(self, other_target: 'StateTarget') -> bool:
        return any(x in other_target.components for x in self.components)

    @property
    def is_neutral(self) -> bool:
        return self.state_parent.is_neutral

    @property
    def is_homodimer(self) -> bool:
        return self.state_parent.is_homodimer

    @property
    def neutral_targets(self) -> List['StateTarget']:
        return [StateTarget(x) for x in self.state_parent.neutral_states]

    @property
    def is_interaction(self) -> bool:
        return isinstance(self.state_parent, InteractionState)

    def is_mutually_exclusive_with(self, other: 'StateTarget') -> bool:
        return self.state_parent.is_mutually_exclusive_with(other.state_parent)

    def complementary_state_targets(self, rxnconsys: RxnConSystem, component: Spec) -> List['StateTarget']:
        others = rxnconsys.complement_states_for_component(component, self.state_parent)
        return [StateTarget(x) for x in others] class UpdateRule:
    """The UpdateRule for a Target gives the Boolean value of the target at time t+1 when we
    evaluate the factor at time t."""
    def __init__(self, target: Target, factor: VennSet[Target]) -> None:
        self.target = target
        self.factor = factor

    def __eq__(self, other: object) -> bool:
        if not isinstance(other, UpdateRule):
            return NotImplemented
        else:
            return self.target == other.target and self.factor.is_equivalent_to(other.factor)

    def __lt__(self, other: object) -> bool:
        if not isinstance(other, UpdateRule):
            return NotImplemented
        else:
            return str(self.target) < str(other.target)

    def __str__(self) -> str:
        return "target: {0}, factors: {1}".format(self.target, self.factor)

    @property
    def factor_targets(self) -> List[Target]:
        return self.factor.values class SmoothingStrategy(Enum):
    """To overcome non biological oscillatory behaviour during the simulation we introduced
    smoothings."""
    no_smoothing = 'no_smoothing'
    smooth_production_sources = 'smooth_production_sources' class KnockoutStrategy(Enum):
    """For which states should knockout rules be generated: none, the neutral ones, or all."""
    no_knockout = 'no_knockout'
    knockout_neutral_states = 'knockout_neutral_states'
    knockout_all_states = 'knockout_all_states' class OverexpressionStrategy(Enum):
    """For which states should overexpression rules be generated: none, the neutral ones, or all."""
    no_overexpression = 'no_overexpression'
    overexpress_neutral_states = 'overexpress_neutral_states'
    overexpress_all_states = 'overexpress_all_states' def boolean_model_from_rxncon(rxncon_sys: RxnConSystem,
                              smoothing_strategy: SmoothingStrategy=SmoothingStrategy.no_smoothing,
                              knockout_strategy: KnockoutStrategy=KnockoutStrategy.no_knockout,
                              overexpression_strategy: OverexpressionStrategy=OverexpressionStrategy.no_overexpression,
                              k_plus_strict: bool=True, k_minus_strict: bool=True) -> BooleanModel:

    def initial_conditions(reaction_targets: List[ReactionTarget], state_targets: List[StateTarget],
                           knockout_targets: List[KnockoutTarget], overexpression_targets: List[OverexpressionTarget]) \
            -> BooleanModelState:
        """As default all the neutral state targets are set to True. All other state targets as well
        as all reaction targets are set to False."""
        conds = {}  # type: Dict[Target, bool]

        for reaction_target in reaction_targets:
            conds[reaction_target] = False

        for state_target in state_targets:
            # Neutral state targets are True.
            if state_target.is_neutral:
                conds[state_target] = True
            # All reaction targets and non-neutral state targets are False.
            else:
                conds[state_target] = False

        for knockout_target in knockout_targets:
            conds[knockout_target] = False

        for overexpression_target in overexpression_targets:
            conds[overexpression_target] = False

        return BooleanModelState(conds)

    def calc_component_presence_factors() -> Tuple[Dict[Spec, VennSet[StateTarget]], List[ComponentStateTarget]]:
        """The form of the component presence factor is:
            (state_a1 | ... | state_an) & (state_b1 | ... | state_bm) & ...

        Mutually exclusive states are combined by boolean OR (state_a1 ... state_an , state_b1 ... state_bm).
        These ORs are then combines with ANDs.

        If a component does not carry states, this will be a ComponentStateTarget.
        """
        component_state_targets = []  # type: List[ComponentStateTarget]
        component_to_factor = {}  # type: Dict[Spec, VennSet[StateTarget]]
        for component in rxncon_sys.components():
            grouped_states = rxncon_sys.states_for_component_grouped(component)
            # component is not part of any state
            if not grouped_states.values():
                component_state_targets.append(ComponentStateTarget(component))
                component_to_factor[component] = ValueSet(ComponentStateTarget(component))
            # component is part of at least one state
            else:
                # mutually exclusive states are combined by OR

                component_to_factor[component] = \
                    Intersection(
                        *(Union(*(ValueSet(StateTarget(x)) for x in group)) for group in grouped_states.values()))

        return component_to_factor, component_state_targets

    def calc_reaction_targets_with_dnf_contingencies(k_plus_strict: bool, k_minus_strict: bool) -> List[ReactionTarget]:
        """Calculates contingency factors for reaction targets.
        Degradation reactions are handled differently then other reactions. An OR contingency will lead to a
        split of the degradation reaction in as many reactions as OR statements. Each OR will be assigned to one
        instance of the reaction."""
        reaction_targets = set()

        for reaction in rxncon_sys.reactions:
            factors = (x.to_venn_set(k_plus_strict=k_plus_strict, k_minus_strict=k_minus_strict, structured=False,
                                     state_wrapper=StateTarget)
                       for x in rxncon_sys.contingencies_for_reaction(reaction))
            cont = Intersection(*factors).to_simplified_set()  # type: VennSet[StateTarget]
            # The reaction is not a degradation reaction or the DNF has just one term.
            if not reaction.degraded_components or len(cont.to_dnf_list()) == 1:
                reaction_targets.add(ReactionTarget(reaction, contingency_factor=cont))
            # The reaction is a degradation reaction
            else:
                # The reaction is split into separated entities according to the number of minterms of the
                # disjunctive normal form (dnf). Each minterm will be assigned to a entity of the degradation reaction.
                for index, factor in enumerate(cont.to_dnf_list()):
                    reaction_targets.add(
                        ReactionTarget(reaction, contingency_variant=index, contingency_factor=factor))

        return list(reaction_targets)

    def update_degs_add_component_states(reaction_targets: List[ReactionTarget],
                                         component_state_targets: List[ComponentStateTarget]) -> List[ReactionTarget]:
        """For degradation reactions, add the stateless components they degrade to the list of targets they degrade."""
        result = deepcopy(reaction_targets)
        for reaction_target in result:
            for degraded_component in reaction_target.degraded_components:
                if ComponentStateTarget(degraded_component) in component_state_targets:
                    reaction_target.degraded_targets.append(ComponentStateTarget(degraded_component))

        return result

    def update_degs_add_contingent_states(reaction_targets: List[ReactionTarget]) -> List[ReactionTarget]:
        """For degradation reactions, add their contingent states to the list of targets they degrade."""
        def degraded_state_targets(component: Spec, soln: Dict[StateTarget, bool]) -> List[StateTarget]:
            # Disregard input states, since they do not influence which states are degraded.
            soln = {k: v for k, v in soln.items() if not k.is_input() and component in k.components}

            # soln evaluates to False if solution is tautology, since when there are no constraints on which
            # states are required to be true/false, soln is an empty dict. Nicely counterintuitive.
            if not soln and ComponentStateTarget(component) in component_state_targets:
                return [ComponentStateTarget(component)]
            elif not soln:
                return [StateTarget(x) for x in rxncon_sys.states_for_component(component)]
            else:
                trues = [target for target, val in soln.items() if val]
                falses = [target for target, val in soln.items() if not val]
                for target in falses:
                    trues += target.complementary_state_targets(rxncon_sys, component)

                return trues

        result = deepcopy(reaction_targets)

        for reaction_target in result:
            solutions = reaction_target.contingency_factor.calc_solutions()

            if reaction_target.degraded_targets and len(solutions) == 1 and not solutions[0]:
                # No contingencies, but targeted degradation. Do not mess with the list of degraded targets.
                continue

            for degraded_component, solution in product(reaction_target.degraded_components, solutions):  # type: ignore
                reaction_target.degraded_targets.extend(degraded_state_targets(degraded_component, solution))

        return result

    def update_degs_add_interaction_state_partner(reaction_targets: List[ReactionTarget]) -> List[ReactionTarget]:
        """Update degradation reactions for interaction states.
        Interaction states are composed out of two components. A degradation reaction degrading an interaction
        state will degrade one of these components. The other component is then released, and thus produced
        by the degrading reaction."""
        result = []

        for reaction_target in reaction_targets:
            appended = False
            degraded_interaction_targets = [x for x in reaction_target.degraded_targets if x.is_interaction]
            for index, interaction_target in enumerate(degraded_interaction_targets):
                empty_partners = [neutral_target for neutral_target in interaction_target.neutral_targets
                                  if not any(component in reaction_target.degraded_components
                                             for component in neutral_target.components)]

                if interaction_target.is_homodimer:
                    assert len(empty_partners) == 0, 'update_degs_add_interaction_state_partner::homodimer error, ' \
                                                     'reaction_target: {}, interaction_target: {}, empty_partners: {}' \
                                                     ''.format(reaction_target, interaction_target, empty_partners)
                    continue

                if len(empty_partners) != 1:
                    raise AssertionError('update_degs_add_interaction_state_partner::empty partners != 1 error\n'
                                         'The full list of degraded targets is {}\n'
                                         'The current reaction target is {}\n'
                                         'The current interaction target is {}\n'
                                         'The current empty partners that have been deduced are: {}\n'
                                         .format(', '.join(str(x) for x in degraded_interaction_targets),
                                                 str(reaction_target), str(interaction_target),
                                                 ', '.join(str(x) for x in empty_partners)))
                new_reaction = deepcopy(reaction_target)
                new_reaction.interaction_variant_index = index if len(degraded_interaction_targets) > 1 else None
                new_reaction.consumed_targets.append(interaction_target)
                new_reaction.produced_targets.append(empty_partners[0])
                result.append(new_reaction)
                appended = True

            if not appended:
                result.append(deepcopy(reaction_target))

        return result

    def update_syns_with_component_states(reaction_targets: List[ReactionTarget],
                                          component_state_targets: List[ComponentStateTarget]) -> List[ReactionTarget]:
        """Update synthesis reaction with component states: stateless components that are synthesised have
        rights too."""
        result = deepcopy(reaction_targets)

        for reaction_target in result:
            for component in reaction_target.synthesised_components:
                if ComponentStateTarget(component) in component_state_targets:
                    reaction_target.synthesised_targets.append(ComponentStateTarget(component))

        return result

    def calc_knockout_targets(knockout_strategy: KnockoutStrategy) -> List[KnockoutTarget]:
        if knockout_strategy == KnockoutStrategy.no_knockout:
            return []
        else:
            return [KnockoutTarget(component) for component in rxncon_sys.components()]

    def calc_overexpression_targets(overexpression_strategy: OverexpressionStrategy) -> List[OverexpressionTarget]:
        if overexpression_strategy == OverexpressionStrategy.no_overexpression:
            return []
        else:
            return [OverexpressionTarget(component) for component in rxncon_sys.components()]

    def calc_reaction_rules() -> None:
        """Calculate the rules of reaction targets: the component factor AND the contingencies."""

        for reaction_target in reaction_targets:
            components = (component_presence_factor[x] for x in reaction_target.components_lhs)
            component_factor = Intersection(*components)  # type: VennSet[StateTarget]
            reaction_rules.append(UpdateRule(reaction_target, Intersection(
                component_factor, reaction_target.contingency_factor).to_simplified_set()))

    def calc_state_rules() -> None:
        """Calculate the rules of the state targets, includes smoothing. For details see our paper."""
        def synthesis_factor(state_target: StateTarget) -> VennSet[Target]:
            fac = EmptySet()  # type: VennSet[Target]
            for rxn in (x for x in reaction_targets if x.synthesises(state_target)):
                fac = Union(fac, ValueSet(rxn))

            for prod_rxn in (x for x in reaction_targets if x.produces(state_target)):
                sources = []
                for source in prod_rxn.consumed_targets:
                    sources.append([source] + [x for x in reaction_targets if x.synthesises(source)])

                for source_combi in product(*sources):
                    # At least one source should be synthesised.
                    if all(isinstance(x, StateTarget) for x in source_combi):
                        continue
                    assert any(isinstance(x, ReactionTarget) and x.synthesised_targets for x in source_combi)

                    fac = Union(fac, Intersection(ValueSet(prod_rxn), *(ValueSet(x) for x in source_combi)))

            return fac

        def component_factor(state_target: StateTarget) -> VennSet[StateTarget]:
            return Intersection(*(component_presence_factor[x] for x in state_target.components))

        def degradation_factor(state_target: StateTarget) -> VennSet[ReactionTarget]:
            return Complement(Union(*(ValueSet(x) for x in reaction_targets if x.degrades(state_target))))

        def pi(state_target: StateTarget, level: int) -> VennSet[Target]:
            res = EmptySet()  # type: VennSet[Target]

            for r in (x for x in reaction_targets if x.produces(state_target)):
                rxn_term = ValueSet(r)  # type: VennSet[Target]
                for s in (x for x in state_targets if r.consumes(x)):
                    if r.degraded_targets:
                        state_term = ValueSet(s)  # type: VennSet[Target]
                    else:
                        state_term = Intersection(ValueSet(s), degradation_factor(s))
                    for l in range(level):
                        state_term = Union(state_term, sigma(s, level - 1))
                    rxn_term = Intersection(rxn_term, state_term)
                res = Union(res, rxn_term)

            return res

        def kappa(state_target: StateTarget, level: int) -> VennSet[Target]:
            res = EmptySet()  # type: VennSet[Target]

            for r in (x for x in reaction_targets if x.consumes(state_target)):
                rxn_term = ValueSet(r)  # type: VennSet[Target]
                for s in (x for x in state_targets if r.consumes(x)):
                    rxn_term = Intersection(rxn_term, ValueSet(s), degradation_factor(s))
                res = Union(res, rxn_term)

            return res

        def sigma(state_target: StateTarget, level: int) -> VennSet[Target]:
            prod_cons_factor = Union(pi(state_target, level),
                                     Intersection(ValueSet(state_target), Complement(kappa(state_target, level))))

            return Union(synthesis_factor(state_target),
                         Intersection(degradation_factor(state_target),
                                      component_factor(state_target),
                                      prod_cons_factor))

        if smoothing_strategy == SmoothingStrategy.no_smoothing:
            level = 0
        elif smoothing_strategy == SmoothingStrategy.smooth_production_sources:
            level = 1
        else:
            raise AssertionError

        for state_target in state_targets:
            state_rules.append(UpdateRule(state_target, sigma(state_target, level).to_simplified_set()))

    def update_state_rules_with_knockouts(knockout_strategy: KnockoutStrategy) -> None:
        if knockout_strategy == KnockoutStrategy.no_knockout:
            return
        elif knockout_strategy in (KnockoutStrategy.knockout_all_states, KnockoutStrategy.knockout_neutral_states):
            for state_rule in state_rules:
                assert isinstance(state_rule.target, StateTarget)

                if knockout_strategy == KnockoutStrategy.knockout_neutral_states and not state_rule.target.is_neutral:
                    continue

                knockout_factor = Complement(
                    Union(*(ValueSet(KnockoutTarget(component)) for component in state_rule.target.components)))
                state_rule.factor = Intersection(knockout_factor, state_rule.factor)

    def update_state_rules_with_overexpressions(overexpression_strategy: OverexpressionStrategy) -> None:
        if overexpression_strategy == OverexpressionStrategy.no_overexpression:
            return
        elif overexpression_strategy in (OverexpressionStrategy.overexpress_all_states,
                                         OverexpressionStrategy.overexpress_neutral_states):
            for state_rule in state_rules:
                assert isinstance(state_rule.target, StateTarget)

                if overexpression_strategy == OverexpressionStrategy.overexpress_neutral_states and not state_rule.target.is_neutral:
                    continue

                overexpression_factor = Intersection(
                    *(ValueSet(OverexpressionTarget(component)) for component in state_rule.target.components))
                state_rule.factor = Union(overexpression_factor, state_rule.factor)

    def calc_knockout_rules() -> None:
        for knockout_target in knockout_targets:
            knockout_rules.append(UpdateRule(knockout_target, ValueSet(knockout_target)))

    def calc_overexpression_rules() -> None:
        for overexpression_target in overexpression_targets:
            overexpression_rules.append(UpdateRule(overexpression_target, ValueSet(overexpression_target)))

    def update_input_output_rules() -> None:
        """If an Input state and an Output reaction share the same name [BLA], they are assumed
        to refer to the same global quantity. Therefore the update rule for the state (which was trivial),
        becomes the update rule for the reaction."""
        to_delete = []  # type: List[UpdateRule]
        for reaction_rule in reaction_rules:
            for state_rule in state_rules:
                if (reaction_rule.target.is_output and state_rule.target.is_input and   # type: ignore
                        str(reaction_rule.target) == str(state_rule.target)):
                    state_rule.factor = reaction_rule.factor
                    to_delete.append(reaction_rule)

        for rule_to_delete in to_delete:
            reaction_targets.remove(rule_to_delete.target)
            reaction_rules.remove(rule_to_delete)

    component_presence_factor, component_state_targets = calc_component_presence_factors()

    state_targets = [StateTarget(x) for x in rxncon_sys.states]  # type: List[StateTarget]
    state_targets += component_state_targets

    reaction_targets = calc_reaction_targets_with_dnf_contingencies(k_plus_strict, k_minus_strict)
    reaction_targets = update_degs_add_component_states(reaction_targets, component_state_targets)
    reaction_targets = update_degs_add_contingent_states(reaction_targets)
    reaction_targets = update_degs_add_interaction_state_partner(reaction_targets)
    reaction_targets = update_syns_with_component_states(reaction_targets, component_state_targets)

    knockout_targets = calc_knockout_targets(knockout_strategy)
    overexpression_targets = calc_overexpression_targets(overexpression_strategy)

    reaction_rules = []  # type: List[UpdateRule]
    state_rules = []  # type: List[UpdateRule]
    knockout_rules = []  # type: List[UpdateRule]
    overexpression_rules = []  # type: List[UpdateRule]

    calc_reaction_rules()
    calc_state_rules()
    update_state_rules_with_knockouts(knockout_strategy)
    update_state_rules_with_overexpressions(overexpression_strategy)
    calc_knockout_rules()
    calc_overexpression_rules()
    update_input_output_rules()

    return BooleanModel(state_targets + reaction_targets + knockout_targets + overexpression_targets,  # type: ignore
                        reaction_rules + state_rules + knockout_rules + overexpression_rules,
                        initial_conditions(reaction_targets, state_targets, knockout_targets, overexpression_targets)) class Set(Generic[T]):
    def calc_solutions(self) -> List[Dict[T, bool]]:
        val_to_sym = self._make_val_to_sym_dict()
        sym_to_val = {sym: val for val, sym in val_to_sym.items()}

        venn_solns = []
        for s in self._to_pyeda_expr(val_to_sym).satisfy_all():
            venn_solns.append({sym_to_val[sym.name]: bool(truth) for sym, truth in s.items()})

        return venn_solns

    def eval_boolean_func(self, vars: Dict[T, bool]) -> bool:
        val_to_sym = self._make_val_to_sym_dict()
        try:
            evaluated = self._to_pyeda_expr(val_to_sym).restrict({expr(sym): vars[val] for val, sym in val_to_sym.items()})
        except KeyError as e:
            raise AssertionError('eval_boolean_func missing variable {}'.format(e.args[0]))

        if evaluated.is_one():
            return True
        elif evaluated.is_zero():
            return False
        else:
            raise AssertionError

    def to_simplified_set(self) -> 'Set[T]':
        val_to_sym = self._make_val_to_sym_dict()
        sym_to_val = {sym: val for val, sym in val_to_sym.items()}
        return venn_from_pyeda(self._to_pyeda_expr(val_to_sym).simplify(), sym_to_val)

    def to_dnf_set(self) -> 'Set[T]':
        val_to_sym = self._make_val_to_sym_dict()
        sym_to_val = {sym: val for val, sym in val_to_sym.items()}
        return venn_from_pyeda(self._to_pyeda_expr(val_to_sym).to_dnf(), sym_to_val)

    def to_dnf_list(self) -> List['Set[T]']:
        val_to_sym = self._make_val_to_sym_dict()
        sym_to_val = {sym: val for val, sym in val_to_sym.items()}
        dnf_set = self._to_pyeda_expr(val_to_sym).to_dnf()

        if dnf_set is One:
            return [UniversalSet()]
        elif dnf_set is Zero:
            return [EmptySet()]
        elif isinstance(dnf_set, Literal):
            return [venn_from_pyeda(dnf_set, sym_to_val)]
        elif isinstance(dnf_set, AndOp):
            return [venn_from_pyeda(dnf_set, sym_to_val)]
        elif isinstance(dnf_set, OrOp):
            return [venn_from_pyeda(term, sym_to_val) for term in dnf_set.xs]
        else:
            raise Exception

    def to_dnf_nested_list(self) -> List[List['Set[T]']]:
        val_to_sym = self._make_val_to_sym_dict()
        sym_to_val = {sym: val for val, sym in val_to_sym.items()}
        dnf_set = self._to_pyeda_expr(val_to_sym).to_dnf()

        if dnf_set is One:
            return [[UniversalSet()]]
        elif dnf_set is Zero:
            return [[EmptySet()]]
        elif isinstance(dnf_set, Literal):
            return [[venn_from_pyeda(dnf_set, sym_to_val)]]
        elif isinstance(dnf_set, AndOp):
            return [[venn_from_pyeda(x, sym_to_val) for x in dnf_set.xs]]

        res = []  # type: List[List[Set[T]]]
        for term in dnf_set.xs:
            if isinstance(term, Literal):
                res.append([venn_from_pyeda(term, sym_to_val)])
            elif isinstance(term, AndOp):
                res.append([venn_from_pyeda(x, sym_to_val) for x in term.xs])
            else:
                raise Exception
        return res

    def is_equivalent_to(self, other: 'Set[T]') -> bool:
        val_to_sym = self._make_val_to_sym_dict()
        val_to_sym = other._make_val_to_sym_dict(val_to_sym)
        return self._to_pyeda_expr(val_to_sym).equivalent(other._to_pyeda_expr(val_to_sym))

    def is_subset_of(self, other: 'Set[T]') -> bool:
        val_to_sym = self._make_val_to_sym_dict()
        val_to_sym = other._make_val_to_sym_dict(val_to_sym)
        return Implies(self._to_pyeda_expr(val_to_sym), other._to_pyeda_expr(val_to_sym)).equivalent(expr(1))

    def is_superset_of(self, other: 'Set[T]') -> bool:
        return other.is_subset_of(self)

    @property
    def values(self) -> List[T]:
        return []

    def _to_pyeda_expr(self, val_to_sym: MutableMapping[T, str]) -> Expression:
        return None

    def _make_val_to_sym_dict(self, existing_dict: Optional[MutableMapping[T, str]]=None) -> MutableMapping[T, str]:
        vals = []  # type: List[T]

        for v in self.values:
            found = False
            for existing in vals:
                try:
                    if v == existing:
                        found = True
                        break
                except AttributeError:
                    pass

            if not found:
                vals.append(v)

        if existing_dict:
            d = existing_dict
            first_new_sym_idx = SYMS.index(next(reversed(d.values()))) + 1  # type: ignore
            for i, val in enumerate(x for x in vals if x not in d.keys()):
                d[val] = SYMS[i + first_new_sym_idx]
        else:
            d = OrderedDict()
            for i, val in enumerate(vals):
                d[val] = SYMS[i]

        return d class ValueSet(UnarySet[T_inv], Generic[T_inv]):
    def __init__(self, value: T_inv) -> None:
        assert value is not None
        assert isinstance(hash(value), int)
        self.value = value

    def __eq__(self, other: object) -> bool:
        if not isinstance(other, Set):
            return NotImplemented
        elif isinstance(other, ValueSet):
            return self.value == other.value
        else:
            return False

    def __hash__(self) -> int:
        return hash('*value-set-{}*'.format(hash(self.value)))

    def __repr__(self) -> str:
        return str(self)

    def __str__(self) -> str:
        if self.value:
            return '{}'.format(str(self.value))
        else:
            return 'UniversalSet'

    @property
    def values(self) -> List[T_inv]:
        return [self.value]

    def _to_pyeda_expr(self, val_to_sym: MutableMapping[Any, str]) -> Expression:
        return expr(val_to_sym[self.value]) class Complement(UnarySet[T], Generic[T]):
    def __init__(self, expr: Set[T]) -> None:
        self.expr = expr

    def __eq__(self, other: object) -> bool:
        if not isinstance(other, Set):
            return NotImplemented
        elif isinstance(other, Complement):
            return self.expr == other.expr
        else:
            return False

    def __hash__(self) -> int:
        return hash('*complement-{}*'.format(self.expr))

    def __repr__(self) -> str:
        return str(self)

    def __str__(self) -> str:
        return '!({})'.format(str(self.expr))

    @property
    def values(self) -> List[T]:
        return self.expr.values

    def _to_pyeda_expr(self, val_to_sym: MutableMapping[T, str]) -> Expression:
        return Not(self.expr._to_pyeda_expr(val_to_sym)) class Intersection(NarySet[T], Generic[T]):
    def __new__(cls, *exprs: Set[T], **kwargs: Any) -> Set[T]:
        if len(exprs) == 0:
            return UniversalSet()
        else:
            return NarySet.__new__(cls, *exprs)  # type: ignore

    def __deepcopy__(self, memodict: Dict) -> Set[T]:
        return Intersection(*(deepcopy(x) for x in self.exprs))

    def __eq__(self, other: object) -> bool:
        if not isinstance(other, Set):
            return NotImplemented
        elif isinstance(other, Intersection):
            return all(my_expr == other_expr for my_expr, other_expr in zip(self.exprs, other.exprs))
        else:
            return False

    def __hash__(self) -> int:
        return hash(str(self))

    def __repr__(self) -> str:
        return str(self)

    def __str__(self) -> str:
        return '({})'.format(' & '.join(str(expr) for expr in self.exprs))

    def _to_pyeda_expr(self, val_to_sym: MutableMapping[T, str]) -> Expression:
        return And(*(expr._to_pyeda_expr(val_to_sym) for expr in self.exprs)) class Union(NarySet[T], Generic[T]):
    def __new__(cls, *exprs: Set[T], **kwargs: Any) -> Set[T]:
        if len(exprs) == 0:
            return EmptySet()
        else:
            return NarySet.__new__(cls, *exprs)  # type: ignore

    def __deepcopy__(self, memodict: Dict) -> Set[T]:
        return Union(*(deepcopy(x) for x in self.exprs))

    def __eq__(self, other: object) -> bool:
        if not isinstance(other, Set):
            return NotImplemented
        elif isinstance(other, Union):
            return all(my_expr == other_expr for my_expr, other_expr in zip(self.exprs, other.exprs))
        else:
            return False

    def __hash__(self) -> int:
        return hash(str(self))

    def __repr__(self) -> str:
        return str(self)

    def __str__(self) -> str:
        return '({})'.format(' | '.join(str(expr) for expr in self.exprs))

    def _to_pyeda_expr(self, val_to_sym: MutableMapping[T, str]) -> Expression:
        return Or(*(expr._to_pyeda_expr(val_to_sym) for expr in self.exprs)) class EmptySet(Set[Any]):
    def __eq__(self, other: object) -> bool:
        if not isinstance(other, Set):
            return NotImplemented
        return isinstance(other, EmptySet)

    def __hash__(self) -> int:
        return hash('*empty-set*')

    def __repr__(self) -> str:
        return str(self)

    def __str__(self) -> str:
        return 'EmptySet'

    def _to_pyeda_expr(self, val_to_sym: MutableMapping[Any, str]) -> Expression:
        return expr(0) class UniversalSet(Set[Any]):
    def __init__(self) -> None:
        pass

    def __eq__(self, other: object) -> bool:
        if not isinstance(other, Set):
            return NotImplemented
        return isinstance(other, UniversalSet)

    def __hash__(self) -> int:
        return hash('UniversalSet')

    def __repr__(self) -> str:
        return 'UniversalSet'

    def __str__(self) -> str:
        return 'UniversalSet'

    def _to_pyeda_expr(self, val_to_sym: MutableMapping[Any, str]) -> Expression:
        return expr(1)
from enum import Enum
from typing import Tuple, Dict
from rxncon.core.rxncon_system import RxnConSystem
from rxncon.simulation.boolean.boolean_model import BooleanModel, Target, ReactionTarget, KnockoutTarget, \
    OverexpressionTarget, StateTarget, UpdateRule, SmoothingStrategy, KnockoutStrategy, OverexpressionStrategy, \
    boolean_model_from_rxncon
from rxncon.venntastic.sets import Set as VennSet, ValueSet, Complement, Intersection, Union, EmptySet, UniversalSet
"""Module containing the functions boolnet_from_boolean_model, boolnet_strs_from_rxncon and
the class QuantitativeContingencyStrategy."""




def boolnet_from_boolean_model(boolean_model: BooleanModel) -> Tuple[str, Dict[str, str], Dict[str, bool]]:
    """Translates the boolean model into BoolNet syntax.

    Returns:
        1. The boolean model in BoolNet syntax,
        2. The (BoolNet name, rxncon name) mapping, since BoolNet is picky about characters,
        3. The initial conditions."""

    def initialize_boolnet_names() -> None:
        nonlocal boolnet_names
        nonlocal reaction_index
        nonlocal state_index
        nonlocal knockout_index
        nonlocal overexpression_index
        for update_rule in boolean_model.update_rules:
            target = update_rule.target
            if target in boolnet_names.keys():
                raise AssertionError
            elif isinstance(target, ReactionTarget):
                name = 'R{}'.format(reaction_index)
                boolnet_names[target] = name
                reaction_index += 1
            elif isinstance(target, KnockoutTarget):
                name = 'K{}'.format(knockout_index)
                boolnet_names[target] = name
                knockout_index += 1
            elif isinstance(target, OverexpressionTarget):
                name = 'O{}'.format(overexpression_index)
                boolnet_names[target] = name
                overexpression_index += 1
            elif isinstance(target, StateTarget):
                name = 'S{}'.format(state_index)
                boolnet_names[target] = name
                state_index += 1
            else:
                raise AssertionError

    def str_from_factor(factor: VennSet) -> str:
        if isinstance(factor, ValueSet):
            return boolnet_names[factor.value]
        elif isinstance(factor, Complement):
            return '!({})'.format(str_from_factor(factor.expr))
        elif isinstance(factor, Intersection):
            return '({})'.format(' & '.join(str_from_factor(x) for x in factor.exprs))
        elif isinstance(factor, Union):
            return '({})'.format(' | '.join(str_from_factor(x) for x in factor.exprs))
        elif isinstance(factor, EmptySet):
            return '0'
        elif isinstance(factor, UniversalSet):
            return '1'
        else:
            raise AssertionError('Could not parse factor {}'.format(factor))

    def str_from_update_rule(update_rule: UpdateRule) -> str:
        return '{0}, {1}'.format(boolnet_names[update_rule.target],
                                 str_from_factor(update_rule.factor))

    boolnet_names = {}  # type: Dict[Target, str]

    reaction_index = 0
    state_index = 0
    knockout_index = 0
    overexpression_index = 0

    initialize_boolnet_names()

    def sort_key(rule_str: str) -> Tuple[str, int]:
        target = rule_str.split(',')[0].strip()
        return target[0], int(target[1:])

    rule_strs = sorted([str_from_update_rule(x) for x in boolean_model.update_rules], key=sort_key)

    return 'targets, factors\n' + '\n'.join(rule for rule in rule_strs) + '\n', \
           {name: str(target) for target, name in boolnet_names.items()}, \
           {boolnet_names[target]: value for target, value in boolean_model.initial_conditions.target_to_value.items()}


class QuantitativeContingencyStrategy(Enum):
    """Strategy for dealing with quantitative contingencies k+ and k-, ignore them or make them into ! resp. x."""
    strict = 'strict'
    ignore = 'ignore'



def boolnet_strs_from_rxncon(rxncon: RxnConSystem, smoothing_strategy: SmoothingStrategy,