Repository: GridTools/gt4py
Pull request #2065 (2025-06-05 22:21:48)

feat[next]: concat where IR passes
Adds infrastructure for the concat_where implementation:
- Introduces itir.InfinityLiteral
- SymbolicDomain operations: union, intersection, promotion, ...
- concat_where specific IR transformations: canonicalize_domain_arguments, expand_tuple_args, transform_to_as_fieldop
- domain inference for concat_where


Patch:
diff --git a/src/gt4py/next/iterator/ir.py b/src/gt4py/next/iterator/ir.py
index d82ef930..e25eaeee 100644
--- a/src/gt4py/next/iterator/ir.py
+++ b/src/gt4py/next/iterator/ir.py
@@ -5,8 +5,10 @@
 #
 # Please, refer to the LICENSE file in the root directory.
 # SPDX-License-Identifier: BSD-3-Clause
+from __future__ import annotations
 
-from typing import ClassVar, List, Optional, Union
+import typing
+from typing import TYPE_CHECKING, ClassVar, List, Optional, Union
 
 import gt4py.eve as eve
 from gt4py.eve import Coerced, SymbolName, SymbolRef
@@ -63,6 +65,22 @@ class NoneLiteral(Expr):
     _none_literal: int = 0
 
 
+class InfinityLiteral(Expr):
+    # TODO(tehrengruber): self referential `ClassVar` not supported in eve.
+    if TYPE_CHECKING:
+        POSITIVE: ClassVar[InfinityLiteral]
+        NEGATIVE: ClassVar[InfinityLiteral]
+
+    name: typing.Literal["POSITIVE", "NEGATIVE"]
+
+    def __str__(self):
+        return f"{type(self).__name__}.{self.name}"
+
+
+InfinityLiteral.NEGATIVE = InfinityLiteral(name="NEGATIVE")
+InfinityLiteral.POSITIVE = InfinityLiteral(name="POSITIVE")
+
+
 class OffsetLiteral(Expr):
     value: Union[int, str]
 
@@ -142,4 +160,5 @@ FunctionDefinition.__hash__ = Node.__hash__  # type: ignore[method-assign]
 Program.__hash__ = Node.__hash__  # type: ignore[method-assign]
 SetAt.__hash__ = Node.__hash__  # type: ignore[method-assign]
 IfStmt.__hash__ = Node.__hash__  # type: ignore[method-assign]
+InfinityLiteral.__hash__ = Node.__hash__  # type: ignore[method-assign]
 Temporary.__hash__ = Node.__hash__  # type: ignore[method-assign]
diff --git a/src/gt4py/next/iterator/ir_utils/domain_utils.py b/src/gt4py/next/iterator/ir_utils/domain_utils.py
index afcf0292..52853899 100644
--- a/src/gt4py/next/iterator/ir_utils/domain_utils.py
+++ b/src/gt4py/next/iterator/ir_utils/domain_utils.py
@@ -10,11 +10,11 @@ from __future__ import annotations
 
 import dataclasses
 import functools
-from typing import Any, Literal, Mapping, Optional
+from typing import Any, Callable, Iterable, Literal, Mapping, Optional
 
 from gt4py.next import common
 from gt4py.next.iterator import builtins, ir as itir
-from gt4py.next.iterator.ir_utils import ir_makers as im
+from gt4py.next.iterator.ir_utils import common_pattern_matcher as cpm, ir_makers as im
 from gt4py.next.iterator.transforms import trace_shifts
 from gt4py.next.iterator.transforms.constant_folding import ConstantFolding
 
@@ -45,41 +45,43 @@ class SymbolicRange:
     start: itir.Expr
     stop: itir.Expr
 
+    def __post_init__(self) -> None:
+        # TODO(havogt): added this defensive checks as code seems to make this reasonable assumption
+        assert self.start is not itir.InfinityLiteral.POSITIVE
+        assert self.stop is not itir.InfinityLiteral.NEGATIVE
+
     def translate(self, distance: int) -> SymbolicRange:
         return SymbolicRange(im.plus(self.start, distance), im.plus(self.stop, distance))
 
 
+_GRID_TYPE_MAPPING = {
+    "unstructured_domain": common.GridType.UNSTRUCTURED,
+    "cartesian_domain": common.GridType.CARTESIAN,
+}
+
+
 @dataclasses.dataclass(frozen=True)
 class SymbolicDomain:
-    grid_type: Literal["unstructured_domain", "cartesian_domain"]
-    ranges: dict[
-        common.Dimension, SymbolicRange
-    ]  # TODO(havogt): remove `AxisLiteral` by `Dimension` everywhere
+    grid_type: common.GridType
+    ranges: dict[common.Dimension, SymbolicRange]
 
     def __hash__(self) -> int:
         return hash((self.grid_type, frozenset(self.ranges.items())))
 
     @classmethod
     def from_expr(cls, node: itir.Node) -> SymbolicDomain:
-        assert isinstance(node, itir.FunCall) and node.fun in [
-            im.ref("unstructured_domain"),
-            im.ref("cartesian_domain"),
-        ]
+        assert cpm.is_call_to(node, ("unstructured_domain", "cartesian_domain"))
 
         ranges: dict[common.Dimension, SymbolicRange] = {}
         for named_range in node.args:
-            assert (
-                isinstance(named_range, itir.FunCall)
-                and isinstance(named_range.fun, itir.SymRef)
-                and named_range.fun.id == "named_range"
-            )
+            assert cpm.is_call_to(named_range, "named_range")
             axis_literal, lower_bound, upper_bound = named_range.args
             assert isinstance(axis_literal, itir.AxisLiteral)
 
             ranges[common.Dimension(value=axis_literal.value, kind=axis_literal.kind)] = (
                 SymbolicRange(lower_bound, upper_bound)
             )
