diff --git a/dwave/optimization/include/dwave-optimization/interval.hpp b/dwave/optimization/include/dwave-optimization/interval.hpp
new file mode 100644
index 00000000..36c9365c
--- /dev/null
+++ b/dwave/optimization/include/dwave-optimization/interval.hpp
@@ -0,0 +1,218 @@
+// Copyright 2026 D-Wave
+//
+//    Licensed under the Apache License, Version 2.0 (the "License");
+//    you may not use this file except in compliance with the License.
+//    You may obtain a copy of the License at
+//
+//        http://www.apache.org/licenses/LICENSE-2.0
+//
+//    Unless required by applicable law or agreed to in writing, software
+//    distributed under the License is distributed on an "AS IS" BASIS,
+//    WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+//    See the License for the specific language governing permissions and
+//    limitations under the License.
+
+#pragma once
+
+#include <algorithm>
+#include <cassert>
+#include <cmath>
+#include <iosfwd>
+#include <type_traits>
+
+#include "dwave-optimization/common.hpp"  // for ssize_t
+#include "dwave-optimization/typing.hpp"
+
+namespace dwave::optimization {
+
+// dev note: true interval arithmetic requires "outward" rounding. For right
+// now I am ignoring this and focusing on getting the API solidified, but
+// we should support it ASAP. I have marked the places we need to make changes
+// with a "numeric issues" comment for easy future reference.
+
+template <DType T>
+requires(not std::is_const_v<T>) class interval {
+ public:
+    // by default we construct an empty interval
+    interval() = default;
+
+    // rule of five stuff all works the way you'd expect
+    interval(const interval&) = default;
+    interval(interval&&) = default;
+    ~interval() = default;
+    interval& operator=(const interval&) = default;
+    interval& operator=(interval&&) = default;
+
+    // todo: note the difference between interval<int>(.5, .5) and interval<int>(interval<float>(.5,
+    // .5))
+    explicit interval(T value) noexcept : interval(value, value) {}
+    interval(T infimum, T supremum) noexcept : infimum_(infimum), supremum_(supremum) {
+        assert(not std::isnan(infimum_) and "infimum must not be nan");
+        assert(not std::isnan(supremum_) and "supremum must not be nan");
+    }
+
+    template <DType U>
+    explicit(std::integral<T> != std::integral<U>) interval(const interval<U>& rhs) : interval() {
+        *this = rhs;
+    }
+    template <DType U>
+    explicit(std::integral<T> != std::integral<U>) interval(interval<U>&& rhs) : interval() {
+        *this = rhs;
+    }
+
+    template <DType U>
+    interval& operator=(const interval<U>& rhs) {
+        // dev note: numeric issues
+        if constexpr (std::floating_point<U> and std::integral<T>) {
+            infimum_ = std::floor(rhs.infimum());
+            supremum_ = std::ceil(rhs.supremum());
+        } else {
+            infimum_ = rhs.infimum();
+            supremum_ = rhs.supremum();
+        }
+        return *this;
+    }
+
+    /// An interval evalutes to `true` if it is not empty.
+    explicit operator bool() const noexcept { return infimum_ <= supremum_; }
+
+    bool operator==(const interval& rhs) const {
+        return infimum_ == rhs.infimum_ and supremum_ == rhs.supremum_;
+    }
+
+    interval& operator|=(DType auto scalar) { return *this |= interval<decltype(scalar)>(scalar); }
+
+    template <DType U>
+    interval& operator|=(const interval<U>& rhs) {
+        // if rhs is an empty interval, do nothing
+        if (not static_cast<bool>(rhs)) return *this;
+
+        // if lhs is an empty interval, then set it to rhs (which we know is not empty)
+        if (not static_cast<bool>(*this)) return *this = rhs;
+
+        // otherwise, adjust our interval to also contain the other.
+        // dev note: numeric issues
+        if constexpr (std::integral<T> and std::floating_point<U>) {
+            infimum_ = std::min<T>(infimum_, std::floor(rhs.infimum()));
+            supremum_ = std::max<T>(supremum_, std::ceil(rhs.supremum()));
+        } else {
+            infimum_ = std::min<T>(infimum_, rhs.infimum());
+            supremum_ = std::max<T>(supremum_, rhs.supremum());
+        }
+        return *this;
+    }
+
+    template <DType U>
+    friend auto operator|(interval lhs, const interval<U>& rhs) {
+        interval<decltype(T() + U())> out(std::move(lhs));
+        out |= rhs;
+        return out;
+    }
+
+    template <DType U>
+    interval& operator+=(const interval<U>& rhs) {
+        // if lhs is an empty interval than the result is also empty
+        if (not static_cast<bool>(*this)) return *this;
+
+        // if rhs is an empty interval, then set lhs to rhs (making lhs empty)
+        if (not static_cast<bool>(rhs)) return *this = rhs;
+
+        // neither are empty! So add everything
+        // dev note: numeric issues
+        if constexpr (std::integral<T> and std::floating_point<U>) {
+            infimum_ += std::floor(rhs.infimum());
+            supremum_ += std::ceil(rhs.supremum());
+        } else {
+            infimum_ += rhs.infimum();
+            supremum_ += rhs.supremum();
+        }
+        return *this;
+    }
+
+    template <DType U>
+    friend auto operator+(interval lhs, const interval<U>& rhs) {
+        interval<decltype(T() + U())> out(std::move(lhs));
+        out += rhs;
+        return out;
+    }
+
+    bool contains(DType auto&& value) const {
+        // dev note: numeric issues
+        return infimum_ <= value and value <= supremum_;
+    }
+
+    // Implement the get() to support the structured binding. It's not OK to
+    // overload std functions so this needs to be found via ADL.
+    // Unfortunately that prevents interval from satisfying the tuple-like
+    // concept. Why can't we have nice things C++?? See note below for
+    // more details.
+    template <std::size_t I>
+    friend const T& get(const interval& in) noexcept requires(I <= 1) {
+        if constexpr (I == 0) return in.infimum_;
+        if constexpr (I == 1) return in.supremum_;
+    }
+    template <std::size_t I>
+    friend T&& get(interval&& in) noexcept requires(I <= 1) {
+        if constexpr (I == 0) return std::move(in.infimum_);
+        if constexpr (I == 1) return std::move(in.supremum_);
+    }
+    template <std::size_t I>
+    friend const T&& get(const interval&& in) noexcept requires(I <= 1) {
+        if constexpr (I == 0) return std::move(in.infimum_);
+        if constexpr (I == 1) return std::move(in.supremum_);
+    }
+
+    const T& infimum() const noexcept { return infimum_; }
+
+    /// The size of the interval is the max(supremum_ - infimum_, 0)
+    double size() const noexcept requires(std::floating_point<T>) {
+        if (supremum_ >= infimum_) return static_cast<double>(supremum_) - infimum_;
+        return 0;
+    }
+    std::ptrdiff_t size() const noexcept requires(std::integral<T>) {
+        if (supremum_ >= infimum_) return static_cast<std::ptrdiff_t>(supremum_) - infimum_;
+        return 0;
+    }
+
+    const T& supremum() const noexcept { return supremum_; }
+
+    static interval unbounded() {
+        using limits = std::numeric_limits<T>;
+        if constexpr (limits::has_infinity) {
+            return interval(-limits::infinity(), limits::infinity());
+        } else {
+            return interval(limits::lowest(), limits::max());
+        }
+    }
+
+ private:
+    T infimum_ = 1;
+    T supremum_ = 0;
+};
+
+std::ostream& operator<<(std::ostream& os, const interval<std::int64_t>& rhs);
+std::ostream& operator<<(std::ostream& os, const interval<double>& rhs);
+
+}  // namespace dwave::optimization
+
+// Do a lot of convoluted stuff to support structured binding.
+// Unfortunately this doesn't support the full tuple-like concept that would let us work
+// out of the box with a lot of nice methods.
+// See https://stackoverflow.com/questions/79660771/stdapply-to-a-custom-tuple-like-type
+// See https://lists.isocpp.org/std-discussion/2021/09/1438.php
+namespace std {
+
+template <dwave::optimization::DType T>
+struct tuple_size<dwave::optimization::interval<T>> : integral_constant<size_t, 2> {};
+
+template <dwave::optimization::DType T>
+struct tuple_element<0, dwave::optimization::interval<T>> {
+    using type = T;
+};
+
+template <dwave::optimization::DType T>
+struct tuple_element<1, dwave::optimization::interval<T>> {
+    using type = T;
+};
+
+}  // namespace std
diff --git a/dwave/optimization/src/interval.cpp b/dwave/optimization/src/interval.cpp
new file mode 100644
index 00000000..cbfb04fa
--- /dev/null
+++ b/dwave/optimization/src/interval.cpp
@@ -0,0 +1,30 @@
+// Copyright 2026 D-Wave
+//
+//    Licensed under the Apache License, Version 2.0 (the "License");
+//    you may not use this file except in compliance with the License.
+//    You may obtain a copy of the License at
+//
+//        http://www.apache.org/licenses/LICENSE-2.0
+//
+//    Unless required by applicable law or agreed to in writing, software
+//    distributed under the License is distributed on an "AS IS" BASIS,
+//    WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+//    See the License for the specific language governing permissions and
+//    limitations under the License.
+
+#include "dwave-optimization/interval.hpp"
+#include "dwave-optimization/common.hpp"  // for ssize_t
+#include "dwave-optimization/typing.hpp"
+
+#include <iostream>
+
+namespace dwave::optimization {
+
+std::ostream& operator<<(std::ostream& os, const interval<std::int64_t>& rhs) {
+    return os << "[" << rhs.infimum() << ", " << rhs.supremum() << "]";
+}
+std::ostream& operator<<(std::ostream& os, const interval<double>& rhs) {
+    return os << "[" << rhs.infimum() << ", " << rhs.supremum() << "]";
+}
+
+}  // namespace dwave::optimization
diff --git a/meson.build b/meson.build
index 690d7165..eb2d252d 100644
--- a/meson.build
+++ b/meson.build
@@ -52,6 +52,7 @@ dwave_optimization_src = [
     'dwave/optimization/src/array.cpp',
     'dwave/optimization/src/fraction.cpp',
     'dwave/optimization/src/graph.cpp',
+    'dwave/optimization/src/interval.cpp',
     'dwave/optimization/src/simplex.cpp',
 ]
 
