quad_form.hpp

Go to the documentation of this file.

/*
 * Copyright © 2025-2026 Akreon LLC
 * SPDX-License-Identifier: Apache-2.0
 */
 
#pragma once
 
#include "akreon/convex/canonical_form_builder.hpp"
#include "akreon/convex/convex_traits.hpp"
#include "akreon/convex/functions/norm.hpp"
#include "akreon/convex/functions/square.hpp"
#include "akreon/convex/sum.hpp"
#include "akreon/linalg/operator/cholesky.hpp"
#include "weighted_sum.hpp"
 
namespace convex {
 
template <linalg::Field T>
struct ScalarQuadForm {
  using Space = linalg::Matrix<T, 1, 1>;
  using Singleton = Space;
 
  template <AtomLike Var>
  constexpr auto operator()(T p, T q, T r, const Var& x) {
    static_assert(
        std::same_as<dom_space_t<meta_codom_t<atom_meta_t<Var>>>, Space>,
        "Scalar quad form requires scalar codomain");
    static_assert(is_fn_type_v<Var> || is_index_type_v<Var>,
                  "Scalar quad form must operate on scalar affine atoms");
 
    const T alpha = linalg::constexpr_sqrt_safe(T(0.5) * p);
    const auto scaled_x = x.mul(alpha);
    const auto quadratic = Square(scaled_x);
    const auto linear = x.mul(q);
 
    const auto sum = Sum(quadratic, linear);
    return sum.add(r);
  }
 
  template <AtomLike Var>
  constexpr auto operator()(T p, T q, const Var& x) {
    return this->operator()(p, q, T(0), x);
  }
};
 
template <linalg::RnSpace Space>
struct QuadForm {
  using T = linalg::space_type_t<Space>;
  static constexpr std::size_t N = linalg::space_dim_v<Space>;
 
  using Mat = linalg::Matrix<T, N, N>;
 
  using Vec = linalg::Vector<T, N>;
 
  using Singleton = linalg::Matrix<T, 1, 1>;
 
  template <AtomLike Var>
  constexpr auto operator()(const Mat& P, const Vec& q, T r, const Var& x) {
    static_assert(N > 1,
                  "Use ScalarQuadForm<T> for scalar quadratic expressions");
    static_assert(is_var_type_v<Var>, "Quad form must operate on a var type");
    static_assert(
        std::same_as<dom_space_t<meta_codom_t<atom_meta_t<Var>>>, Space>,
        "Invalid quad operator mismatch on vector types");
 
    const auto chol = linalg::Cholesky<Mat>(T(0.5) * P);
    AK_ASSERT(chol.is_valid() && chol.is_psd() && "Chol(P) is not PSD");
 
    const auto LT_x = x.mul(linalg::transpose(chol.G()));
    const auto l2 =
        Norm<NormMeta<Space, Singleton>, Var, linalg::L2<Vec>>(LT_x);
    const auto quadratic = Square(l2);
    const auto weighted_sum = WeightedSum<Space>()(q, x);
 
    const auto sum = Sum(quadratic, weighted_sum);
    return sum.add(r);
  }
 
  template <AtomLike Var>
  constexpr auto operator()(const Mat& P, const Vec& q, const Var& x) {
    return this->operator()(P, q, T(0), x);
  }
};
 
template <AtomLike Var>
constexpr auto make_scalar_quad_form(
    const dom_type_t<meta_codom_t<atom_meta_t<Var>>>& p,
    const dom_type_t<meta_codom_t<atom_meta_t<Var>>>& q,
    const dom_type_t<meta_codom_t<atom_meta_t<Var>>>& r, const Var& x) {
  using T = dom_type_t<meta_codom_t<atom_meta_t<Var>>>;
  return ScalarQuadForm<T>()(p, q, r, x);
}
 
template <AtomLike Var>
constexpr auto make_scalar_quad_form(
    const dom_type_t<meta_codom_t<atom_meta_t<Var>>>& p,
    const dom_type_t<meta_codom_t<atom_meta_t<Var>>>& q, const Var& x) {
  using T = dom_type_t<meta_codom_t<atom_meta_t<Var>>>;
  return ScalarQuadForm<T>()(p, q, x);
}
 
template <AtomLike Var>
constexpr auto make_quad_form(
    const linalg::Matrix<dom_type_t<meta_codom_t<atom_meta_t<Var>>>,
                         dom_size_v<meta_codom_t<atom_meta_t<Var>>>,
                         dom_size_v<meta_codom_t<atom_meta_t<Var>>>>& P,
    const linalg::Vector<dom_type_t<meta_codom_t<atom_meta_t<Var>>>,
                         dom_size_v<meta_codom_t<atom_meta_t<Var>>>>& q,
    const dom_type_t<meta_codom_t<atom_meta_t<Var>>>& r, const Var& x) {
  using Space = dom_space_t<meta_codom_t<atom_meta_t<Var>>>;
  return QuadForm<Space>()(P, q, r, x);
}
 
template <AtomLike Var>
constexpr auto make_quad_form(
    const linalg::Matrix<dom_type_t<meta_codom_t<atom_meta_t<Var>>>,
                         dom_size_v<meta_codom_t<atom_meta_t<Var>>>,
                         dom_size_v<meta_codom_t<atom_meta_t<Var>>>>& P,
    const linalg::Vector<dom_type_t<meta_codom_t<atom_meta_t<Var>>>,
                         dom_size_v<meta_codom_t<atom_meta_t<Var>>>>& q,
    const Var& x) {
  using Space = dom_space_t<meta_codom_t<atom_meta_t<Var>>>;
  return QuadForm<Space>()(P, q, x);
}
 
} // namespace convex