diff --git a/baybe/acquisition/_builder.py b/baybe/acquisition/_builder.py
index 8c92185752..61979d492b 100644
--- a/baybe/acquisition/_builder.py
+++ b/baybe/acquisition/_builder.py
@@ -23,6 +23,7 @@
     _ExpectedHypervolumeImprovement,
     qExpectedHypervolumeImprovement,
     qLogExpectedHypervolumeImprovement,
+    qLogNoisyExpectedImprovement,
     qNegIntegratedPosteriorVariance,
     qThompsonSampling,
 )
@@ -75,6 +76,7 @@ class BotorchAcquisitionArgs:
     # Optional, depending on the specific acquisition function being used
     best_f: float | None = _OPT_FIELD
     beta: float | None = _OPT_FIELD
+    constraints: list | None = _OPT_FIELD
     maximize: bool | None = _OPT_FIELD
     mc_points: Tensor | None = _OPT_FIELD
     num_fantasies: int | None = _OPT_FIELD
@@ -197,6 +199,7 @@ def build(self) -> BoAcquisitionFunction:
         # Set context-specific parameters
         self._set_best_f()
         self._set_target_transformation()
+        self._set_constraints()
         self._set_X_baseline()
         self._set_X_pending()
         self._set_mc_points()
@@ -222,6 +225,18 @@ def _set_target_transformation(self) -> None:
                 return
 
         if self.acqf.is_analytic:
+            # TODO: Certain analytic acquisition functions (e.g. analytic EI with
+            #   constraints) do support outcome constraints and will be added to BayBE
+            #   in the future. Once available, this guard should be scoped to only
+            #   those analytic acqfs that do NOT support constraints, and
+            #   `to_botorch_posterior_transform()` must be fixed to pad the weight
+            #   vector to length `n_models`.
+            if self.objective.outcome_constraints:
+                raise IncompatibilityError(
+                    f"Analytical acquisition function '{type(self.acqf).__name__}' "
+                    f"does not support outcome constraints. Use an MC-based "
+                    f"acquisition function instead."
+                )
             try:
                 transform = self.objective.to_botorch_posterior_transform()
             except NonGaussianityError as ex:
@@ -253,17 +268,111 @@ def _set_target_transformation(self) -> None:
 
             self._args.objective = self.objective.to_botorch()
 
+    def _set_constraints(self) -> None:
+        """Set BoTorch's ``constraints`` argument from outcome constraints.
+
+        Outcome constraint compatibility check — Layer 2 (acquisition function level).
+        Raises IncompatibilityError if the acqf's BoTorch __init__ signature does not
+        include a ``constraints`` parameter.
+        """
+        if not self.objective.outcome_constraints:
+            return
+
+        if flds.constraints.name not in self._signature:
+            raise IncompatibilityError(
+                f"The selected acquisition function "
+                f"'{type(self.acqf).__name__}' does not support outcome "
+                f"constraints. Use a compatible acquisition function such as "
+                f"'{qLogNoisyExpectedImprovement.__name__}' instead."
+            )
+        constraints = self.objective.to_botorch_constraints()
+        if constraints:
+            self._args.constraints = constraints
+
     def _set_best_f(self) -> None:
-        """Set BoTorch's ``best_f`` argument."""
+        """Set BoTorch's ``best_f`` argument.
+
+        best_f is a constant reference value (not differentiable). When outcome
+        constraints are present, only feasible training points are considered.
+        """
         if flds.best_f.name not in self._signature:
             return
 
         match self.objective:
             case SingleTargetObjective() | DesirabilityObjective():
-                self._args.best_f = self._posterior_mean_comp.max().item()
+                if not (constraints := self.objective.to_botorch_constraints()):
+                    self._args.best_f = self._posterior_mean_comp.max().item()
+                else:
+                    self._args.best_f = self._compute_best_f_with_constraints(
+                        constraints
+                    )
             case _:
                 raise NotImplementedError("This line should be impossible to reach.")
 
+    def _compute_best_f_with_constraints(
+        self, constraints: list[Callable[[Tensor], Tensor]]
+    ) -> float:
+        """Compute the best objective value considering outcome constraints.
+
+        Falls back to the global maximum if no feasible training point exists.
+
+        Args:
+            constraints: Constraint functions from
+                :meth:`~baybe.objectives.base.Objective.to_botorch_constraints`.
+
+        Returns:
+            The best feasible objective value, or the global maximum as fallback.
+        """
+        # Get objective values for all training points
+        objective_values = self._posterior_mean_comp
+
+        # Get raw model predictions for constraint evaluation
+        batched = to_tensor(self._train_x).unsqueeze(-2)
+        posterior = self._botorch_surrogate.posterior(batched)
+        model_predictions = posterior.mean.squeeze(-2)
+
+        # Apply constraint functions to filter feasible points
+        feasible_mask = self._compute_feasible_mask(model_predictions, constraints)
+
+        if not feasible_mask.any():
+            # TODO: other mechanisms, e.g. steer towards feasible region?
+            # No feasible training points - fall back to global maximum
+            return objective_values.max().item()
+
+        # Return maximum among feasible points
+        feasible_objectives = objective_values[feasible_mask]
+        return feasible_objectives.max().item()
+
+    def _compute_feasible_mask(
+        self,
+        model_predictions: Tensor,
+        constraints: list[Callable[[Tensor], Tensor]],
+    ) -> Tensor:
+        """Compute boolean mask indicating which points satisfy all constraints.
+
+        Uses hard thresholding (feasible when constraint value <= 0) combined
+        via boolean AND across all constraints.
+
+        Args:
+            model_predictions: Raw model predictions [n_points, n_outputs]
+            constraints: Constraint functions from to_botorch_constraints()
+
+        Returns:
+            Boolean mask [n_points] where True = feasible, False = infeasible
+        """
+        n_points = model_predictions.shape[0]
+        feasible_mask = torch.ones(n_points, dtype=torch.bool)
+
+        for constraint_func in constraints:
+            # Constraint func: [batch, q, m] -> [batch, q]; we insert q=1
+            # via unsqueeze(-2), so output is [n_points, 1]; squeeze q dim.
+            constraint_violations = constraint_func(
+                model_predictions.unsqueeze(-2)
+            ).squeeze(-1)
+            feasible_mask &= constraint_violations <= 0
+
+        return feasible_mask
+
     def set_default_sample_shape(self, acqf: BoAcquisitionFunction, /):
         """Apply temporary workaround for Thompson sampling."""
         # TODO: Needs redesign once bandits are supported more generally