diff --git a/releasenotes/notes/feature-cpp-interval-232c1ed8d06f7fec.yaml b/releasenotes/notes/feature-cpp-interval-232c1ed8d06f7fec.yaml
new file mode 100644
index 00000000..d0369db2
--- /dev/null
+++ b/releasenotes/notes/feature-cpp-interval-232c1ed8d06f7fec.yaml
@@ -0,0 +1,3 @@
+---
+features:
+  - Add C++ ``interval<T>()`` class for interval arithmetic.
diff --git a/tests/cpp/meson.build b/tests/cpp/meson.build
index ef23f764..f9cb3925 100644
--- a/tests/cpp/meson.build
+++ b/tests/cpp/meson.build
@@ -35,6 +35,7 @@ tests_all = executable(
         'test_functional.cpp',
         'test_functional_.cpp',
         'test_graph.cpp',
+        'test_interval.cpp',
         'test_iterators.cpp',
         'test_simplex.cpp',
         'test_type_list.cpp',
diff --git a/tests/cpp/test_interval.cpp b/tests/cpp/test_interval.cpp
new file mode 100644
index 00000000..95dcd9dd
--- /dev/null
+++ b/tests/cpp/test_interval.cpp
@@ -0,0 +1,121 @@
+// Copyright 2026 D-Wave
+//
+//    Licensed under the Apache License, Version 2.0 (the "License");
+//    you may not use this file except in compliance with the License.
+//    You may obtain a copy of the License at
+//
+//        http://www.apache.org/licenses/LICENSE-2.0
+//
+//    Unless required by applicable law or agreed to in writing, software
+//    distributed under the License is distributed on an "AS IS" BASIS,
+//    WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+//    See the License for the specific language governing permissions and
+//    limitations under the License.
+
+#include <catch2/catch_template_test_macros.hpp>
+#include <catch2/catch_test_macros.hpp>
+#include <catch2/matchers/catch_matchers_all.hpp>
+#include <catch2/matchers/catch_matchers_range_equals.hpp>
+
+#include "dwave-optimization/interval.hpp"
+
+using Catch::Matchers::RangeEquals;
+
+namespace dwave::optimization {
+
+TEMPLATE_TEST_CASE(
+    "interval",
+    "",
+    float,
+    double,
+    bool,
+    std::int8_t,
+    std::int16_t,
+    std::int32_t,
+    std::int64_t
+) {
+    GIVEN("a default-constructed interval") {
+        auto in = interval<TestType>();
+
+        CHECK(not static_cast<bool>(in));
+        CHECK(in.size() == 0);
+    }
+
+    GIVEN("an interval containing only the value 0") {
+        auto in = interval<TestType>(0);
+
+        CHECK(static_cast<bool>(in));
+        CHECK(in.size() == 0);
+        CHECK(in.contains(0));
+        CHECK(not in.contains(1));
+        CHECK(not in.contains(std::numeric_limits<double>::epsilon()));
+    }
+
+    SECTION("get<...>(interval)") {
+        SECTION("const rvalue") {
+            const interval<TestType> in{0, 1};
+            const TestType&& inf = get<0>(std::move(in));
+            CHECK(inf == 0);
+        }
+    }
+
+    SECTION("printing") {
+        std::stringstream ss;
+        ss << interval<TestType>(0, 1);
+        CHECK(ss.str() == "[0, 1]");
+    }
+
+    SECTION("structured binding") {
+        SECTION("const reference") {
+            auto in = interval<TestType>(0, 1);
+            const auto& [inf, sup] = in;
+            STATIC_REQUIRE(std::same_as<decltype(inf), const TestType>);
+            STATIC_REQUIRE(std::same_as<decltype(sup), const TestType>);
+            CHECK(inf == 0);
+            CHECK(sup == 1);
+        }
+
+        SECTION("rvalue") {
+            auto in = interval<TestType>(0, 1);
+            auto [inf, sup] = in;
+            STATIC_REQUIRE(std::same_as<decltype(inf), TestType>);
+            STATIC_REQUIRE(std::same_as<decltype(sup), TestType>);
+            CHECK(inf == 0);
+            CHECK(sup == 1);
+        }
+
+        SECTION("const rvalue") {
+            const auto in = interval<TestType>(0, 1);
+
+            auto&& [inf, sup] = in;
+            STATIC_REQUIRE(std::same_as<decltype(inf), const TestType>);
+            STATIC_REQUIRE(std::same_as<decltype(sup), const TestType>);
+            CHECK(inf == 0);
+            CHECK(sup == 1);
+        }
+    }
+}
+
+TEST_CASE("interval") {
+    SECTION("addition") {
+        CHECK((interval<int>(0, 3) + interval<float>(.5, .5)) == interval<float>(.5, 3.5));
+    }
+
+    SECTION("union") {
+        auto in = interval<int>();
+        in |= .5;
+        CHECK(in == interval<int>(0, 1));
+        in |= -5;
+        CHECK(in == interval<int>(-5, 1));
+
+        CHECK((interval<int>(0, 1) | interval<double>(.5, 1.5)) == interval<double>(0, 1.5));
+        CHECK((interval<float>(0, 1) | interval<double>(.5, 1.5)) == interval<double>(0, 1.5));
+        CHECK((interval<float>(0, 1.5) | interval<bool>(0, 1)) == interval<float>(0, 1.5));
+    }
+
+    SECTION("Copy construction") {
+        CHECK(interval<int>(interval<float>(.2, .8)) == interval<int>(0, 1));
+    }
+}
+
+}  // namespace dwave::optimization