-        return cls(node.fun.id, ranges)  # type: ignore[attr-defined]  # ensure by assert above
+        return cls(_GRID_TYPE_MAPPING[node.fun.id], ranges)
 
     def as_expr(self) -> itir.FunCall:
         converted_ranges: dict[common.Dimension, tuple[itir.Expr, itir.Expr]] = {
@@ -164,22 +166,101 @@ class SymbolicDomain:
             raise AssertionError("Number of shifts must be a multiple of 2.")
 
 
-def domain_union(*domains: SymbolicDomain) -> SymbolicDomain:
-    """Return the (set) union of a list of domains."""
-    new_domain_ranges = {}
+def _reduce_ranges(
+    *ranges: SymbolicRange,
+    start_reduce_op: Callable[[itir.Expr, itir.Expr], itir.Expr],
+    stop_reduce_op: Callable[[itir.Expr, itir.Expr], itir.Expr],
+) -> SymbolicRange:
+    """Uses start_op and stop_op to fold the start and stop of a list of ranges."""
+    start = functools.reduce(
+        lambda current_expr, el_expr: start_reduce_op(current_expr, el_expr),
+        [range_.start for range_ in ranges],
+    )
+    stop = functools.reduce(
+        lambda current_expr, el_expr: stop_reduce_op(current_expr, el_expr),
+        [range_.stop for range_ in ranges],
+    )
+    # constant fold expression to keep the tree small
+    start, stop = ConstantFolding.apply(start), ConstantFolding.apply(stop)  # type: ignore[assignment]  # always an itir.Expr
+    return SymbolicRange(start, stop)
+
+
+_range_union = functools.partial(
+    _reduce_ranges, start_reduce_op=im.minimum, stop_reduce_op=im.maximum
+)
+_range_intersection = functools.partial(
+    _reduce_ranges, start_reduce_op=im.maximum, stop_reduce_op=im.minimum
+)
+
+
+def _reduce_domains(
+    *domains: SymbolicDomain,
+    range_reduce_op: Callable[..., SymbolicRange],
+) -> SymbolicDomain:
+    """
+    Applies range_op to the ranges of a list of domains with same dimensions and grid_type.
+    """
     assert all(domain.grid_type == domains[0].grid_type for domain in domains)
     assert all(domain.ranges.keys() == domains[0].ranges.keys() for domain in domains)
-    for dim in domains[0].ranges.keys():
-        start = functools.reduce(
-            lambda current_expr, el_expr: im.minimum(current_expr, el_expr),
-            [domain.ranges[dim].start for domain in domains],
-        )
-        stop = functools.reduce(
-            lambda current_expr, el_expr: im.maximum(current_expr, el_expr),
-            [domain.ranges[dim].stop for domain in domains],
-        )
-        # constant fold expression to keep the tree small
-        start, stop = ConstantFolding.apply(start), ConstantFolding.apply(stop)  # type: ignore[assignment]  # always an itir.Expr
-        new_domain_ranges[dim] = SymbolicRange(start, stop)
+
+    dims = domains[0].ranges.keys()
+    new_domain_ranges = {dim: range_reduce_op(*(d.ranges[dim] for d in domains)) for dim in dims}
 
     return SymbolicDomain(domains[0].grid_type, new_domain_ranges)
+
+
+domain_union = functools.partial(_reduce_domains, range_reduce_op=_range_union)
+"""Return the (set) union of a list of domains."""
+
+domain_intersection = functools.partial(_reduce_domains, range_reduce_op=_range_intersection)
+"""Return the intersection of a list of domains."""
+
+
+def domain_complement(domain: SymbolicDomain) -> SymbolicDomain:
+    """
+    Return the (set) complement of a half-infinite domain.
+
+    Note: after canonicalization of concat_where, the domain is always half-infinite,
+    i.e. it has ranges of the form `]-inf, a[` or `[a, inf[`.
+    """
+    dims_dict = {}
+    for dim in domain.ranges.keys():
+        lb, ub = domain.ranges[dim].start, domain.ranges[dim].stop
+        assert (lb == itir.InfinityLiteral.NEGATIVE) != (ub == itir.InfinityLiteral.POSITIVE)
+        # `]-inf, a[` -> `[a, inf[`
+        if lb == itir.InfinityLiteral.NEGATIVE:
+            dims_dict[dim] = SymbolicRange(start=ub, stop=itir.InfinityLiteral.POSITIVE)
+        # `[a, inf]` -> `]-inf, a]`
+        else:  # ub == itir.InfinityLiteral.POSITIVE:
+            dims_dict[dim] = SymbolicRange(start=itir.InfinityLiteral.NEGATIVE, stop=lb)
+    return SymbolicDomain(domain.grid_type, dims_dict)
+
+
+def promote_domain(
+    domain: SymbolicDomain, target_dims: Iterable[common.Dimension]
+) -> SymbolicDomain:
+    """Return a domain that is extended with the dimensions of target_dims."""
+    assert set(domain.ranges.keys()).issubset(target_dims)
+    dims_dict = {
+        dim: domain.ranges[dim]
+        if dim in domain.ranges
+        else SymbolicRange(itir.InfinityLiteral.NEGATIVE, itir.InfinityLiteral.POSITIVE)
+        for dim in target_dims
+    }
+    return SymbolicDomain(domain.grid_type, dims_dict)
+
+
+def is_finite(range_or_domain: SymbolicRange | SymbolicDomain) -> bool:
+    """
+    Return whether a range is unbounded in (at least) one direction.
+
+    The expression is required to be constant folded before for the result to be reliable.
+    """
+    match range_or_domain:
+        case SymbolicRange() as range_:
+            infinity_literals = (itir.InfinityLiteral.POSITIVE, itir.InfinityLiteral.NEGATIVE)
+            return not (range_.start in infinity_literals or range_.stop in infinity_literals)
+        case SymbolicDomain() as domain:
+            return all(is_finite(range_) for range_ in domain.ranges.values())
+        case _:
+            raise ValueError("Expected a 'SymbolicRange' or 'SymbolicDomain'.")
diff --git a/src/gt4py/next/iterator/ir_utils/ir_makers.py b/src/gt4py/next/iterator/ir_utils/ir_makers.py
index 7acdeb2f..739aa5d9 100644
--- a/src/gt4py/next/iterator/ir_utils/ir_makers.py
+++ b/src/gt4py/next/iterator/ir_utils/ir_makers.py
@@ -243,6 +243,12 @@ def if_(cond, true_val, false_val):
     return call("if_")(cond, true_val, false_val)
 
 
+def concat_where(cond, true_field, false_field):
+    """Create a concat_where FunCall, shorthand for ``call("concat_where")(expr)``."""
+
+    return call("concat_where")(cond, true_field, false_field)
+
+
 def lift(expr):
     """Create a lift FunCall, shorthand for ``call(call("lift")(expr))``."""
     return call(call("lift")(expr))
@@ -437,18 +443,18 @@ def domain(
     """
     if isinstance(grid_type, common.GridType):
         grid_type = f"{grid_type!s}_domain"
-    return call(grid_type)(
+    expr = call(grid_type)(
         *[
             call("named_range")(
-                itir.AxisLiteral(value=d.value, kind=d.kind)
-                if isinstance(d, common.Dimension)
-                else itir.AxisLiteral(value=d),
+                axis_literal(d),
                 r[0],
                 r[1],
             )
             for d, r in ranges.items()
         ]
     )
+    expr.type = ts.DomainType(dims=list(ranges.keys()))
+    return expr
 
 
 def as_fieldop(expr: itir.Expr | str, domain: Optional[itir.Expr] = None) -> Callable:
@@ -478,7 +484,8 @@ def as_fieldop(expr: itir.Expr | str, domain: Optional[itir.Expr] = None) -> Cal
     def _populate_domain_annex_wrapper(*args, **kwargs):
         node = result(*args, **kwargs)
         # note: if the domain is not a direct construction, e.g. because it is only a reference
-        # to a domain defined in a let, don't populate the annex
+        # to a domain defined in a let, don't populate the annex, since we can not create a
+        # symbolic domain for it.
         if domain and cpm.is_call_to(domain, ("cartesian_domain", "unstructured_domain")):
             node.annex.domain = domain_utils.SymbolicDomain.from_expr(domain)
         return node
@@ -515,6 +522,10 @@ def op_as_fieldop(
     return _impl
 
 
+def axis_literal(dim: common.Dimension) -> itir.AxisLiteral:
+    return itir.AxisLiteral(value=dim.value, kind=dim.kind)
+
+
 def cast_as_fieldop(type_: str, domain: Optional[itir.FunCall] = None):
     """
     Promotes the function `cast_` to a field_operator.
diff --git a/src/gt4py/next/iterator/ir_utils/misc.py b/src/gt4py/next/iterator/ir_utils/misc.py
index 0c75679c..00ff9abb 100644
--- a/src/gt4py/next/iterator/ir_utils/misc.py
+++ b/src/gt4py/next/iterator/ir_utils/misc.py
@@ -12,6 +12,7 @@ from typing import Callable, Iterable, TypeVar
 
 from gt4py import eve
 from gt4py.eve import utils as eve_utils
+from gt4py.next import common
 from gt4py.next.iterator import ir as itir
 from gt4py.next.iterator.ir_utils import common_pattern_matcher as cpm, ir_makers as im
 from gt4py.next.iterator.transforms import inline_lambdas
@@ -218,6 +219,24 @@ def extract_projector(
         return extract_projector(expr, projector, _depth + 1)
 
 
+def grid_type_from_domain(domain: itir.FunCall) -> common.GridType:
+    if cpm.is_call_to(domain, "cartesian_domain"):
+        return common.GridType.CARTESIAN
+    else:
+        assert cpm.is_call_to(domain, "unstructured_domain")
+        return common.GridType.UNSTRUCTURED
+
+
+def grid_type_from_program(program: itir.Program) -> common.GridType:
+    domains = program.walk_values().if_isinstance(itir.SetAt).getattr("domain").to_set()
+    grid_types = {grid_type_from_domain(d) for d in domains}
+    if len(grid_types) != 1:
+        raise ValueError(
+            f"Found 'set_at' with more than one 'GridType': '{grid_types}'. This is currently not supported."
+        )
+    return grid_types.pop()
+
+
 SymOrStr = TypeVar("SymOrStr", itir.Sym, str)
 
 
diff --git a/src/gt4py/next/iterator/pretty_printer.py b/src/gt4py/next/iterator/pretty_printer.py
index 7acbf5d2..5063e263 100644
--- a/src/gt4py/next/iterator/pretty_printer.py
+++ b/src/gt4py/next/iterator/pretty_printer.py
@@ -133,6 +133,13 @@ class PrettyPrinter(NodeTranslator):
     def visit_Literal(self, node: ir.Literal, *, prec: int) -> list[str]:
         return [str(node.value)]
 
+    def visit_InfinityLiteral(self, node: ir.InfinityLiteral, *, prec: int) -> list[str]:
+        if node == ir.InfinityLiteral.POSITIVE:
+            return ["∞"]
+        elif node == ir.InfinityLiteral.NEGATIVE:
+            return ["-∞"]
+        raise AssertionError()
+
     def visit_OffsetLiteral(self, node: ir.OffsetLiteral, *, prec: int) -> list[str]:
         return [str(node.value) + "ₒ"]
 
diff --git a/src/gt4py/next/iterator/transforms/collapse_tuple.py b/src/gt4py/next/iterator/transforms/collapse_tuple.py
index eae97b51..69881dd3 100644
--- a/src/gt4py/next/iterator/transforms/collapse_tuple.py
+++ b/src/gt4py/next/iterator/transforms/collapse_tuple.py
@@ -112,7 +112,7 @@ def _flattened_as_fieldop_param_el_name(param: str, idx: int) -> str:
 #  should revisit the pattern here and try to find a more general mechanism.
 @dataclasses.dataclass(frozen=True, kw_only=True)
 class CollapseTuple(
-    fixed_point_transformation.FixedPointTransformation, eve.PreserveLocationVisitor
+    fixed_point_transformation.CombinedFixedPointTransform, eve.PreserveLocationVisitor
 ):
     """
     Simplifies `make_tuple`, `tuple_get` calls.
diff --git a/src/gt4py/next/iterator/transforms/concat_where/__init__.py b/src/gt4py/next/iterator/transforms/concat_where/__init__.py
new file mode 100644
index 00000000..31f6872a
--- /dev/null
+++ b/src/gt4py/next/iterator/transforms/concat_where/__init__.py
@@ -0,0 +1,18 @@
+# GT4Py - GridTools Framework
+#
+# Copyright (c) 2014-2024, ETH Zurich
+# All rights reserved.
+#
+# Please, refer to the LICENSE file in the root directory.
+# SPDX-License-Identifier: BSD-3-Clause
+
+from gt4py.next.iterator.transforms.concat_where.canonicalize_domain_argument import (
+    canonicalize_domain_argument,
+)
+from gt4py.next.iterator.transforms.concat_where.expand_tuple_args import expand_tuple_args
+from gt4py.next.iterator.transforms.concat_where.transform_to_as_fieldop import (
+    transform_to_as_fieldop,
+)
+
+
+__all__ = ["canonicalize_domain_argument", "expand_tuple_args", "transform_to_as_fieldop"]
diff --git a/src/gt4py/next/iterator/transforms/concat_where/canonicalize_domain_argument.py b/src/gt4py/next/iterator/transforms/concat_where/canonicalize_domain_argument.py
new file mode 100644
index 00000000..37146604
--- /dev/null
+++ b/src/gt4py/next/iterator/transforms/concat_where/canonicalize_domain_argument.py
@@ -0,0 +1,97 @@
+# GT4Py - GridTools Framework
+#
+# Copyright (c) 2014-2024, ETH Zurich
+# All rights reserved.
+#
+# Please, refer to the LICENSE file in the root directory.
+# SPDX-License-Identifier: BSD-3-Clause
+
+from typing import Optional
+
+from gt4py.eve import PreserveLocationVisitor
+from gt4py.next.iterator import ir as itir
+from gt4py.next.iterator.ir_utils import (
+    common_pattern_matcher as cpm,
+    domain_utils,
+    ir_makers as im,
+)
+from gt4py.next.iterator.ir_utils.domain_utils import SymbolicDomain
+from gt4py.next.iterator.transforms import fixed_point_transformation
+
+
+def _range_complement(
+    range_: domain_utils.SymbolicRange,
+) -> tuple[domain_utils.SymbolicRange, domain_utils.SymbolicRange]:
+    # `[a, b[` -> `[-inf, a[` ∪ `[b, inf[`  # noqa: RUF003
+    assert not any(isinstance(b, itir.InfinityLiteral) for b in [range_.start, range_.stop])
+    return (
+        domain_utils.SymbolicRange(itir.InfinityLiteral.NEGATIVE, range_.start),
+        domain_utils.SymbolicRange(range_.stop, itir.InfinityLiteral.POSITIVE),
+    )
+
+
+class _CanonicalizeDomainArgument(
+    PreserveLocationVisitor, fixed_point_transformation.FixedPointTransformation
+):
+    """
+    TODO(tehrengruber): Explain why this is the canonical form.
+    """
+
+    @classmethod
+    def apply(cls, node: itir.Node):
+        return cls().visit(node)
+
+    def transform(self, node: itir.Node) -> Optional[itir.Node]:  # type: ignore[override] # ignore kwargs for simplicity
+        if cpm.is_call_to(node, "concat_where"):
+            cond_expr, field_a, field_b = node.args
+            if cpm.is_call_to(cond_expr, "and_"):
+                conds = cond_expr.args
+                return im.let(("__cwcda_field_a", field_a), ("__cwcda_field_b", field_b))(
+                    self.fp_transform(
+                        im.concat_where(
+                            conds[0],
+                            self.fp_transform(
+                                im.concat_where(conds[1], "__cwcda_field_a", "__cwcda_field_b")
+                            ),
+                            "__cwcda_field_b",
+                        )
+                    )
+                )
+            if cpm.is_call_to(cond_expr, "or_"):
+                conds = cond_expr.args
+                return im.let(("__cwcda_field_a", field_a), ("__cwcda_field_b", field_b))(
+                    self.fp_transform(
+                        im.concat_where(
+                            conds[0],
+                            "__cwcda_field_a",
+                            self.fp_transform(
+                                im.concat_where(conds[1], "__cwcda_field_a", "__cwcda_field_b")
+                            ),
+                        )
+                    )
+                )
+
+            # concat_where([1, 2[, a, b) -> concat_where([-inf, 1] | [2, inf[, b, a)
+            if cpm.is_call_to(cond_expr, ("cartesian_domain", "unstructured_domain")):
+                domain = SymbolicDomain.from_expr(cond_expr)
+                if len(domain.ranges) == 1:
+                    dim, range_ = next(iter(domain.ranges.items()))
+                    if domain_utils.is_finite(range_):
+                        complement = _range_complement(range_)
+                        new_domains = [
+                            im.domain(domain.grid_type, {dim: (cr.start, cr.stop)})
+                            for cr in complement
+                        ]
+                        return self.fp_transform(
+                            im.concat_where(im.call("or_")(*new_domains), field_b, field_a)
+                        )
+                else:
+                    # TODO(tehrengruber): Implement. Note that this case can not be triggered by
+                    #  the frontend yet since domains can only be created by expressions like
+                    #  `IDim < 10`.
+                    raise NotImplementedError()
+
+        return None
+
+
+canonicalize_domain_argument = _CanonicalizeDomainArgument.apply
diff --git a/src/gt4py/next/iterator/transforms/concat_where/expand_tuple_args.py b/src/gt4py/next/iterator/transforms/concat_where/expand_tuple_args.py
new file mode 100644
index 00000000..ea408697
--- /dev/null
+++ b/src/gt4py/next/iterator/transforms/concat_where/expand_tuple_args.py
@@ -0,0 +1,62 @@
+# GT4Py - GridTools Framework
+#
+# Copyright (c) 2014-2024, ETH Zurich
+# All rights reserved.
+#
+# Please, refer to the LICENSE file in the root directory.
+# SPDX-License-Identifier: BSD-3-Clause
+
+from gt4py.eve import NodeTranslator, PreserveLocationVisitor
+from gt4py.next import common
+from gt4py.next.iterator import ir as itir
+from gt4py.next.iterator.ir_utils import common_pattern_matcher as cpm, ir_makers as im
+from gt4py.next.iterator.type_system import inference as type_inference
+from gt4py.next.type_system import type_specifications as ts
+
+
+class _ExpandTupleArgs(PreserveLocationVisitor, NodeTranslator):
+    PRESERVED_ANNEX_ATTRS = (
+        "type",
+        "domain",
+    )
+
+    @classmethod
+    def apply(
+        cls,
+        node: itir.Node,
+        *,
+        offset_provider_type: common.OffsetProviderType,
+        allow_undeclared_symbols: bool = False,
+    ) -> itir.Node:
+        node = type_inference.infer(
+            node,
+            offset_provider_type=offset_provider_type,
+            allow_undeclared_symbols=allow_undeclared_symbols,
+        )
+        return cls().visit(node)
+
+    def visit_FunCall(self, node: itir.FunCall) -> itir.FunCall:
+        node = self.generic_visit(node)
+
+        # `concat_where(cond, {a, b}, {c, d})`
+        # -> `{concat_where(cond, a, c), concat_where(cond, a, c)}`
+        if cpm.is_call_to(node, "concat_where") and isinstance(
+            type_inference.reinfer(node.args[1]).type, ts.TupleType
+        ):
+            cond, true_branch, false_branch = node.args
+            new_els = []
+            assert isinstance(true_branch.type, ts.TupleType)
+            for i in range(len(true_branch.type.types)):
+                new_els.append(
+                    im.concat_where(cond, im.tuple_get(i, "__tb"), im.tuple_get(i, "__fb"))
+                )
+
+            new_node = im.let(("__tb", true_branch), ("__fb", false_branch))(
+                im.make_tuple(*new_els)
+            )
+            return new_node
+
+        return node
+
+
+expand_tuple_args = _ExpandTupleArgs.apply
diff --git a/src/gt4py/next/iterator/transforms/concat_where/transform_to_as_fieldop.py b/src/gt4py/next/iterator/transforms/concat_where/transform_to_as_fieldop.py
new file mode 100644
index 00000000..a693770a
--- /dev/null
+++ b/src/gt4py/next/iterator/transforms/concat_where/transform_to_as_fieldop.py
@@ -0,0 +1,93 @@
+# GT4Py - GridTools Framework
+#
+# Copyright (c) 2014-2024, ETH Zurich
+# All rights reserved.
+#
+# Please, refer to the LICENSE file in the root directory.
+# SPDX-License-Identifier: BSD-3-Clause
+
+import functools
+
+from gt4py.eve import NodeTranslator, PreserveLocationVisitor
+from gt4py.next import utils
+from gt4py.next.iterator import ir as itir
+from gt4py.next.iterator.ir_utils import (
+    common_pattern_matcher as cpm,
+    domain_utils,
+    ir_makers as im,
+)
+from gt4py.next.iterator.transforms import symbol_ref_utils
+from gt4py.next.iterator.type_system import inference as type_inference
+from gt4py.next.type_system import type_specifications as ts
+
+
+def _in(pos: itir.Expr, domain: itir.Expr) -> itir.Expr:
+    """
+    Given a position and a domain return an expression that evaluates to `True` if the position is inside the domain.
+
+    `in_({i, j, k}, u⟨ Iₕ: [i0, i1[, Iₕ: [j0, j1[, Iₕ: [k0, k1[ ⟩`
+    -> `i0 <= i < i1 & j0 <= j < j1 & k0 <= k < k1`
+    """
+    ret = [
+        im.and_(
+            im.less_equal(v.start, im.tuple_get(i, pos)),
+            im.less(im.tuple_get(i, pos), v.stop),
+        )
+        for i, v in enumerate(domain_utils.SymbolicDomain.from_expr(domain).ranges.values())
+    ]
+    return functools.reduce(im.and_, ret)
+
+
+class _TransformToAsFieldop(PreserveLocationVisitor, NodeTranslator):
+    PRESERVED_ANNEX_ATTRS = (
+        "type",
+        "domain",
+    )
+
+    @classmethod
+    def apply(cls, node: itir.Node):
+        """
+        Transform `concat_where` expressions into equivalent `as_fieldop` expressions.
+
+        Note that (backward) domain inference may not be executed after this pass as it can not
+        correctly infer the accessed domains when the value selection is represented as an `if_`
+        inside the `as_fieldop.
+        """
+        node = cls().visit(node)
+        node = type_inference.SanitizeTypes().visit(node)
+        return node
+
+    def visit_FunCall(self, node: itir.FunCall) -> itir.FunCall:
+        node = self.generic_visit(node)
+        if cpm.is_call_to(node, "concat_where"):
+            cond, true_branch, false_branch = node.args
+            assert isinstance(cond.type, ts.DomainType)
+            position = [im.index(dim) for dim in cond.type.dims]
+            refs = symbol_ref_utils.collect_symbol_refs(cond)
+
+            domains: tuple[domain_utils.SymbolicDomain, ...] = utils.flatten_nested_tuple(
+                node.annex.domain
+            )
+            assert all(domain == domains[0] for domain in domains), (
+                "At this point all `concat_where` arguments should be posed on the same domain."
+            )
+            assert isinstance(domains[0], domain_utils.SymbolicDomain)
+            domain_expr = domains[0].as_expr()
+
+            return im.as_fieldop(
+                im.lambda_("__tcw_pos", "__tcw_arg0", "__tcw_arg1", *refs)(
+                    im.let(*zip(refs, map(im.deref, refs), strict=True))(
+                        im.if_(
+                            _in(im.deref("__tcw_pos"), cond),
+                            im.deref("__tcw_arg0"),
+                            im.deref("__tcw_arg1"),
+                        )
+                    )
+                ),
+                domain_expr,
+            )(im.make_tuple(*position), true_branch, false_branch, *refs)
+
+        return node
+
+
+transform_to_as_fieldop = _TransformToAsFieldop.apply
diff --git a/src/gt4py/next/iterator/transforms/constant_folding.py b/src/gt4py/next/iterator/transforms/constant_folding.py
index 774e7a67..4a16122a 100644
--- a/src/gt4py/next/iterator/transforms/constant_folding.py
+++ b/src/gt4py/next/iterator/transforms/constant_folding.py
@@ -56,7 +56,7 @@ _COMMUTATIVE_OPS = ("plus", "multiplies", "minimum", "maximum")
 
 @dataclasses.dataclass(frozen=True, kw_only=True)
 class ConstantFolding(
-    fixed_point_transformation.FixedPointTransformation, eve.PreserveLocationVisitor
+    fixed_point_transformation.CombinedFixedPointTransform, eve.PreserveLocationVisitor
 ):
     PRESERVED_ANNEX_ATTRS = (
         "type",
diff --git a/src/gt4py/next/iterator/transforms/fixed_point_transformation.py b/src/gt4py/next/iterator/transforms/fixed_point_transformation.py
index b48e0cc3..3818f386 100644
--- a/src/gt4py/next/iterator/transforms/fixed_point_transformation.py
+++ b/src/gt4py/next/iterator/transforms/fixed_point_transformation.py
@@ -15,20 +15,11 @@ from gt4py.next.iterator import ir
 from gt4py.next.iterator.type_system import inference as itir_type_inference
 
 
-@dataclasses.dataclass(frozen=True, kw_only=True)
 class FixedPointTransformation(eve.NodeTranslator):
     """
-    Transformation pass that transforms until no transformation is applicable anymore.
+    Base class for iterative transformations that converge when a fixed-point is reached.
     """
 
-    #: Enum of all transformation (names). The transformations need to be defined as methods
-    #: named `transform_<NAME>`.
-    Transformation: ClassVar[Type[enum.Flag]]
-
-    #: All transformations enabled in this instance, e.g. `Transformation.T1 & Transformation.T2`.
-    #: Usually the default value is chosen to be all transformations.
-    enabled_transformations: enum.Flag
-
     REINFER_TYPES: ClassVar[bool] = False
 
     def visit(self, node, **kwargs):
@@ -43,18 +34,44 @@ class FixedPointTransformation(eve.NodeTranslator):
             new_node = self.transform(node, **kwargs)
             if new_node is None:
                 break
+            else:
+                new_node = self._post_transform(node, new_node)
             assert new_node != node
             node = new_node
         return node
 
+    def _post_transform(self, node: ir.Node, new_node: ir.Node) -> ir.Node:
+        if self.REINFER_TYPES:
+            itir_type_inference.reinfer(new_node)
+        self._preserve_annex(node, new_node)
+        return new_node
+
     def transform(self, node: ir.Node, **kwargs) -> Optional[ir.Node]:
         """
         Transform node once.
 
-        Execute transformations until one is applicable. As soon as a transformation occured
-        the function will return the transformed node. Note that the transformation itself
-        may call other transformations on child nodes again.
+        Execute transformation if applicable. When a transformation occurred the function will return
+        the transformed node. Note that the transformation itself may call other transformations on
+        child nodes again.
         """
+        ...
+
+
+@dataclasses.dataclass(frozen=True, kw_only=True)
+class CombinedFixedPointTransform(FixedPointTransformation):
+    """
+    Base class for a set of iterative transformations that converge when a fixed-point is reached.
+    """
+
+    #: Enum of all transformation (names). The transformations need to be defined as methods
+    #: named `transform_<NAME>`.
+    Transformation: ClassVar[Type[enum.Flag]]
+
+    #: All transformations enabled in this instance, e.g. `Transformation.T1 & Transformation.T2`.
+    #: Usually the default value is chosen to be all transformations.
+    enabled_transformations: enum.Flag
+
+    def transform(self, node: ir.Node, **kwargs) -> Optional[ir.Node]:
         for transformation in self.Transformation:
             if self.enabled_transformations & transformation:
                 assert isinstance(transformation.name, str)
@@ -64,8 +81,5 @@ class FixedPointTransformation(eve.NodeTranslator):
                     assert result is not node, (
                         f"Transformation {transformation.name.lower()} should have returned None, since nothing changed."
                     )
-                    if self.REINFER_TYPES:
-                        itir_type_inference.reinfer(result)
-                    self._preserve_annex(node, result)
                     return result
         return None
diff --git a/src/gt4py/next/iterator/transforms/fuse_as_fieldop.py b/src/gt4py/next/iterator/transforms/fuse_as_fieldop.py
index 26a8bcad..0e01dafe 100644
--- a/src/gt4py/next/iterator/transforms/fuse_as_fieldop.py
+++ b/src/gt4py/next/iterator/transforms/fuse_as_fieldop.py
@@ -231,7 +231,7 @@ def _make_tuple_element_inline_predicate(node: itir.Expr):
 
 @dataclasses.dataclass(frozen=True, kw_only=True)
 class FuseAsFieldOp(
-    fixed_point_transformation.FixedPointTransformation, eve.PreserveLocationVisitor
+    fixed_point_transformation.CombinedFixedPointTransform, eve.PreserveLocationVisitor
 ):
     """
     Merge multiple `as_fieldop` calls into one.
diff --git a/src/gt4py/next/iterator/transforms/infer_domain.py b/src/gt4py/next/iterator/transforms/infer_domain.py
index 3c6fbedf..c22b7754 100644
--- a/src/gt4py/next/iterator/transforms/infer_domain.py
+++ b/src/gt4py/next/iterator/transforms/infer_domain.py
@@ -21,6 +21,7 @@ from gt4py.next.iterator.ir_utils import (
     ir_makers as im,
     misc as ir_misc,
 )
+from gt4py.next.iterator.ir_utils.domain_utils import SymbolicDomain
 from gt4py.next.iterator.transforms import constant_folding, trace_shifts
 from gt4py.next.iterator.type_system import inference as itir_type_inference
 from gt4py.next.type_system import type_info, type_specifications as ts
@@ -56,6 +57,7 @@ class InferenceOptions(typing.TypedDict):
     offset_provider: common.OffsetProvider | common.OffsetProviderType
     symbolic_domain_sizes: Optional[dict[str, str]]
     allow_uninferred: bool
+    keep_existing_domains: bool
 
 
 class DomainAnnexDebugger(eve.NodeVisitor):
@@ -182,11 +184,16 @@ def _infer_as_fieldop(
     offset_provider: common.OffsetProvider | common.OffsetProviderType,
     symbolic_domain_sizes: Optional[dict[str, str]],
     allow_uninferred: bool,
+    keep_existing_domains: bool,
 ) -> tuple[itir.FunCall, AccessedDomains]:
     assert isinstance(applied_fieldop, itir.FunCall)
     assert cpm.is_call_to(applied_fieldop.fun, "as_fieldop")
     if not allow_uninferred and target_domain is DomainAccessDescriptor.NEVER:
         raise ValueError("'target_domain' cannot be 'NEVER' unless `allow_uninferred=True`.")
+
+    if len(applied_fieldop.fun.args) == 2 and keep_existing_domains:
+        target_domain = SymbolicDomain.from_expr(applied_fieldop.fun.args[1])
+
     # FIXME[#1582](tehrengruber): Temporary solution for `tuple_get` on scan result. See `test_solve_triag`.
     if isinstance(target_domain, tuple):
         target_domain = _domain_union(*flatten_nested_tuple(target_domain))  # type: ignore[arg-type]  # mypy not smart enough
@@ -226,6 +233,7 @@ def _infer_as_fieldop(
             offset_provider=offset_provider,
             symbolic_domain_sizes=symbolic_domain_sizes,
             allow_uninferred=allow_uninferred,
+            keep_existing_domains=keep_existing_domains,
         )
         transformed_inputs.append(transformed_input)
 
@@ -348,6 +356,40 @@ def _infer_if(
     return result_expr, actual_domains
 
 
+def _infer_concat_where(
+    expr: itir.Expr,
+    domain: DomainAccess,
+    **kwargs: Unpack[InferenceOptions],
+) -> tuple[itir.Expr, AccessedDomains]:
+    assert cpm.is_call_to(expr, "concat_where")
+    infered_args_expr = []
+    actual_domains: AccessedDomains = {}
+    cond, true_field, false_field = expr.args
+    symbolic_cond = domain_utils.SymbolicDomain.from_expr(cond)
+    cond_complement = domain_utils.domain_complement(symbolic_cond)
+
+    for arg in [true_field, false_field]:
+
+        @tree_map
+        def mapper(d: NonTupleDomainAccess):
+            if isinstance(d, DomainAccessDescriptor):
+                return d
+            promoted_cond = domain_utils.promote_domain(
+                symbolic_cond if arg == true_field else cond_complement,  # noqa: B023 # function is never used outside the loop
+                d.ranges.keys(),
+            )
+            return domain_utils.domain_intersection(d, promoted_cond)
+
+        domain_ = mapper(domain)
+
+        infered_arg_expr, actual_domains_arg = infer_expr(arg, domain_, **kwargs)
+        infered_args_expr.append(infered_arg_expr)
+        actual_domains = _merge_domains(actual_domains, actual_domains_arg)
+
+    result_expr = im.call(expr.fun)(cond, *infered_args_expr)
+    return result_expr, actual_domains
+
+
 def _infer_broadcast(
     expr: itir.Expr,
     domain: DomainAccess,
@@ -380,6 +422,8 @@ def _infer_expr(
         return _infer_tuple_get(expr, domain, **kwargs)
     elif cpm.is_call_to(expr, "if_"):
         return _infer_if(expr, domain, **kwargs)
+    elif cpm.is_call_to(expr, "concat_where"):
+        return _infer_concat_where(expr, domain, **kwargs)
     elif cpm.is_call_to(expr, "broadcast"):
         return _infer_broadcast(expr, domain, **kwargs)
     elif (
@@ -399,6 +443,7 @@ def infer_expr(
     offset_provider: common.OffsetProvider | common.OffsetProviderType,
     symbolic_domain_sizes: Optional[dict[str, str]] = None,
     allow_uninferred: bool = False,
+    keep_existing_domains: bool = False,
 ) -> tuple[itir.Expr, AccessedDomains]:
     """
     Infer the domain of all field subexpressions of `expr`.
@@ -409,10 +454,18 @@ def infer_expr(
     Arguments:
     - expr: The expression to be inferred.
     - domain: The domain `expr` is read at.
+
+    Keyword Arguments:
     - symbolic_domain_sizes: A dictionary mapping axes names, e.g., `I`, `Vertex`, to a symbol
       name that evaluates to the length of that axis.
     - allow_uninferred: Allow `as_fieldop` expressions whose domain is either unknown (e.g.
       because of a dynamic shift) or never accessed.
+    - keep_existing_domains: If `True`, keep existing domains in `as_fieldop` expressions and
+      use them to propagate the domain further. This is useful in cases where after a
+      transformation some nodes are missing domain information that needs to be repopulated,
+      but we can't reinfer everything because some domain access information has been lost.
+      For example when a `concat_where` is transformed into an `as_fieldop` with an if we lose
+      some information that could lead to unnecessary overcomputation and out-of-bounds accesses.
 
     Returns:
       A tuple containing the inferred expression with all applied `as_fieldop` (that are accessed)
@@ -458,8 +511,10 @@ def infer_expr(
         offset_provider=offset_provider,
         symbolic_domain_sizes=symbolic_domain_sizes,
         allow_uninferred=allow_uninferred,
+        keep_existing_domains=keep_existing_domains,
     )
-    expr.annex.domain = domain
+    if not keep_existing_domains or not hasattr(expr.annex, "domain"):
+        expr.annex.domain = domain
 
     return expr, accessed_domains
 
@@ -497,6 +552,7 @@ def infer_program(
     offset_provider: common.OffsetProvider | common.OffsetProviderType,
     symbolic_domain_sizes: Optional[dict[str, str]] = None,
     allow_uninferred: bool = False,
+    keep_existing_domains: bool = False,
 ) -> itir.Program:
     """
     Infer the domain of all field subexpressions inside a program.
@@ -522,6 +578,7 @@ def infer_program(
                 offset_provider=offset_provider,
                 symbolic_domain_sizes=symbolic_domain_sizes,
                 allow_uninferred=allow_uninferred,
+                keep_existing_domains=keep_existing_domains,
             )
             for stmt in program.body
         ],
diff --git a/src/gt4py/next/iterator/transforms/infer_domain_ops.py b/src/gt4py/next/iterator/transforms/infer_domain_ops.py
new file mode 100644
index 00000000..6447ceca
--- /dev/null
+++ b/src/gt4py/next/iterator/transforms/infer_domain_ops.py
@@ -0,0 +1,107 @@
+# GT4Py - GridTools Framework
+#
+# Copyright (c) 2014-2024, ETH Zurich
+# All rights reserved.
+#
+# Please, refer to the LICENSE file in the root directory.
+# SPDX-License-Identifier: BSD-3-Clause
+import dataclasses
+
+from gt4py.eve import NodeTranslator, PreserveLocationVisitor
+from gt4py.next import common
+from gt4py.next.iterator import builtins, ir as itir
+from gt4py.next.iterator.ir_utils import (
+    common_pattern_matcher as cpm,
+    domain_utils,
+    ir_makers as im,
+    misc as ir_misc,
+)
+from gt4py.next.iterator.type_system import inference
+from gt4py.next.type_system import type_specifications as ts
+
+
+@dataclasses.dataclass
+class InferDomainOps(PreserveLocationVisitor, NodeTranslator):
+    grid_type: common.GridType
+
+    @classmethod
+    def apply(cls, program: itir.Program):
+        return cls(grid_type=ir_misc.grid_type_from_program(program)).visit(program, recurse=True)
+
+    def visit_FunCall(self, node: itir.FunCall, **kwargs) -> itir.Node:
+        if kwargs["recurse"]:
+            node = self.generic_visit(node, **kwargs)
+
+        # e.g. `IDim < a`
+        if cpm.is_call_to(node, builtins.BINARY_MATH_COMPARISON_BUILTINS) and any(
+            isinstance(arg, itir.AxisLiteral) for arg in node.args
+        ):
+            arg1, arg2 = node.args
+            if isinstance(arg2, itir.AxisLiteral):
+                # take complementary operation if we have e.g. `0 < IDim` use `IDim > 0`
+                complementary_op = {
+                    "less": "greater",
+                    "less_equal": "greater_equal",
+                    "greater": "less",
+                    "greater_equal": "less_equal",
+                    "eq": "eq",
+                    "not_eq": "not_eq",
+                }
+                return self.visit(
+                    im.call(complementary_op[node.fun.id])(arg2, arg1),
+                    **{**kwargs, "recurse": False},
+                )
+
+            inference.reinfer(arg1)
+            assert isinstance(arg1.type, ts.DimensionType)
+            dim: common.Dimension = arg1.type.dim
+            value: itir.Expr = arg2
+
+            if cpm.is_call_to(node, ("less", "less_equal", "greater", "greater_equal", "eq")):
+                min_: itir.Expr
+                max_: itir.Expr
+
+                # `IDim < 1`
+                if cpm.is_call_to(node, "less"):
+                    min_ = itir.InfinityLiteral.NEGATIVE
+                    max_ = value
+                # `IDim <= 1`
+                elif cpm.is_call_to(node, "less_equal"):
+                    min_ = itir.InfinityLiteral.NEGATIVE
+                    max_ = im.plus(value, 1)
+                # `IDim > 1`
+                elif cpm.is_call_to(node, "greater"):
+                    min_ = im.plus(value, 1)
+                    max_ = itir.InfinityLiteral.POSITIVE
+                # `IDim >= 1`
+                elif cpm.is_call_to(node, "greater_equal"):
+                    min_ = value
+                    max_ = itir.InfinityLiteral.POSITIVE
+                # `IDim == 1`
+                elif cpm.is_call_to(node, "eq"):
+                    min_ = value
+                    max_ = im.plus(value, 1)
+
+                domain = domain_utils.SymbolicDomain(
+                    self.grid_type,
+                    ranges={dim: domain_utils.SymbolicRange(start=min_, stop=max_)},
+                )
+
+                return domain.as_expr()
+            elif cpm.is_call_to(node, "not_eq"):
+                # `IDim != a` -> `IDim < a & IDim > a`
+                return self.visit(
+                    im.call("and_")(
+                        self.visit(
+                            im.less(im.axis_literal(dim), value), **(kwargs | {"recurse": False})
+                        ),
+                        self.visit(
+                            im.greater(im.axis_literal(dim), value), **(kwargs | {"recurse": False})
+                        ),
+                    ),
+                    **(kwargs | {"recurse": False}),
+                )
+            else:
+                raise AssertionError()
+
+        return node
diff --git a/src/gt4py/next/iterator/transforms/inline_fundefs.py b/src/gt4py/next/iterator/transforms/inline_fundefs.py
index 03b20d14..2b8767e4 100644
--- a/src/gt4py/next/iterator/transforms/inline_fundefs.py
+++ b/src/gt4py/next/iterator/transforms/inline_fundefs.py
@@ -32,6 +32,7 @@ def prune_unreferenced_fundefs(program: itir.Program) -> itir.Program:
     """
     Remove all function declarations that are never called.
 
+    >>> from gt4py.next import common
     >>> from gt4py.next.iterator.ir_utils import ir_makers as im
     >>> fun1 = itir.FunctionDefinition(
     ...     id="fun1",
@@ -43,6 +44,7 @@ def prune_unreferenced_fundefs(program: itir.Program) -> itir.Program:
     ...     params=[im.sym("a")],
     ...     expr=im.deref("a"),
     ... )
+    >>> IDim = common.Dimension(value="IDim", kind=common.DimensionKind.HORIZONTAL)
     >>> program = itir.Program(
     ...     id="testee",
     ...     function_definitions=[fun1, fun2],
@@ -51,7 +53,7 @@ def prune_unreferenced_fundefs(program: itir.Program) -> itir.Program:
     ...     body=[
     ...         itir.SetAt(
     ...             expr=im.call("fun1")("inp"),
-    ...             domain=im.domain("cartesian_domain", {"IDim": (0, 10)}),
+    ...             domain=im.domain("cartesian_domain", {IDim: (0, 10)}),
     ...             target=im.ref("out"),
     ...         )
     ...     ],
diff --git a/src/gt4py/next/iterator/type_system/inference.py b/src/gt4py/next/iterator/type_system/inference.py
index 237ad351..c33e3a71 100644
--- a/src/gt4py/next/iterator/type_system/inference.py
+++ b/src/gt4py/next/iterator/type_system/inference.py
@@ -511,6 +511,12 @@ class ITIRTypeInference(eve.NodeTranslator):
         assert isinstance(node.type, ts.ScalarType)
         return node.type
 
+    def visit_InfinityLiteral(self, node: itir.InfinityLiteral, **kwargs) -> ts.ScalarType:
+        return ts.ScalarType(kind=ts.ScalarKind.INT32)
+
+    def visit_NegInfinityLiteral(self, node: itir.InfinityLiteral, **kwargs) -> ts.ScalarType:
+        return ts.ScalarType(kind=ts.ScalarKind.INT32)
+
     def visit_SymRef(
         self, node: itir.SymRef, *, ctx: dict[str, ts.TypeSpec]
     ) -> ts.TypeSpec | type_synthesizer.TypeSynthesizer:
diff --git a/src/gt4py/next/iterator/type_system/type_synthesizer.py b/src/gt4py/next/iterator/type_system/type_synthesizer.py
index 8f683c2f..ce995326 100644
--- a/src/gt4py/next/iterator/type_system/type_synthesizer.py
+++ b/src/gt4py/next/iterator/type_system/type_synthesizer.py
@@ -15,7 +15,7 @@ from typing import TypeVar, cast, overload
 
 from gt4py.eve import utils as eve_utils
 from gt4py.eve.extended_typing import Callable, Iterable, Optional, Union
-from gt4py.next import common
+from gt4py.next import common, utils
 from gt4py.next.ffront import fbuiltins
 from gt4py.next.iterator import builtins, ir as itir
 from gt4py.next.iterator.type_system import type_specifications as it_ts
@@ -149,13 +149,31 @@ def _(arg: ts.ScalarType) -> ts.ScalarType:
     return ts.ScalarType(kind=ts.ScalarKind.BOOL)
 
 
-@_register_builtin_type_synthesizer(
-    fun_names=builtins.BINARY_MATH_COMPARISON_BUILTINS | builtins.BINARY_LOGICAL_BUILTINS
-)
-def _(lhs: ts.ScalarType, rhs: ts.ScalarType) -> ts.ScalarType | ts.TupleType:
+def synthesize_binary_math_comparison_builtins(
+    lhs, rhs
+) -> ts.ScalarType | ts.TupleType | ts.DomainType:
+    if isinstance(lhs, ts.ScalarType) and isinstance(rhs, ts.DimensionType):
+        return ts.DomainType(dims=[rhs.dim])
+    if isinstance(lhs, ts.DimensionType) and isinstance(rhs, ts.ScalarType):
+        return ts.DomainType(dims=[lhs.dim])
+    assert all(isinstance(lhs, (ts.ScalarType, ts.DeferredType)) for arg in (lhs, rhs))
     return ts.ScalarType(kind=ts.ScalarKind.BOOL)
 
 
+@_register_builtin_type_synthesizer(fun_names=builtins.BINARY_MATH_COMPARISON_BUILTINS)
+def _(lhs, rhs) -> ts.ScalarType | ts.TupleType | ts.DomainType:
+    return synthesize_binary_math_comparison_builtins(lhs, rhs)
+
+
+@_register_builtin_type_synthesizer(fun_names=builtins.BINARY_LOGICAL_BUILTINS)
+def _(lhs, rhs) -> ts.ScalarType | ts.TupleType | ts.DomainType:
+    if isinstance(lhs, ts.DomainType) and isinstance(rhs, ts.DomainType):
+        assert lhs.dims != "unknown" and rhs.dims != "unknown"
+        return ts.DomainType(dims=common.promote_dims(lhs.dims, rhs.dims))
+    else:
+        return synthesize_binary_math_comparison_builtins(lhs, rhs)
+
+
 @_register_builtin_type_synthesizer
 def deref(it: it_ts.IteratorType | ts.DeferredType) -> ts.DataType | ts.DeferredType:
     if isinstance(it, ts.DeferredType):
@@ -245,6 +263,39 @@ def index(arg: ts.DimensionType) -> ts.FieldType:
     )
 
 
+@_register_builtin_type_synthesizer
+def concat_where(
+    domain: ts.DomainType,
+    true_field: ts.FieldType | ts.TupleType | ts.DeferredType,
+    false_field: ts.FieldType | ts.TupleType | ts.DeferredType,
+) -> ts.FieldType | ts.TupleType | ts.DeferredType:
+    if isinstance(true_field, ts.DeferredType) or isinstance(false_field, ts.DeferredType):
+        return ts.DeferredType(constraint=None)
+
+    @utils.tree_map(
+        collection_type=ts.TupleType,
+        result_collection_constructor=lambda el: ts.TupleType(types=list(el)),
+    )
+    def deduce_return_type(tb: ts.FieldType | ts.ScalarType, fb: ts.FieldType | ts.ScalarType):
+        if any(isinstance(b, ts.DeferredType) for b in [tb, fb]):
+            return ts.DeferredType(constraint=ts.FieldType)
+
+        tb_dtype, fb_dtype = (type_info.extract_dtype(b) for b in [tb, fb])
+
+        assert tb_dtype == fb_dtype, (
+            f"Field arguments must be of same dtype, got '{tb_dtype}' != '{fb_dtype}'."
+        )
+        dtype = tb_dtype
+
+        return_dims = common.promote_dims(
+            domain.dims, type_info.extract_dims(type_info.promote(tb, fb))
+        )
+        return_type = ts.FieldType(dims=return_dims, dtype=dtype)
+        return return_type
+
+    return deduce_return_type(true_field, false_field)
+
+
 @_register_builtin_type_synthesizer
 def broadcast(
     arg: ts.FieldType | ts.ScalarType | ts.DeferredType, dims: tuple[ts.DimensionType]
@@ -321,11 +372,7 @@ def _collect_and_check_dimensions(input_: ts.TypeSpec) -> list[common.Dimension]
         .filter(lambda dims: len(dims) > 0)
         .to_list()
     )
-    if all_input_dims:
-        assert all(cur_input_dims == all_input_dims[0] for cur_input_dims in all_input_dims)
-        return all_input_dims[0]
-
-    return []
+    return common.promote_dims(*all_input_dims)
 
 
 def _convert_as_fieldop_input_to_iterator(
diff --git a/src/gt4py/next/program_processors/codegens/gtfn/itir_to_gtfn_ir.py b/src/gt4py/next/program_processors/codegens/gtfn/itir_to_gtfn_ir.py
index 696cfc62..a4453905 100644
--- a/src/gt4py/next/program_processors/codegens/gtfn/itir_to_gtfn_ir.py
+++ b/src/gt4py/next/program_processors/codegens/gtfn/itir_to_gtfn_ir.py
@@ -84,24 +84,6 @@ def _get_domains(nodes: Iterable[itir.Stmt]) -> Iterable[itir.FunCall]:
     return result
 
 
-def _extract_grid_type(domain: itir.FunCall) -> common.GridType:
-    if domain.fun == itir.SymRef(id="cartesian_domain"):
-        return common.GridType.CARTESIAN
-    else:
-        assert domain.fun == itir.SymRef(id="unstructured_domain")
-        return common.GridType.UNSTRUCTURED
-
-
-def _get_gridtype(body: list[itir.Stmt]) -> common.GridType:
-    domains = _get_domains(body)
-    grid_types = {_extract_grid_type(d) for d in domains}
-    if len(grid_types) != 1:
-        raise ValueError(
-            f"Found 'set_at' with more than one 'GridType': '{grid_types}'. This is currently not supported."
-        )
-    return grid_types.pop()
-
-
 def _name_from_named_range(named_range_call: itir.FunCall) -> str:
     assert isinstance(named_range_call, itir.FunCall) and named_range_call.fun == itir.SymRef(
         id="named_range"
@@ -342,7 +324,7 @@ class GTFN_lowering(eve.NodeTranslator, eve.VisitorWithSymbolTableTrait):
             raise TypeError(f"Expected a 'Program', got '{type(node).__name__}'.")
 
         node = itir_type_inference.infer(node, offset_provider_type=offset_provider_type)
-        grid_type = _get_gridtype(node.body)
+        grid_type = ir_utils_misc.grid_type_from_program(node)
         if grid_type == common.GridType.UNSTRUCTURED:
             node = _CannonicalizeUnstructuredDomain.apply(node)
         return cls(
diff --git a/src/gt4py/next/program_processors/runners/dace/gtir_builtin_translators.py b/src/gt4py/next/program_processors/runners/dace/gtir_builtin_translators.py
index ac7df22c..3de4bd87 100644
--- a/src/gt4py/next/program_processors/runners/dace/gtir_builtin_translators.py
+++ b/src/gt4py/next/program_processors/runners/dace/gtir_builtin_translators.py
@@ -523,7 +523,6 @@ def extract_domain(node: gtir.Node) -> FieldopDomain:
             domain.append((dim, lower_bound, upper_bound))
 
     elif isinstance(node, domain_utils.SymbolicDomain):
-        assert str(node.grid_type) in {"cartesian_domain", "unstructured_domain"}
         for dim, drange in node.ranges.items():
             domain.append(
                 (dim, parse_range_boundary(drange.start), parse_range_boundary(drange.stop))


Test patch:
diff --git a/tests/next_tests/unit_tests/iterator_tests/ir_utils_test.py/test_domain_utils.py b/tests/next_tests/unit_tests/iterator_tests/ir_utils_test.py/test_domain_utils.py
new file mode 100644
index 00000000..69a2ed77
--- /dev/null
+++ b/tests/next_tests/unit_tests/iterator_tests/ir_utils_test.py/test_domain_utils.py
@@ -0,0 +1,182 @@
+# GT4Py - GridTools Framework
+#
+# Copyright (c) 2014-2024, ETH Zurich
+# All rights reserved.
+#
+# Please, refer to the LICENSE file in the root directory.
+# SPDX-License-Identifier: BSD-3-Clause
+
+import pytest
+
+from gt4py.next import common
+from gt4py.next.iterator import ir as itir
+from gt4py.next.iterator.ir_utils import domain_utils, ir_makers as im
+
+I = common.Dimension("I")
+J = common.Dimension("J")
+
+a_range = domain_utils.SymbolicRange(0, 10)
+another_range = domain_utils.SymbolicRange(5, 15)
+infinity_range = domain_utils.SymbolicRange(
+    itir.InfinityLiteral.NEGATIVE, itir.InfinityLiteral.POSITIVE
+)
+# the two next ranges are complement of each other (which is used in tests)
+right_infinity_range = domain_utils.SymbolicRange(0, itir.InfinityLiteral.POSITIVE)
+left_infinity_range = domain_utils.SymbolicRange(itir.InfinityLiteral.NEGATIVE, 0)
+
+
+def _make_domain(i: int):
+    return domain_utils.SymbolicDomain(
+        grid_type=common.GridType.CARTESIAN,
+        ranges={
+            I: domain_utils.SymbolicRange(im.ref(f"start_I_{i}"), im.ref(f"end_I_{i}")),
+            J: domain_utils.SymbolicRange(im.ref(f"start_J_{i}"), im.ref(f"end_J_{i}")),
+        },
+    )
+
+
+def test_symbolic_range():
+    with pytest.raises(AssertionError):
+        domain_utils.SymbolicRange(itir.InfinityLiteral.POSITIVE, 0)
+    with pytest.raises(AssertionError):
+        domain_utils.SymbolicRange(0, itir.InfinityLiteral.NEGATIVE)
+
+
+def test_domain_op_preconditions():
+    domain_a = domain_utils.SymbolicDomain(
+        grid_type=common.GridType.CARTESIAN,
+        ranges={I: domain_utils.SymbolicRange(0, 10)},
+    )
+    domain_b = domain_utils.SymbolicDomain(
+        grid_type=common.GridType.CARTESIAN,
+        ranges={J: domain_utils.SymbolicRange(5, 15)},
+    )
+    with pytest.raises(AssertionError):
+        domain_utils._reduce_domains(domain_a, domain_b, range_reduce_op=domain_utils._range_union)
+    domain_c = domain_utils.SymbolicDomain(
+        grid_type=common.GridType.UNSTRUCTURED, ranges={I: domain_utils.SymbolicRange(0, 10)}
+    )
+    with pytest.raises(AssertionError):
+        domain_utils._reduce_domains(domain_a, domain_c, range_reduce_op=domain_utils._range_union)
+
+
+def test_domain_union():
+    domain0 = _make_domain(0)
+    domain1 = _make_domain(1)
+    domain2 = _make_domain(2)
+
+    expected = domain_utils.SymbolicDomain(
+        grid_type=common.GridType.CARTESIAN,
+        ranges={
+            I: domain_utils.SymbolicRange(
+                im.minimum(
+                    im.minimum(im.ref("start_I_0"), im.ref("start_I_1")), im.ref("start_I_2")
+                ),
+                im.maximum(im.maximum(im.ref("end_I_0"), im.ref("end_I_1")), im.ref("end_I_2")),
+            ),
+            J: domain_utils.SymbolicRange(
+                im.minimum(
+                    im.minimum(im.ref("start_J_0"), im.ref("start_J_1")), im.ref("start_J_2")
+                ),
+                im.maximum(im.maximum(im.ref("end_J_0"), im.ref("end_J_1")), im.ref("end_J_2")),
+            ),
+        },
+    )
+    assert expected == domain_utils.domain_union(domain0, domain1, domain2)
+
+
+def test_domain_intersection():
+    domain0 = _make_domain(0)
+    domain1 = _make_domain(1)
+    domain2 = _make_domain(2)
+
+    expected = domain_utils.SymbolicDomain(
+        grid_type=common.GridType.CARTESIAN,
+        ranges={
+            I: domain_utils.SymbolicRange(
+                im.maximum(
+                    im.maximum(im.ref("start_I_0"), im.ref("start_I_1")), im.ref("start_I_2")
+                ),
+                im.minimum(im.minimum(im.ref("end_I_0"), im.ref("end_I_1")), im.ref("end_I_2")),
+            ),
+            J: domain_utils.SymbolicRange(
+                im.maximum(
+                    im.maximum(im.ref("start_J_0"), im.ref("start_J_1")), im.ref("start_J_2")
+                ),
+                im.minimum(im.minimum(im.ref("end_J_0"), im.ref("end_J_1")), im.ref("end_J_2")),
+            ),
+        },
+    )
+    assert expected == domain_utils.domain_intersection(domain0, domain1, domain2)
+
+
+@pytest.mark.parametrize(
+    "ranges, expected",
+    [
+        ({I: a_range}, None),
+        ({I: infinity_range}, None),
+        ({I: a_range, J: right_infinity_range}, None),
+        (
+            {I: right_infinity_range, J: left_infinity_range},
+            {I: left_infinity_range, J: right_infinity_range},
+        ),
+    ],
+)
+def test_domain_complement(ranges, expected):
+    if expected is None:
+        with pytest.raises(AssertionError):
+            domain_utils.domain_complement(
+                domain_utils.SymbolicDomain(grid_type=common.GridType.CARTESIAN, ranges=ranges)
+            )
+    else:
+        assert domain_utils.domain_complement(
+            domain_utils.SymbolicDomain(grid_type=common.GridType.CARTESIAN, ranges=ranges)
+        ) == domain_utils.SymbolicDomain(grid_type=common.GridType.CARTESIAN, ranges=expected)
+
+
+@pytest.mark.parametrize(
+    "testee_ranges, dimensions, expected_ranges",
+    [
+        ({I: a_range}, [I, J], {I: a_range, J: infinity_range}),
+        ({I: a_range}, [J, I], {I: a_range, J: infinity_range}),
+        ({I: a_range}, [I], {I: a_range}),
+        ({I: a_range}, [J], None),
+        ({I: a_range, J: another_range}, [J], None),
+    ],
+)
+def test_promote_domain(testee_ranges, dimensions, expected_ranges):
+    testee = domain_utils.SymbolicDomain(grid_type=common.GridType.CARTESIAN, ranges=testee_ranges)
+    if expected_ranges is None:
+        with pytest.raises(AssertionError):
+            domain_utils.promote_domain(testee, dimensions)
+    else:
+        expected = domain_utils.SymbolicDomain(
+            grid_type=common.GridType.CARTESIAN, ranges=expected_ranges
+        )
+        promoted = domain_utils.promote_domain(testee, dimensions)
+        assert promoted == expected
+
+
+def test_is_finite_symbolic_range():
+    assert not domain_utils.is_finite(infinity_range)
+    assert not domain_utils.is_finite(left_infinity_range)
+    assert not domain_utils.is_finite(right_infinity_range)
+    assert domain_utils.is_finite(a_range)
+
+
+@pytest.mark.parametrize(
+    "ranges, expected",
+    [
+        ({I: a_range, J: a_range}, True),
+        ({I: a_range, J: another_range}, True),
+        ({I: right_infinity_range, J: a_range}, False),
+        ({I: a_range, J: right_infinity_range}, False),
+    ],
+)
+def test_is_finite_symbolic_domain(ranges, expected):
+    assert (
+        domain_utils.is_finite(
+            domain_utils.SymbolicDomain(grid_type=common.GridType.CARTESIAN, ranges=ranges)
+        )
+        == expected
+    )
diff --git a/tests/next_tests/unit_tests/iterator_tests/test_type_inference.py b/tests/next_tests/unit_tests/iterator_tests/test_type_inference.py
index 4a574f25..a0361e7b 100644
--- a/tests/next_tests/unit_tests/iterator_tests/test_type_inference.py
+++ b/tests/next_tests/unit_tests/iterator_tests/test_type_inference.py
@@ -259,6 +259,31 @@ def expression_test_cases():
             ),
             ts.TupleType(types=[float_i_field, float_i_field]),
         ),
+        # concat_where
+        (
+            im.concat_where(
+                im.domain(common.GridType.CARTESIAN, {IDim: (0, 1)}),
+                im.ref("a", float_i_field),
+                im.ref("b", float_ij_field),
+            ),
+            float_ij_field,
+        ),
+        (
+            im.concat_where(
+                im.domain(common.GridType.CARTESIAN, {IDim: (0, 1)}),
+                im.ref("a", ts.TupleType(types=[float_i_field] * 2)),
+                im.ref("b", ts.TupleType(types=[float_i_field] * 2)),
+            ),
+            ts.TupleType(types=[float_i_field] * 2),
+        ),
+        (
+            im.concat_where(
+                im.domain(common.GridType.CARTESIAN, {IDim: (0, 1)}),
+                im.ref("a", ts.TupleType(types=[float_i_field, float_ij_field])),
+                im.ref("b", ts.TupleType(types=[float_i_field] * 2)),
+            ),
+            ts.TupleType(types=[float_i_field, float_ij_field]),
+        ),
     )
 
 
diff --git a/tests/next_tests/unit_tests/iterator_tests/transforms_tests/test_concat_where_canonicalize_domain_args.py b/tests/next_tests/unit_tests/iterator_tests/transforms_tests/test_concat_where_canonicalize_domain_args.py
new file mode 100644
index 00000000..58086670
--- /dev/null
+++ b/tests/next_tests/unit_tests/iterator_tests/transforms_tests/test_concat_where_canonicalize_domain_args.py
@@ -0,0 +1,68 @@
+# GT4Py - GridTools Framework
+#
+# Copyright (c) 2014-2024, ETH Zurich
+# All rights reserved.
+#
+# Please, refer to the LICENSE file in the root directory.
+# SPDX-License-Identifier: BSD-3-Clause
+from gt4py.next import common
+import pytest
+from gt4py.next.iterator import ir as itir
+from gt4py.next.iterator.ir_utils import ir_makers as im
+from gt4py.next.iterator.transforms import (
+    concat_where,
+    inline_lambdas,
+)
+from gt4py.next.type_system import type_specifications as ts
+
+int_type = ts.ScalarType(kind=ts.ScalarKind.INT32)
+IDim = common.Dimension(value="IDim", kind=common.DimensionKind.HORIZONTAL)
+field_type = ts.FieldType(dims=[IDim], dtype=int_type)
+
+
+def cases():
+    return [
+        # testee, expected
+        (
+            im.concat_where(im.and_("cond1", "cond2"), "a", "b"),
+            im.concat_where("cond1", im.concat_where("cond2", "a", "b"), "b"),
+        ),
+        (
+            im.concat_where(im.or_("cond1", "cond2"), "a", "b"),
+            im.concat_where("cond1", "a", im.concat_where("cond2", "a", "b")),
+        ),
+        (
+            im.concat_where(im.domain(common.GridType.CARTESIAN, {IDim: (0, 1)}), "a", "b"),
+            im.concat_where(
+                im.domain(common.GridType.CARTESIAN, {IDim: (itir.InfinityLiteral.NEGATIVE, 0)}),
+                "b",
+                im.concat_where(
+                    im.domain(
+                        common.GridType.CARTESIAN, {IDim: (1, itir.InfinityLiteral.POSITIVE)}
+                    ),
+                    "b",
+                    "a",
+                ),
+            ),
+        ),
+        (  # not transformed
+            im.concat_where(
+                im.domain(common.GridType.CARTESIAN, {IDim: (0, itir.InfinityLiteral.POSITIVE)}),
+                "a",
+                "b",
+            ),
+            im.concat_where(
+                im.domain(common.GridType.CARTESIAN, {IDim: (0, itir.InfinityLiteral.POSITIVE)}),
+                "a",
+                "b",
+            ),
+        ),
+    ]
+
+
+@pytest.mark.parametrize("testee, expected", cases())
+def test_nested_concat_where(testee, expected):
+    actual = concat_where.canonicalize_domain_argument(testee)
+    actual = inline_lambdas.InlineLambdas.apply(actual, opcount_preserving=True)
+
+    assert actual == expected
diff --git a/tests/next_tests/unit_tests/iterator_tests/transforms_tests/test_concat_where_expand_tuple_args.py b/tests/next_tests/unit_tests/iterator_tests/transforms_tests/test_concat_where_expand_tuple_args.py
new file mode 100644
index 00000000..42ad2920
--- /dev/null
+++ b/tests/next_tests/unit_tests/iterator_tests/transforms_tests/test_concat_where_expand_tuple_args.py
@@ -0,0 +1,53 @@
+# GT4Py - GridTools Framework
+#
+# Copyright (c) 2014-2024, ETH Zurich
+# All rights reserved.
+#
+# Please, refer to the LICENSE file in the root directory.
+# SPDX-License-Identifier: BSD-3-Clause
+from gt4py.next import common
+import pytest
+from gt4py.next.iterator import ir as itir
+from gt4py.next.iterator.ir_utils import ir_makers as im, domain_utils
+from gt4py.next.iterator.transforms import (
+    concat_where,
+    inline_lambdas,
+    infer_domain,
+    collapse_tuple,
+)
+from gt4py.next.type_system import type_specifications as ts
+from gt4py.next.iterator.type_system import type_specifications as it_ts
+
+int_type = ts.ScalarType(kind=ts.ScalarKind.INT32)
+IDim = common.Dimension(value="IDim", kind=common.DimensionKind.HORIZONTAL)
+field_type = ts.FieldType(dims=[IDim], dtype=int_type)
+
+
+def test_trivial():
+    cond = im.domain(common.GridType.CARTESIAN, {IDim: (itir.InfinityLiteral.NEGATIVE, 1)})
+    domain = im.domain(common.GridType.CARTESIAN, {IDim: (0, 2)})
+    symbolic_domain = domain_utils.SymbolicDomain.from_expr(domain)
+
+    testee = im.concat_where(
+        cond,
+        im.make_tuple(im.ref("a", field_type), im.ref("c", field_type)),
+        im.make_tuple(im.ref("b", field_type), im.ref("d", field_type)),
+    )
+    testee, _ = infer_domain.infer_expr(
+        testee,
+        (symbolic_domain, symbolic_domain),
+        keep_existing_domains=True,
+        offset_provider={},
+    )
+
+    expected = im.make_tuple(im.concat_where(cond, "a", "b"), im.concat_where(cond, "c", "d"))
+
+    actual = concat_where.expand_tuple_args(
+        testee, offset_provider_type={}, allow_undeclared_symbols=True
+    )
+
+    actual = collapse_tuple.CollapseTuple.apply(
+        actual, allow_undeclared_symbols=True, within_stencil=False
+    )
+
+    assert actual == expected
diff --git a/tests/next_tests/unit_tests/iterator_tests/transforms_tests/test_concat_where_transform_to_as_fieldop.py b/tests/next_tests/unit_tests/iterator_tests/transforms_tests/test_concat_where_transform_to_as_fieldop.py
new file mode 100644
index 00000000..2517c7ab
--- /dev/null
+++ b/tests/next_tests/unit_tests/iterator_tests/transforms_tests/test_concat_where_transform_to_as_fieldop.py
@@ -0,0 +1,89 @@
+# GT4Py - GridTools Framework
+#
+# Copyright (c) 2014-2024, ETH Zurich
+# All rights reserved.
+#
+# Please, refer to the LICENSE file in the root directory.
+# SPDX-License-Identifier: BSD-3-Clause
+from gt4py.next import common
+
+from gt4py.next.iterator.ir_utils import ir_makers as im, domain_utils
+from gt4py.next.iterator.transforms import concat_where, inline_lambdas
+from gt4py.next.iterator.transforms.concat_where import transform_to_as_fieldop
+from gt4py.next.iterator.transforms.concat_where.transform_to_as_fieldop import _in
+from gt4py.next.type_system import type_specifications as ts
+
+int_type = ts.ScalarType(kind=ts.ScalarKind.INT32)
+IDim = common.Dimension(value="IDim", kind=common.DimensionKind.HORIZONTAL)
+JDim = common.Dimension(value="JDim", kind=common.DimensionKind.HORIZONTAL)
+
+
+def test_in_helper():
+    pos = im.make_tuple(0, 1)
+    bounds = {
+        IDim: (3, 4),
+        JDim: (5, 6),
+    }
+    expected = im.and_(
+        im.and_(
+            im.less_equal(bounds[IDim][0], im.tuple_get(0, pos)),
+            im.less(im.tuple_get(0, pos), bounds[IDim][1]),
+        ),
+        im.and_(
+            im.less_equal(bounds[JDim][0], im.tuple_get(1, pos)),
+            im.less(im.tuple_get(1, pos), bounds[JDim][1]),
+        ),
+    )
+    actual = _in(pos, im.domain(common.GridType.CARTESIAN, bounds))
+    assert actual == expected
+
+
+def test_trivial():
+    domain = im.domain(common.GridType.CARTESIAN, {IDim: (0, 2)})
+
+    cond = im.domain(common.GridType.CARTESIAN, {IDim: (0, 1)})
+    testee = im.concat_where(cond, "true_branch", "false_branch")
+    testee.annex.domain = domain_utils.SymbolicDomain.from_expr(domain)
+    expected = im.as_fieldop(
+        im.lambda_("__tcw_pos", "__tcw_arg0", "__tcw_arg1")(
+            im.if_(
+                _in(im.deref("__tcw_pos"), cond),
+                im.deref("__tcw_arg0"),
+                im.deref("__tcw_arg1"),
+            )
+        ),
+        domain,
+    )(im.make_tuple(im.index(IDim)), "true_branch", "false_branch")
+
+    actual = concat_where.transform_to_as_fieldop(testee)
+    actual = inline_lambdas.InlineLambdas.apply(actual)  # simplify
+
+    assert actual == expected
+
+
+def test_capturing_cond():
+    domain = im.domain(common.GridType.CARTESIAN, {IDim: (0, 1)})
+
+    cond = im.domain(common.GridType.CARTESIAN, {IDim: ("start", "stop")})
+    testee = im.concat_where(cond, "true_branch", "false_branch")
+    testee.annex.domain = domain_utils.SymbolicDomain.from_expr(domain)
+    expected = im.as_fieldop(
+        im.lambda_("__tcw_pos", "__tcw_arg0", "__tcw_arg1", "start", "stop")(
+            im.if_(
+                _in(
+                    im.deref("__tcw_pos"),
+                    im.domain(
+                        common.GridType.CARTESIAN, {IDim: (im.deref("start"), im.deref("stop"))}
+                    ),
+                ),
+                im.deref("__tcw_arg0"),
+                im.deref("__tcw_arg1"),
+            )
+        ),
+        domain,
+    )(im.make_tuple(im.index(IDim)), "true_branch", "false_branch", "start", "stop")
+
+    actual = concat_where.transform_to_as_fieldop(testee)
+    actual = inline_lambdas.InlineLambdas.apply(actual)  # simplify
+
+    assert actual == expected
diff --git a/tests/next_tests/unit_tests/iterator_tests/transforms_tests/test_infer_domain_ops.py b/tests/next_tests/unit_tests/iterator_tests/transforms_tests/test_infer_domain_ops.py
new file mode 100644
index 00000000..77ba3719
--- /dev/null
+++ b/tests/next_tests/unit_tests/iterator_tests/transforms_tests/test_infer_domain_ops.py
@@ -0,0 +1,71 @@
+# GT4Py - GridTools Framework
+#
+# Copyright (c) 2014-2024, ETH Zurich
+# All rights reserved.
+#
+# Please, refer to the LICENSE file in the root directory.
+# SPDX-License-Identifier: BSD-3-Clause
+
+import pytest
+import textwrap
+
+from gt4py.eve.utils import UIDGenerator
+from gt4py.next import common
+from gt4py.next.iterator import ir as itir
+from gt4py.next.iterator.ir_utils import ir_makers as im
+from gt4py.next.iterator.transforms.infer_domain_ops import InferDomainOps
+from gt4py.next.iterator.transforms.constant_folding import ConstantFolding
+
+from next_tests.integration_tests.cases import IDim, JDim, KDim
+
+
+def test_data():
+    return [
+        (
+            im.less(im.axis_literal(IDim), 1),
+            im.domain(common.GridType.CARTESIAN, {IDim: (itir.InfinityLiteral.NEGATIVE, 1)}),
+        ),
+        (
+            im.less_equal(im.axis_literal(IDim), 1),
+            im.domain(common.GridType.CARTESIAN, {IDim: (itir.InfinityLiteral.NEGATIVE, 2)}),
+        ),
+        (
+            im.greater(im.axis_literal(IDim), 1),
+            im.domain(common.GridType.CARTESIAN, {IDim: (2, itir.InfinityLiteral.POSITIVE)}),
+        ),
+        (
+            im.greater_equal(im.axis_literal(IDim), 1),
+            im.domain(common.GridType.CARTESIAN, {IDim: (1, itir.InfinityLiteral.POSITIVE)}),
+        ),
+        (
+            im.less(1, im.axis_literal(IDim)),
+            im.domain(common.GridType.CARTESIAN, {IDim: (2, itir.InfinityLiteral.POSITIVE)}),
+        ),
+        (
+            im.less_equal(1, im.axis_literal(IDim)),
+            im.domain(common.GridType.CARTESIAN, {IDim: (1, itir.InfinityLiteral.POSITIVE)}),
+        ),
+        (
+            im.greater(1, im.axis_literal(IDim)),
+            im.domain(common.GridType.CARTESIAN, {IDim: (itir.InfinityLiteral.NEGATIVE, 1)}),
+        ),
+        (
+            im.greater_equal(1, im.axis_literal(IDim)),
+            im.domain(common.GridType.CARTESIAN, {IDim: (itir.InfinityLiteral.NEGATIVE, 2)}),
+        ),
+        (im.eq(1, im.axis_literal(IDim)), im.domain(common.GridType.CARTESIAN, {IDim: (1, 2)})),
+        (
+            im.not_eq(1, im.axis_literal(IDim)),
+            im.and_(
+                im.domain(common.GridType.CARTESIAN, {IDim: (itir.InfinityLiteral.NEGATIVE, 1)}),
+                im.domain(common.GridType.CARTESIAN, {IDim: (2, itir.InfinityLiteral.POSITIVE)}),
+            ),
+        ),
+    ]
+
+
+@pytest.mark.parametrize("testee,expected", test_data())
+def test_trivial(testee, expected):
+    actual = InferDomainOps(grid_type=common.GridType.CARTESIAN).visit(testee, recurse=True)
+    actual = ConstantFolding.apply(actual)  # simplify expr to get simpler expected expressions
+    assert actual == expected