diff --git a/baybe/constraints/__init__.py b/baybe/constraints/__init__.py
index 8b92ecd6fe..32644fc1aa 100644
--- a/baybe/constraints/__init__.py
+++ b/baybe/constraints/__init__.py
@@ -21,6 +21,7 @@
     DiscreteProductConstraint,
     DiscreteSumConstraint,
 )
+from baybe.constraints.outcome import OutcomeConstraint
 from baybe.constraints.validation import validate_constraints
 
 __all__ = [
@@ -42,6 +43,8 @@
     "DiscretePermutationInvarianceConstraint",
     "DiscreteProductConstraint",
     "DiscreteSumConstraint",
+    # --- Outcome constraints ---#
+    "OutcomeConstraint",
     # --- Other --- #
     "validate_constraints",
     "DISCRETE_CONSTRAINTS_FILTERING_ORDER",
diff --git a/baybe/constraints/outcome.py b/baybe/constraints/outcome.py
new file mode 100644
index 0000000000..2f5ec1b0a2
--- /dev/null
+++ b/baybe/constraints/outcome.py
@@ -0,0 +1,89 @@
+"""Functionality for outcome constraints."""
+
+from __future__ import annotations
+
+from collections.abc import Callable
+from typing import Literal
+
+import pandas as pd
+import torch
+from attrs import define, field
+from attrs.validators import in_, instance_of
+
+from baybe.serialization.mixin import SerialMixin
+from baybe.targets.base import Target
+
+
+@define(frozen=True, slots=False)
+class OutcomeConstraint(SerialMixin):
+    """A constraint applied to target outcomes in the output space.
+
+    Outcome constraints restrict the feasible region based on target predictions,
+    different from parameter constraints which restrict the input space.
+    """
+
+    target: Target = field(validator=instance_of(Target))
+    """The target to be constrained."""
+
+    operator: Literal["<=", ">=", "=="] = field(validator=in_(["<=", ">=", "=="]))
+    """The constraint operator."""
+
+    threshold: float = field(validator=instance_of((int, float)), converter=float)
+    """The constraint threshold value in experimental units."""
+
+    def __str__(self) -> str:
+        """Return string representation."""
+        return f"{self.target.name} {self.operator} {self.threshold}"
+
+    def get_computational_threshold(self) -> float:
+        """Convert experimental threshold to computational units.
+
+        Returns:
+            The threshold value in computational units.
+        """
+        # Create dummy series with threshold value in experimental units
+        experimental_series = pd.Series([self.threshold], name=self.target.name)
+
+        # Apply the same transformations as the target
+        computational_series = self.target.transform(experimental_series)
+
+        return computational_series.iloc[0]
+
+    def to_botorch_constraint_func(
+        self, target_idx: int
+    ) -> Callable[[torch.Tensor], torch.Tensor]:
+        """Create a botorch-compatible constraint function.
+
+        Args:
+            target_idx: Index of the target in model output.
+
+        Returns:
+            A constraint function that returns <= 0 for feasible region.
+        """
+        computational_threshold = self.get_computational_threshold()
+
+        def constraint_func(samples: torch.Tensor) -> torch.Tensor:
+            """Constraint function operating on computational level.
+
+            Args:
+                samples: Model output samples in computational units.
+
+            Returns:
+                Constraint values where <= 0 indicates feasible region.
+
+            Raises:
+                ValueError: If the constraint operator is not supported.
+            """
+            if self.operator == "<=":
+                return samples[..., target_idx] - computational_threshold
+            elif self.operator == ">=":
+                return computational_threshold - samples[..., target_idx]
+            elif self.operator == "==":
+                # Equality constraint with small tolerance
+                return (
+                    torch.abs(samples[..., target_idx] - computational_threshold) - 1e-6
+                )
+            else:
+                raise ValueError(f"Unsupported constraint operator: {self.operator}")
+
+        return constraint_func
diff --git a/baybe/exceptions.py b/baybe/exceptions.py
index 0be2273341..afbc7b291c 100644
--- a/baybe/exceptions.py
+++ b/baybe/exceptions.py
@@ -44,6 +44,10 @@ class MinimumCardinalityViolatedWarning(UserWarning):
     """Minimum cardinality constraints are violated."""
 
 
+class OutcomeConstraintIgnoredWarning(UserWarning):
+    """Outcome constraints are present but cannot be enforced by the recommender."""
+
+
 ##### Exceptions #####
 
 
diff --git a/baybe/objectives/base.py b/baybe/objectives/base.py
index eee36ee007..f3a5955534 100644
--- a/baybe/objectives/base.py
+++ b/baybe/objectives/base.py
@@ -5,15 +5,18 @@
 import gc
 import warnings
 from abc import ABC, abstractmethod
+from collections.abc import Callable
 from typing import TYPE_CHECKING, ClassVar
 
 import pandas as pd
 from attrs import define, field
+from attrs.validators import deep_iterable, instance_of
 
+from baybe.constraints.outcome import OutcomeConstraint
 from baybe.serialization.mixin import SerialMixin
 from baybe.targets.base import Target
 from baybe.targets.numerical import NumericalTarget
-from baybe.utils.basic import is_all_instance
+from baybe.utils.basic import is_all_instance, to_tuple
 from baybe.utils.dataframe import get_transform_objects, to_tensor
 from baybe.utils.dataframe import (
     handle_missing_values as df_handle_missing_values,
@@ -22,6 +25,7 @@
 from baybe.utils.validation import validate_target_input
 
 if TYPE_CHECKING:
+    import torch
     from botorch.acquisition.objective import MCAcquisitionObjective, PosteriorTransform
 
 
@@ -43,15 +47,56 @@ class Objective(ABC, SerialMixin):
     )
     """Optional metadata containing description and other information."""
 
+    outcome_constraints: tuple[OutcomeConstraint, ...] = field(
+        default=(),
+        converter=to_tuple,
+        validator=deep_iterable(member_validator=instance_of(OutcomeConstraint)),
+        kw_only=True,
+    )
+    """Outcome constraints applied to the optimization problem."""
+
     @property
     def description(self) -> str | None:
         """The description of the objective."""
         return self.metadata.description
 
+    def __attrs_post_init__(self) -> None:
+        """Validate outcome constraints against optimization targets."""
+        # Disallow overlap between optimization and constraint targets
+        optimization_names = {t.name for t in self._optimization_targets}
+        constraint_names = {t.name for t in self.constraint_targets}
+        if overlap := optimization_names & constraint_names:
+            raise ValueError(
+                f"Targets cannot be both optimized and constrained: {sorted(overlap)}."
+            )
+
+        # Constraint-only targets must have minimize=None
+        if invalid := [
+            t.name
+            for t in self.constraint_targets
+            if isinstance(t, NumericalTarget) and t.minimize is not None
+        ]:
+            raise ValueError(
+                f"Constraint-only targets must have minimize=None, "
+                f"but the following do not: {invalid}."
+            )
+
     @property
     @abstractmethod
+    def _optimization_targets(self) -> tuple[Target, ...]:
+        """The targets being optimized."""
+
+    @property
+    def constraint_targets(self) -> tuple[Target, ...]:
+        """Targets referenced in constraints (guaranteed disjoint from optimization)."""
+        return tuple(
+            {c.target.name: c.target for c in self.outcome_constraints}.values()
+        )
+
+    @property
     def targets(self) -> tuple[Target, ...]:
-        """The targets included in the objective."""
+        """All targets requiring data processing and modeling (opt + constraints)."""
+        return (*self._optimization_targets, *self.constraint_targets)
 
     @property
     def _modeled_quantities(self) -> tuple[Target, ...]:
@@ -82,9 +127,9 @@ def _is_multi_model(self) -> bool:
         return self._n_models > 1
 
     @property
-    @abstractmethod
     def output_names(self) -> tuple[str, ...]:
         """The names of the outputs of the objective."""
+        return tuple(t.name for t in self._optimization_targets)
 
     @property
     def n_outputs(self) -> int:
@@ -98,10 +143,10 @@ def supports_partial_measurements(self) -> bool:
 
     @property
     def _oriented_targets(self) -> tuple[Target, ...]:
-        """The targets with optional negation transformation for minimization."""
+        """The targets to optimize with optional negation transform for minimization."""
         return tuple(
             t.negate() if isinstance(t, NumericalTarget) and t.minimize else t
-            for t in self.targets
+            for t in self._optimization_targets
         )
 
     @property
@@ -149,6 +194,23 @@ def to_botorch(self) -> MCAcquisitionObjective:
             )
         )
 
