halfspace.hpp

Go to the documentation of this file.

/*
 * Copyright © 2025-2026 Akreon LLC
 * SPDX-License-Identifier: Apache-2.0
 */
 
#pragma once
 
#include "convex_set.hpp"
 
namespace convex {
 
template <linalg::RnSpace Space>
using HalfspaceDom = ConvexDomain<Space, VarType::FREE>;
 
template <linalg::RnSpace Space>
using HalfspaceCodom = HalfspaceDom<Space>;
 
template <linalg::RnSpace Space>
using HalfspaceMeta = ConvexAtomMeta<HalfspaceDom<Space>, HalfspaceCodom<Space>,
                                     CurvatureType::CONVEX>;
 
template <typename Meta, typename Child, HalfspaceType Type>
struct Halfspace : ConvexSet<Meta, Halfspace<Meta, Child, Type>, Child> {
  using ConvexBase = ConvexSet<Meta, Halfspace<Meta, Child, Type>, Child>;
  using typename ConvexBase::PointT;
  using typename ConvexBase::T;
 
  using ConvexBase::N;
 
  constexpr explicit Halfspace(const PointT& a_in, T b_in,
                               const Child& child_in)
      : ConvexBase(child_in), a(a_in), b(b_in) {}
 
  constexpr bool contains_impl(const PointT& x, T /* tol */) const {
    switch (Type) {
      case HalfspaceType::UH:
        return linalg::inner_product(a, x) >= b;
      case HalfspaceType::LH:
      default:
        return linalg::inner_product(a, x) <= b;
    }
  }
 
  constexpr bool is_well_formed() const {
    if (linalg::all_zero(a)) {
      if constexpr (Type == HalfspaceType::LH) {
        if (b < T(0)) return false;
      } else {
        if (b > T(0)) return false;
      }
    }
    return true;
  }
 
  template <linalg::SquareMatrix Mat>
    requires(std::same_as<T, linalg::matrix_scalar_t<Mat>> &&
             linalg::mat_cols_v<Mat> == N)
  constexpr auto transform_impl(const affine::AffineMap<Mat>& map) const {
    const auto [U_inv, event] = linalg::MatrixInverse<Mat>()(map.linear_part());
    static_cast<void>(event);
    const auto aT_U = transpose(as_col_vector(a)) * U_inv;
    return Halfspace<Meta, Child, Type>(as_vector(aT_U),
                                        b + (aT_U * map.translation())(0),
                                        std::get<0>(this->children));
  }
 
  template <typename Other>
  constexpr bool is_equal_impl(const Other& other) const {
    return other.self().a == a && other.self().b == b;
  }
 
  template <AtomType AType, conic::BuilderLike Builder>
  constexpr auto canonicalize_set(const Builder& builder) const {
    const auto child_builder =
        std::get<0>(this->children)
            .template canonicalize<AtomType::SET>(builder);
    const auto child_coord = child_builder.coord_map;
 
    static constexpr std::size_t ChildCodom =
        linalg::mat_cols_v<decltype(child_coord.linear_part())>;
    static_assert(ChildCodom == N,
                  "Convex set domain does not equal the child coordinate "
                  "linear map's codomain in size");
 
    linalg::RowVector<T, N> F;
    F = linalg::set_block<0, 1, 0, N>(F, T(-1) * linalg::as_row_vector(a));
    linalg::Vector<T, 1> g{b};
    const auto map = affine::AffineMap<decltype(F)>(F, g);
    const auto rplus = conic::NonnegOrthantCone<linalg::Vector<T, 1>>();
    const auto cones = std::make_tuple(rplus);
    const auto emission =
        conic::LocalEmission(map, cones, child_builder.eq_map);
 
    return child_builder.template realize<AType>(emission);
  }
 
  PointT a; //< n-dim vector of type T
  T b;      //< bias scalar of type T
};
 
} // namespace convex