+    def _check_posterior_transform_constraint_support(self) -> None:
+        """Raise error if outcome constraints are present.
+
+        Analytic acquisition functions do not yet support outcome constraints.
+        This guard prevents a shape mismatch in the posterior transform by
+        raising early when constraints are detected.
+
+        Raises:
+            NotImplementedError: If outcome constraints are present.
+        """
+        if self.constraint_targets:
+            raise NotImplementedError(
+                f"'{type(self).__name__}.to_botorch_posterior_transform()' does not "
+                f"yet support objectives with outcome constraints. Use an MC-based "
+                f"acquisition function instead."
+            )
+
     @abstractmethod
     def to_botorch_posterior_transform(self) -> PosteriorTransform:
         """Convert to BoTorch posterior transform, if possible.
@@ -278,13 +340,54 @@ def identify_non_dominated_configurations(
         """
         from botorch.utils.multi_objective.pareto import is_non_dominated
 
-        validate_target_input(configurations, self.targets)
+        validate_target_input(configurations, self._optimization_targets)
 
         targets = self.transform(configurations, allow_extra=True)
         non_dominated = is_non_dominated(Y=to_tensor(targets), deduplicate=False)
 
         return pd.Series(non_dominated.numpy(), name="is_non_dominated")
 
+    def to_botorch_constraints(
+        self,
+    ) -> list[Callable[[torch.Tensor], torch.Tensor]]:
+        """Convert outcome constraints to BoTorch constraint callables.
+
+        Returns a list of callables, each mapping samples of shape
+        ``sample_shape x batch_shape x q x m`` to constraint values of shape
+        ``sample_shape x batch_shape x q``. A constraint is satisfied when
+        the return value is <= 0. Multiple constraints on the same target are
+        combined via ``torch.max`` (most-violated semantics).
+        """
+        if not self.outcome_constraints:
+            return []
+
+        from collections import defaultdict
+
+        import torch
+
+        # Group constraints by target index in model output
+        all_quantities = list(self._modeled_quantity_names)
+        constraints_by_target: dict[int, list[OutcomeConstraint]] = defaultdict(list)
+        for constraint in self.outcome_constraints:
+            if constraint.target.name in all_quantities:
+                idx = all_quantities.index(constraint.target.name)
+                constraints_by_target[idx].append(constraint)
+
+        # Build one callable per target; combine multiple via max (most-violated)
+        constraint_callables: list[Callable[[torch.Tensor], torch.Tensor]] = []
+        for target_idx, constraints in constraints_by_target.items():
+            funcs = [c.to_botorch_constraint_func(target_idx) for c in constraints]
+            if len(funcs) == 1:
+                constraint_callables.append(funcs[0])
+            else:
+                constraint_callables.append(
+                    lambda s, f=funcs: torch.max(
+                        torch.stack([fn(s) for fn in f]), dim=0
+                    )[0]
+                )
+
+        return constraint_callables
+
 
 def to_objective(x: Target | Objective, /) -> Objective:
     """Convert a target into an objective (with objective passthrough)."""
diff --git a/baybe/objectives/desirability.py b/baybe/objectives/desirability.py
index 5e371c6819..2b1eb46f57 100644
--- a/baybe/objectives/desirability.py
+++ b/baybe/objectives/desirability.py
@@ -111,7 +111,7 @@ class DesirabilityObjective(Objective):
     @weights.default
     def _default_weights(self) -> tuple[float, ...]:
         """Create unit weights for all targets."""
-        return tuple(1.0 for _ in range(len(self.targets)))
+        return tuple(1.0 for _ in range(len(self._optimization_targets)))
 
     @_targets.validator
     def _validate_targets(self, _, targets) -> None:  # noqa: DOC101, DOC103
@@ -147,7 +147,7 @@ def _validate_targets(self, _, targets) -> None:  # noqa: DOC101, DOC103
 
     @weights.validator
     def _validate_weights(self, _, weights) -> None:  # noqa: DOC101, DOC103
-        if (lw := len(weights)) != (lt := len(self.targets)):
+        if (lw := len(weights)) != (lt := len(self._optimization_targets)):
             raise ValueError(
                 f"If custom weights are specified, there must be one for each target. "
                 f"Specified number of targets: {lt}. Specified number of weights: {lw}."
@@ -155,14 +155,16 @@ def _validate_weights(self, _, weights) -> None:  # noqa: DOC101, DOC103
 
     @override
     @property
-    def targets(self) -> tuple[NumericalTarget, ...]:
+    def _optimization_targets(self) -> tuple[Target, ...]:
+        """Only the targets being optimized."""
         return self._targets
 
     @override
     @property
     def _modeled_quantities(self) -> tuple[Target, ...]:
         if self.as_pre_transformation:
-            return (NumericalTarget(_OUTPUT_NAME),)
+            # TODO: add constraint_targets
+            return (NumericalTarget(_OUTPUT_NAME),) + self.constraint_targets
         else:
             return self.targets
 
@@ -170,7 +172,9 @@ def _modeled_quantities(self) -> tuple[Target, ...]:
     @property
     def _model_quantities_to_target_names(self) -> dict[str, list[str]]:
         if self.as_pre_transformation:
-            return {_OUTPUT_NAME: [t.name for t in self.targets]}
+            result = {_OUTPUT_NAME: [t.name for t in self._optimization_targets]}
+            result.update({t.name: [t.name] for t in self.constraint_targets})
+            return result
         else:
             return {t.name: [t.name] for t in self.targets}
 
@@ -191,7 +195,7 @@ def normalized_weights(self) -> tuple[float, ...]:
 
     @override
     def __str__(self) -> str:
-        targets_list = [target.summary() for target in self.targets]
+        targets_list = [target.summary() for target in self._optimization_targets]
         targets_df = pd.DataFrame(targets_list)
         targets_df["Weights"] = self.normalized_weights
 
@@ -211,7 +215,7 @@ def _oriented_targets(self) -> tuple[Target, ...]:
         # to enable geometric averaging.
         return tuple(
             t.negate() + 1 if isinstance(t, NumericalTarget) and t.minimize else t
-            for t in self.targets
+            for t in self._optimization_targets
         )
 
     @override
@@ -265,10 +269,11 @@ def _to_botorch_full(self) -> MCAcquisitionObjective:
 
     @override
     def to_botorch_posterior_transform(self) -> ScalarizedPosteriorTransform:
+        self._check_posterior_transform_constraint_support()
         if self.as_pre_transformation:
             return IdentityTransformation().to_botorch_posterior_transform()
 
-        targets = self.targets
+        targets = self._optimization_targets
         transformations = [t.transformation for t in targets]
         if (
             self.scalarizer is not Scalarizer.MEAN
diff --git a/baybe/objectives/pareto.py b/baybe/objectives/pareto.py
index 987c2a981d..628bf145dc 100644
--- a/baybe/objectives/pareto.py
+++ b/baybe/objectives/pareto.py
@@ -11,6 +11,7 @@
 
 from baybe.objectives.base import Objective
 from baybe.objectives.validation import validate_target_names
+from baybe.targets.base import Target
 from baybe.targets.numerical import NumericalTarget
 from baybe.utils.basic import to_tuple
 
@@ -35,14 +36,10 @@ class ParetoObjective(Objective):
 
     @override
     @property
-    def targets(self) -> tuple[NumericalTarget, ...]:
+    def _optimization_targets(self) -> tuple[Target, ...]:
+        """Only the targets being optimized."""
         return self._targets
 
-    @override
-    @property
-    def output_names(self) -> tuple[str, ...]:
-        return tuple(target.name for target in self.targets)
-
     @override
     @property
     def supports_partial_measurements(self) -> bool:
diff --git a/baybe/objectives/single.py b/baybe/objectives/single.py
index d0f76b0f22..ef6438f3d8 100644
--- a/baybe/objectives/single.py
+++ b/baybe/objectives/single.py
@@ -45,18 +45,18 @@ def __str__(self) -> str:
             to_string("Targets", pretty_print_df(targets_df)),
         ]
 
+        if self.outcome_constraints:
+            constraints_str = "\n".join(str(c) for c in self.outcome_constraints)
+            fields.append(to_string("Outcome Constraints", constraints_str))
+
         return to_string("Objective", *fields)
 
     @override
     @property
-    def targets(self) -> tuple[Target, ...]:
+    def _optimization_targets(self) -> tuple[Target, ...]:
+        """Only the optimization targets."""
         return (self._target,)
 
-    @override
-    @property
-    def output_names(self) -> tuple[str, ...]:
-        return (self._target.name,)
-
     @override
     @property
     def supports_partial_measurements(self) -> bool:
@@ -75,6 +75,7 @@ def to_botorch(self) -> MCAcquisitionObjective:
 
     @override
     def to_botorch_posterior_transform(self) -> ScalarizedPosteriorTransform:
+        self._check_posterior_transform_constraint_support()
         if not (
             isinstance((t := self._target), NumericalTarget)
             and isinstance(
diff --git a/baybe/recommenders/naive.py b/baybe/recommenders/naive.py
index 8039755443..fece22ab8e 100644
--- a/baybe/recommenders/naive.py
+++ b/baybe/recommenders/naive.py
@@ -1,12 +1,14 @@
 """Naive recommender for hybrid spaces."""
 
 import gc
+import warnings
 from typing import ClassVar
 
 import pandas as pd
 from attrs import define, field
 from typing_extensions import override
 
+from baybe.exceptions import OutcomeConstraintIgnoredWarning
 from baybe.objectives.base import Objective
 from baybe.recommenders.pure.base import PureRecommender
 from baybe.recommenders.pure.bayesian.base import BayesianRecommender
@@ -85,6 +87,25 @@ def recommend(
 
         # We are in a hybrid setting now
 
+        # Outcome constraint compatibility check — Layer 1 gap fill.
+        # In this hybrid path, _recommend_discrete() is called directly on the
+        # disc_recommender, bypassing its recommend() and therefore the general
+        # Layer 1 warning in PureRecommender.recommend().
+        # See notes_compability_check_outcome_constraints.md for full design.
+        if (
+            not disc_is_bayesian
+            and objective is not None
+            and objective.outcome_constraints
+        ):
+            warnings.warn(
+                f"'{type(self.disc_recommender).__name__}' does not support outcome "
+                f"constraints. The constraints will be ignored for the discrete "
+                f"subspace during this recommendation step. Use a Bayesian "
+                f"recommender for the discrete subspace to enforce outcome "
+                f"constraints.",
+                OutcomeConstraintIgnoredWarning,
+            )
+
         # We will attach continuous parts to discrete parts and the other way round.
         # To make things simple, we sample a single point in the continuous space which
         # will then be attached to every discrete point when the acquisition function
diff --git a/baybe/recommenders/pure/base.py b/baybe/recommenders/pure/base.py
index 16eefe1016..b04522f093 100644
--- a/baybe/recommenders/pure/base.py
+++ b/baybe/recommenders/pure/base.py
@@ -1,6 +1,7 @@
 """Base classes for all pure recommenders."""
 
 import gc
+import warnings
 from abc import ABC
 from collections.abc import Callable
 from typing import Any, ClassVar, NoReturn
@@ -11,7 +12,11 @@
 from cattrs.gen import make_dict_unstructure_fn
 from typing_extensions import override
 
-from baybe.exceptions import DeprecationError, NotEnoughPointsLeftError
+from baybe.exceptions import (
+    DeprecationError,
+    NotEnoughPointsLeftError,
+    OutcomeConstraintIgnoredWarning,
+)
 from baybe.objectives.base import Objective
 from baybe.recommenders.base import RecommenderProtocol
 from baybe.searchspace import SearchSpace
@@ -102,6 +107,26 @@ def recommend(
         measurements: pd.DataFrame | None = None,
         pending_experiments: pd.DataFrame | None = None,
     ) -> pd.DataFrame:
+        # Outcome constraint compatibility check — Layer 1 (recommender level).
+        # Outcome constraints are only enforced by BayesianRecommender (via the acqf
+        # builder's _set_constraints, which is Layer 2). All other recommenders silently
+        # ignore them, so we warn the user here. BayesianRecommender is excluded via
+        # isinstance because it delegates enforcement to Layer 2.
+        from baybe.recommenders.pure.bayesian.base import BayesianRecommender
+
+        if (
+            not isinstance(self, BayesianRecommender)
+            and objective is not None
+            and objective.outcome_constraints
+        ):
+            warnings.warn(
+                f"'{self.__class__.__name__}' does not support outcome constraints. "
+                f"The constraints will be ignored during this recommendation step. "
+                f"Use a Bayesian recommender with a compatible acquisition function "
+                f"to enforce outcome constraints.",
+                OutcomeConstraintIgnoredWarning,
+            )
+
         if objective is not None:
             validate_object_names(searchspace.parameters + objective.targets)
 
diff --git a/baybe/targets/numerical.py b/baybe/targets/numerical.py
index 5aee308d53..25a5ee368d 100644
--- a/baybe/targets/numerical.py
+++ b/baybe/targets/numerical.py
@@ -122,8 +122,17 @@ class NumericalTarget(Target, SerialMixin):
     )
     """An optional target transformation."""
 
-    minimize: bool = field(default=False, validator=instance_of(bool), kw_only=True)
-    """Boolean flag indicating if the target is to be minimized."""
+    minimize: bool | None = field(
+        default=False,
+        validator=instance_of((bool, type(None))),
+        kw_only=True,
+    )
+    """Boolean flag indicating if the target is to be minimized.
+
+    - ``False``: maximize the target (default)
+    - ``True``: minimize the target
+    - ``None``: no optimization direction (constraint-only target)
+    """
 
     _constructor_info: dict[str, Any] | None = field(default=None, init=False, eq=False)
     """Helper to keep track of the target's construction details."""
diff --git a/tests/constraints/test_outcome_constraints.py b/tests/constraints/test_outcome_constraints.py
new file mode 100644
index 0000000000..fc6e5d914a
--- /dev/null
+++ b/tests/constraints/test_outcome_constraints.py
@@ -0,0 +1,290 @@
+"""Tests for outcome constraints functionality."""
+
+import numpy as np
+import pandas as pd
+import pytest
+import torch
+from pytest import param
+
+from baybe.acquisition.acqfs import qLogExpectedImprovement
+from baybe.constraints.outcome import OutcomeConstraint
+from baybe.objectives.desirability import DesirabilityObjective
+from baybe.objectives.pareto import ParetoObjective
+from baybe.objectives.single import SingleTargetObjective
+from baybe.parameters.numerical import NumericalDiscreteParameter
+from baybe.recommenders.pure.bayesian.botorch import BotorchRecommender
+from baybe.searchspace import SearchSpace
+from baybe.targets.numerical import NumericalTarget
+from baybe.transformations.basic import AffineTransformation
+
+# ---------------------------------------------------------------------------
+# Module-level test data
+# ---------------------------------------------------------------------------
+
+_temp_target = NumericalTarget("temperature", minimize=None)
+_yield_target = NumericalTarget("yield", minimize=False)
+_pressure_target = NumericalTarget("pressure", minimize=None)
+_pressure_constraint = OutcomeConstraint(_pressure_target, ">=", 2.0)
+_temp_constraint_le = OutcomeConstraint(_temp_target, "<=", 60.0)
+_temp_constraint_ge = OutcomeConstraint(_temp_target, ">=", 100.0)
+
+# Shared training data: yield and temperature both increase with x
+_rng = np.random.default_rng(42)
+_xs = np.linspace(1, 20, 15)
+_noise = _rng.normal(0, 1, 15)
+_measurements_temp = pd.DataFrame(
+    {"x": _xs, "yield": 5 * _xs + _noise, "temperature": 6 * _xs + 20 + _noise}
+)
+# For >=: yield decreases with x, temperature increases → constraint pushes high
+_measurements_temp_ge = pd.DataFrame(
+    {"x": _xs, "yield": -5 * _xs + 100 + _noise, "temperature": 6 * _xs + 20 + _noise}
+)
+
+# Targets for multi-target objectives
+_yield_target_normalized = NumericalTarget.normalized_ramp("yield", cutoffs=(0, 120))
+_purity_target_normalized = NumericalTarget.normalized_ramp("purity", cutoffs=(0, 1))
+_yield_target_pareto = NumericalTarget("yield", minimize=False)
+_purity_target_pareto = NumericalTarget("purity", minimize=False)
+
+# Measurements including yield, purity, and temperature for multi-target tests
+_measurements_multi = pd.DataFrame(
+    {
+        "x": _xs,
+        "yield": 5 * _xs + _noise,
+        "purity": 0.04 * _xs + 0.2 + _noise * 0.01,
+        "temperature": 6 * _xs + 20 + _noise,
+    }
+)
+
+# ---------------------------------------------------------------------------
+# Unit tests
+# ---------------------------------------------------------------------------
+
+
+@pytest.mark.parametrize(
+    ("operator", "sample_val", "threshold", "feasible"),
+    [
+        param("<=", 0.7, 150.0, True, id="le_feasible"),
+        param("<=", 0.8, 150.0, False, id="le_infeasible"),
+        param(">=", 0.8, 110.0, True, id="ge_feasible"),
+        param(">=", 0.4, 110.0, False, id="ge_infeasible"),
+    ],
+)
+def test_outcome_constraint(operator, sample_val, threshold, feasible):
+    """Threshold maps to computational space; feasible iff constraint_value <= 0."""
+    target = NumericalTarget(
+        "temperature",
+        transformation=AffineTransformation(factor=1 / 180, shift=-20 / 180),
+        minimize=None,
+    )
+    constraint = OutcomeConstraint(target, operator, threshold)
+
+    # Threshold transformation: experimental -> computational space
+    expected_comp = (threshold - 20) / 180
+    assert constraint.get_computational_threshold() == pytest.approx(
+        expected_comp, abs=1e-5
+    )
+
+    # Operator semantics: feasible iff constraint_value <= 0
+    constraint_func = constraint.to_botorch_constraint_func(0)
+    result = constraint_func(torch.tensor([[sample_val]]))
+    assert (result.item() <= 0) == feasible
+
+
+# ---------------------------------------------------------------------------
+# Objective with/without outcome constraint
+# ---------------------------------------------------------------------------
+
+
+@pytest.mark.parametrize(
+    ("outcome_constraints", "n_constraints", "n_targets", "has_str"),
+    [
+        param((), 0, 1, False, id="no_constraints"),
+        param((_temp_constraint_le,), 1, 2, True, id="single_constraint"),
+        param(
+            (_temp_constraint_le, _pressure_constraint),
+            2,
+            3,
+            True,
+            id="multiple_constraints",
+        ),
+    ],
+)
+def test_objective_with_outcome_constraints(
+    outcome_constraints, n_constraints, n_targets, has_str
+):
+    """SingleTargetObjective correctly integrates outcome constraints."""
+    objective = SingleTargetObjective(
+        target=_yield_target, outcome_constraints=outcome_constraints
+    )
+
+    assert len(objective.outcome_constraints) == n_constraints
+    assert len(objective.to_botorch_constraints()) == n_constraints
+    assert len(objective.targets) == n_targets
+    assert len(objective.constraint_targets) == n_targets - 1
+    assert objective._optimization_targets == (_yield_target,)
+    assert ("Outcome Constraints" in str(objective)) == has_str
+
+
+# ---------------------------------------------------------------------------
+# Behavioral tests
+# ---------------------------------------------------------------------------
+@pytest.fixture(name="wide_searchspace")
+def fixture_wide_searchspace():
+    """A discrete search space with x in [1, 20]."""
+    parameter = NumericalDiscreteParameter(name="x", values=list(range(1, 21)))
+    return SearchSpace.from_product(parameters=[parameter])
+
+
+@pytest.mark.parametrize(
+    (
+        "constrained_obj",
+        "unconstrained_obj",
+        "measurements",
+        "assert_direction",
+        "n_outputs",
+    ),
+    [
+        param(
+            SingleTargetObjective(
+                target=_yield_target,
+                outcome_constraints=(_temp_constraint_le,),
+            ),
+            SingleTargetObjective(target=_yield_target),
+            _measurements_temp,
+            "low",
+            2,
+            id="single_target_le",
+        ),
+        param(
+            SingleTargetObjective(
+                target=_yield_target,
+                outcome_constraints=(_temp_constraint_ge,),
+            ),
+            SingleTargetObjective(target=_yield_target),
+            _measurements_temp_ge,
+            "high",
+            2,
+            id="single_target_ge",
+        ),
+        param(
+            DesirabilityObjective(
+                targets=[_yield_target_normalized, _purity_target_normalized],
+                outcome_constraints=(_temp_constraint_le,),
+            ),
+            DesirabilityObjective(
+                targets=[_yield_target_normalized, _purity_target_normalized],
+            ),
+            _measurements_multi,
+            "low",
+            3,
+            id="desirability_le",
+        ),
+        param(
+            ParetoObjective(
+                targets=[_yield_target_pareto, _purity_target_pareto],
+                outcome_constraints=(_temp_constraint_le,),
+            ),
+            ParetoObjective(
+                targets=[_yield_target_pareto, _purity_target_pareto],
+            ),
+            _measurements_multi,
+            "low",
+            3,
+            id="pareto_le",
+        ),
+    ],
+)
+def test_constraint_effectiveness(
+    wide_searchspace,
+    constrained_obj,
+    unconstrained_obj,
+    measurements,
+    assert_direction,
+    n_outputs,
+):
+    """Outcome constraints steer recommendations toward the feasible region.
+
+    Compares constrained vs unconstrained recommendations to verify that
+    constraints shift the mean recommendation toward the feasible region.
+    """
+    torch.manual_seed(42)
+
+    means = []
+    recommenders = []
+    unconstrained_cols = [c for c in measurements.columns if c != "temperature"]
+    for objective, train_data in zip(
+        [constrained_obj, unconstrained_obj],
+        [measurements, measurements[unconstrained_cols]],
+    ):
+        recommender = BotorchRecommender()
+        result = recommender.recommend(
+            batch_size=5,
+            searchspace=wide_searchspace,
+            objective=objective,
+            measurements=train_data,
+        )
+        means.append(result["x"].mean())
+        recommenders.append(recommender)
+
+    mean_constrained, mean_unconstrained = means
+
+    # Constraint steers recommendations toward the feasible region
+    mean_small, mean_large = (
+        (mean_constrained, mean_unconstrained)
+        if assert_direction == "low"
+        else (mean_unconstrained, mean_constrained)
+    )
+    assert mean_small < mean_large, (
+        f"Constraint should steer toward {assert_direction} x: "
+        f"constrained mean={mean_constrained}, "
+        f"unconstrained mean={mean_unconstrained}"
+    )
+
+    # Constrained surrogate predicts all modeled targets
+    posterior = recommenders[0]._surrogate_model.posterior(
+        pd.DataFrame({"x": [3.0, 10.0, 17.0]}), joint=False
+    )
+    assert posterior.mean.shape[-1] == n_outputs, (
+        f"Surrogate should predict {n_outputs} outputs, got {posterior.mean.shape[-1]}"
+    )
+
+
+@pytest.mark.parametrize(
+    ("measurements", "best_f_range"),
+    [
+        pytest.param(
+            _measurements_temp,
+            (50.0, 80.0),
+            id="mixed_feasibility",
+        ),
+        pytest.param(
+            # All infeasible: falls back to global max from surrogate posterior mean
+            _measurements_temp.assign(
+                temperature=_measurements_temp["temperature"] + 200
+            ),
+            (80.0, 110.0),
+            id="all_infeasible",
+        ),
+    ],
+)
+def test_best_f_feasibility(wide_searchspace, measurements, best_f_range):
+    """best_f reflects best feasible value or falls back to global max."""
+    torch.manual_seed(42)
+
+    constraint = OutcomeConstraint(target=_temp_target, operator="<=", threshold=100.0)
+    objective = SingleTargetObjective(
+        target=_yield_target, outcome_constraints=(constraint,)
+    )
+
+    recommender = BotorchRecommender(acquisition_function=qLogExpectedImprovement())
+    recommender.get_acquisition_function(
+        searchspace=wide_searchspace,
+        objective=objective,
+        measurements=measurements,
+    )
+
+    best_f = recommender._botorch_acqf.best_f.item()
+    assert best_f_range[0] < best_f < best_f_range[1], (
+        f"best_f ({best_f}) should be in range {best_f_range}"
+    )
diff --git a/tests/validation/test_outcome_constraint_validation.py b/tests/validation/test_outcome_constraint_validation.py
new file mode 100644
index 0000000000..5290d97047
--- /dev/null
+++ b/tests/validation/test_outcome_constraint_validation.py
@@ -0,0 +1,138 @@
+"""Validation tests for outcome constraints."""
+
+import pandas as pd
+import pytest
+from pytest import param
+
+from baybe.acquisition.acqfs import qUpperConfidenceBound
+from baybe.constraints.outcome import OutcomeConstraint
+from baybe.exceptions import IncompatibilityError, OutcomeConstraintIgnoredWarning
+from baybe.objectives.single import SingleTargetObjective
+from baybe.parameters.numerical import NumericalDiscreteParameter
+from baybe.recommenders.pure.bayesian.botorch import BotorchRecommender
+from baybe.recommenders.pure.nonpredictive.sampling import RandomRecommender
+from baybe.searchspace import SearchSpace
+from baybe.targets.numerical import NumericalTarget
+
+# ---------------------------------------------------------------------------
+# Module-level test data
+# ---------------------------------------------------------------------------
+
+_temp_target = NumericalTarget("temperature", minimize=None)
+_yield_target = NumericalTarget("yield", minimize=False)
+
+_parameter = NumericalDiscreteParameter(name="x", values=list(range(1, 11)))
+_searchspace = SearchSpace.from_product(parameters=[_parameter])
+_measurements = pd.DataFrame(
+    {
+        "x": [1, 2, 3, 4, 5],
+        "yield": [10.0, 20.0, 30.0, 40.0, 50.0],
+        "temperature": [80.0, 90.0, 95.0, 105.0, 110.0],
+    }
+)
+
+# ---------------------------------------------------------------------------
+# OutcomeConstraint creation validation
+# ---------------------------------------------------------------------------
+
+
+@pytest.mark.parametrize(
+    ("kwargs", "error", "match"),
+    [
+        param(
+            {"target": "temperature", "operator": "<=", "threshold": 100.0},
+            TypeError,
+            None,
+            id="non_target_object",
+        ),
+        param(
+            {"target": _temp_target, "operator": ">", "threshold": 100.0},
+            ValueError,
+            "'operator' must be in",
+            id="invalid_operator",
+        ),
+        param(
+            {"target": _temp_target, "operator": "<=", "threshold": "high"},
+            ValueError,
+            "could not convert string to float",
+            id="non_numeric_threshold",
+        ),
+    ],
+)
+def test_invalid_outcome_constraint(kwargs, error, match):
+    """Invalid OutcomeConstraint arguments raise expected errors."""
+    with pytest.raises(error, match=match):
+        OutcomeConstraint(**kwargs)
+
+
+# ---------------------------------------------------------------------------
+# Objective + constraint combination validation
+# ---------------------------------------------------------------------------
+
+
+@pytest.mark.parametrize(
+    ("constraint_target", "match"),
+    [
+        param(
+            _yield_target,
+            "both optimized and constrained",
+            id="same_target",
+        ),
+        param(
+            NumericalTarget("yield", minimize=None),
+            "both optimized and constrained",
+            id="same_name_different_object",
+        ),
+        param(
+            NumericalTarget("temperature", minimize=False),
+            "must have minimize=None",
+            id="minimize_not_none",
+        ),
+    ],
+)
+def test_invalid_objective_with_outcome_constraint(constraint_target, match):
+    """Invalid outcome constraint + objective combinations raise expected errors."""
+    constraint = OutcomeConstraint(constraint_target, "<=", 100.0)
+    with pytest.raises(ValueError, match=match):
+        SingleTargetObjective(target=_yield_target, outcome_constraints=(constraint,))
+
+
+# ---------------------------------------------------------------------------
+# Recommender incompatibility
+# ---------------------------------------------------------------------------
+
+
+@pytest.mark.parametrize(
+    ("recommender", "expected"),
+    [
+        param(
+            BotorchRecommender(acquisition_function=qUpperConfidenceBound()),
+            IncompatibilityError,
+            id="incompatible_acqf",
+        ),
+        param(
+            RandomRecommender(),
+            OutcomeConstraintIgnoredWarning,
+            id="nonpredictive_recommender",
+        ),
+    ],
+)
+def test_recommender_constraint_incompatibility(recommender, expected):
+    """Incompatible recommenders raise errors or emit warnings."""
+    constraint = OutcomeConstraint(target=_temp_target, operator="<=", threshold=100.0)
+    objective = SingleTargetObjective(
+        target=_yield_target, outcome_constraints=(constraint,)
+    )
+
+    context = (
+        pytest.warns(expected)
+        if issubclass(expected, Warning)
+        else pytest.raises(expected, match="does not support outcome")
+    )
+    with context:
+        recommender.recommend(
+            batch_size=2,
+            searchspace=_searchspace,
+            objective=objective,
+            measurements=_measurements,
+        )