.. _program_listing_file_SeQuant_domain_mbpt_op.hpp:

Program Listing for File op.hpp
===============================

|exhale_lsh| :ref:`Return to documentation for file <file_SeQuant_domain_mbpt_op.hpp>` (``SeQuant/domain/mbpt/op.hpp``)

.. |exhale_lsh| unicode:: U+021B0 .. UPWARDS ARROW WITH TIP LEFTWARDS

.. code-block:: cpp

   //
   // Created by Eduard Valeyev on 2019-03-26.
   //
   
   #ifndef SEQUANT_DOMAIN_MBPT_OP_HPP
   #define SEQUANT_DOMAIN_MBPT_OP_HPP
   
   #include <SeQuant/domain/mbpt/fwd.hpp>
   
   #include <SeQuant/domain/mbpt/spin.hpp>
   
   #include <SeQuant/core/attr.hpp>
   #include <SeQuant/core/container.hpp>
   #include <SeQuant/core/context.hpp>
   #include <SeQuant/core/expr.hpp>
   #include <SeQuant/core/hash.hpp>
   #include <SeQuant/core/index.hpp>
   #include <SeQuant/core/interval.hpp>
   #include <SeQuant/core/math.hpp>
   #include <SeQuant/core/op.hpp>
   #include <SeQuant/core/rational.hpp>
   #include <SeQuant/core/space.hpp>
   #include <SeQuant/core/utility/strong.hpp>
   
   #include <range/v3/iterator/basic_iterator.hpp>
   #include <range/v3/range/conversion.hpp>
   #include <range/v3/range/primitives.hpp>
   #include <range/v3/view/map.hpp>
   #include <range/v3/view/transform.hpp>
   #include <range/v3/view/view.hpp>
   #include <range/v3/view/zip.hpp>
   
   #include <algorithm>
   #include <array>
   #include <cassert>
   #include <cstddef>
   #include <cstdint>
   #include <functional>
   #include <initializer_list>
   #include <iterator>
   #include <map>
   #include <optional>
   #include <stdexcept>
   #include <string>
   #include <string_view>
   #include <type_traits>
   #include <utility>
   #include <vector>
   
   namespace sequant {
   namespace mbpt {
   
   DEFINE_STRONG_TYPE_FOR_INTEGER(nₚ, std::int64_t);  // define nₚ
   DEFINE_STRONG_TYPE_FOR_INTEGER(nₕ, std::int64_t);  // define nₕ
   
   #ifndef DEFINE_SINGLE_SIGNED_ARGUMENT_OP_VARIANT
   #define DEFINE_SINGLE_SIGNED_ARGUMENT_OP_VARIANT(OP)                      \
     inline ExprPtr OP(std::int64_t Rank) { return OP(nₚ(Rank), nₕ(Rank)); } \
     inline ExprPtr OP(nₚ Rank) { return OP(Rank, nₕ(Rank)); }
   #endif  // DEFINE_SINGLE_SIGNED_ARGUMENT_OP_VARIANT
   
   template <typename QuantumNumbers>
   bool is_vacuum(QuantumNumbers qns);
   
   inline IndexSpace make_space(const IndexSpace::Type& type) {
     return get_default_context().index_space_registry()->retrieve(type,
                                                                   Spin::any);
   }
   
   enum class OpType {
     h,    
     f,    
     f̃,    
     θ,    
     g,    
     t,    
     λ,    
     A,    
     S,    
     L,    
     R,    
     R12,  
     GR,   
     C,    
     RDM,  
     RDMCumulant,  
     δ,            
     h_1,          
     t_1,          
     λ_1,          
     invalid       
   };
   
   inline const std::map<OpType, std::wstring> optype2label{
       {OpType::h, L"h"},
       {OpType::f, L"f"},
       {OpType::f̃, L"f̃"},
       {OpType::g, L"g"},
       {OpType::θ, L"θ"},
       {OpType::t, L"t"},
       {OpType::λ, L"λ"},
       {OpType::A, L"A"},
       {OpType::S, L"S"},
       {OpType::L, L"L"},
       {OpType::R, L"R"},
       {OpType::R12, L"F"},
       {OpType::GR, L"GR"},
       {OpType::C, L"C"},
       {OpType::RDM, L"γ"},
       // see https://en.wikipedia.org/wiki/Cumulant
       {OpType::RDMCumulant, L"κ"},
       {OpType::δ, L"δ"},
       {OpType::h_1, L"h¹"},
       {OpType::t_1, L"t¹"},
       {OpType::λ_1, L"λ¹"}};
   
   inline const std::map<std::wstring, OpType> label2optype =
       ranges::views::zip(ranges::views::values(optype2label),
                          ranges::views::keys(optype2label)) |
       ranges::to<std::map<std::wstring, OpType>>();
   
   enum class OpClass { ex, deex, gen };
   
   std::vector<std::wstring> cardinal_tensor_labels();
   
   std::wstring to_wstring(OpType op);
   OpClass to_class(OpType op);
   
   
   
   template <typename QuantumNumbers, Statistics S = Statistics::FermiDirac>
   class Operator;
   
   template <typename QuantumNumbers>
   using FOperator = Operator<QuantumNumbers, Statistics::FermiDirac>;
   template <typename QuantumNumbers>
   using BOperator = Operator<QuantumNumbers, Statistics::BoseEinstein>;
   
   template <Statistics S>
   class Operator<void, S> : public Expr, public Labeled {
    protected:
     Operator() = default;
     Operator(std::function<std::wstring_view()> label_generator,
              std::function<ExprPtr()> tensor_form_generator)
         : label_generator_(std::move(label_generator)),
           tensor_form_generator_(tensor_form_generator) {}
   
    public:
     virtual ~Operator() = default;
   
     std::wstring_view label() const override {
       assert(label_generator_);
       return label_generator_();
     }
   
     virtual ExprPtr tensor_form() const {
       assert(tensor_form_generator_);
       return tensor_form_generator_();
     }
   
     bool is_cnumber() const override { return false; }
   
    private:
     std::function<std::wstring_view()> label_generator_;
     std::function<ExprPtr()> tensor_form_generator_;
   };
   
   using FOperatorBase = FOperator<void>;
   using BOperatorBase = BOperator<void>;
   
   struct default_qns_tag {};
   
   // clang-format off
   
   // clang-format on
   template <typename QNV = std::int64_t, typename Tag = default_qns_tag>
   class QuantumNumberChange
       : public container::svector<
             boost::numeric::interval<std::make_signed_t<QNV>>, 8> {
    public:
     using QNC = std::make_signed_t<QNV>;  // change in quantum numbers
     using interval_t = boost::numeric::interval<QNC>;
     using base_type =
         container::svector<boost::numeric::interval<std::make_signed_t<QNV>>, 8>;
     using this_type = QuantumNumberChange<QNV, Tag>;
   
     const std::size_t size() const {
       if (get_default_context().vacuum() == Vacuum::Physical) {
         return 2;
       } else if (get_default_context().vacuum() == Vacuum::SingleProduct) {
         auto isr = get_default_context().index_space_registry();
         const auto& isr_base_spaces = isr->base_spaces();
         assert(isr_base_spaces.size() > 0);
         return isr_base_spaces.size() * 2;
       } else {
         throw std::logic_error("unknown Vacuum type");
       }
     }
   
     QuantumNumberChange() {
       this->resize(this->size());
       assert(this->base().size() != 0);
       std::fill(this->begin(), this->end(), interval_t{});
     }
   
     template <typename I, typename Range,
               typename = std::enable_if_t<
                   meta::is_range_v<std::remove_reference_t<Range>> &&
                   std::is_convertible_v<I, interval_t>>>
     explicit QuantumNumberChange(Range&& i) : QuantumNumberChange() {
       assert(i.size() == size());
       std::copy(i.begin(), i.end(), this->begin());
     }
   
     template <typename I, typename = std::enable_if_t<std::is_convertible_v<
                               std::initializer_list<I>, interval_t>>>
     explicit QuantumNumberChange(
         std::initializer_list<std::initializer_list<I>> i)
         : QuantumNumberChange() {
       assert(i.size() == size());
   #ifndef NDEBUG
       if (std::find_if(i.begin(), i.end(),
                        [](const auto& ii) { return ii.size() != 2; }) != i.end())
         throw std::invalid_argument(
             "QuantumNumberChange<N>(initializer_list<initializer_list> i): each "
             "element of i must contain 2 elements");
   #endif
       for (std::size_t c = 0; c != size(); ++c) {
         this->operator[](c) = interval_t{*((i.begin() + c)->begin()),
                                          *((i.begin() + c)->begin() + 1)};
       }
     }
   
     QuantumNumberChange& operator+=(const QuantumNumberChange& other) {
       for (std::size_t c = 0; c != size(); ++c) this->operator[](c) += other[c];
       return *this;
     }
   
     bool operator==(const this_type& b) const {
       return std::equal(
           this->begin(), this->end(), b.begin(),
           [](const auto& ia, const auto& ib) { return equal(ia, ib); });
     }
     bool operator!=(const this_type& b) const { return !this->operator==(b); }
   
     // determines the number of physical vacuum creators and annihilators for the
     // active particle and hole space from the Context. for general operators this
     // is not defined. for example: O_{e_1}^{i_1 m_1} a_{i_1 m_1}^{e_1} asking for
     // the active particle annihilators in this example is nonsense and will
     // return -1.
   
     interval_t ncre_particles() {
       const auto& qnvec = this->base();
       auto isr = get_default_context().index_space_registry();
       const auto& base_spaces = isr->base_spaces();
       interval_t result = 0;
       for (unsigned int i = 0; i < base_spaces.size(); i++) {
         const auto& base_space = base_spaces[i];
         const auto intersect_type =
             base_space.attr()
                 .intersection(isr->particle_space(base_space.qns()).attr())
                 .type();
         if (IndexSpace::Type{} != intersect_type) {
           result += qnvec[2 * i];
         }
       }
       return result;
     }
   
     interval_t nann_particles() {
       const auto& qnvec = this->base();
       auto isr = get_default_context().index_space_registry();
       const auto& base_spaces = isr->base_spaces();
       interval_t result = 0;
       for (unsigned int i = 0; i < base_spaces.size(); i++) {
         const auto& base_space = base_spaces[i];
         const auto intersect_type =
             base_space.attr()
                 .intersection(isr->particle_space(base_space.qns()).attr())
                 .type();
         if (IndexSpace::Type{} != intersect_type) {
           result += qnvec[2 * i + 1];
         }
       }
       return result;
     }
   
     interval_t ncre_holes() {
       const auto& qnvec = this->base();
       auto isr = get_default_context().index_space_registry();
       const auto& base_spaces = isr->base_spaces();
       interval_t result = 0;
       for (unsigned int i = 0; i < base_spaces.size(); i++) {
         const auto& base_space = base_spaces[i];
         const auto intersect_type =
             base_space.attr()
                 .intersection(isr->hole_space(base_space.qns()).attr())
                 .type();
         if (IndexSpace::Type{} != intersect_type) {
           result += qnvec[2 * i];
         }
       }
       return result;
     }
   
     interval_t nann_holes() {
       const auto& qnvec = this->base();
       auto isr = get_default_context().index_space_registry();
       const auto& base_spaces = isr->base_spaces();
       interval_t result = 0;
       for (unsigned int i = 0; i < base_spaces.size(); i++) {
         const auto& base_space = base_spaces[i];
         const auto intersect_type =
             base_space.attr()
                 .intersection(isr->hole_space(base_space.qns()).attr())
                 .type();
         if (IndexSpace::Type{} != intersect_type) {
           result += qnvec[2 * i + 1];
         }
       }
       return result;
     }
   
     template <typename I>
     bool in(I i) {
       return boost::numeric::in(static_cast<int64_t>(i), this->front());
     }
   
     template <typename I, typename = std::enable_if_t<std::is_integral_v<I>>>
     bool in(std::initializer_list<I> i) {
       assert(i.size() == size());
       std::array<I, 4> i_arr;
       std::copy(i.begin(), i.end(), i_arr.begin());
       return this->in(i_arr);
     }
   
     bool overlaps_with(base_type i) {
       for (std::size_t c = 0; c != this->size(); ++c) {
         if (!boost::numeric::overlap(i[c], this->operator[](c))) {
           return false;
         }
       }
       return true;
     }
   
     auto hash_value() const {
       assert(size() > 0);
       auto val = sequant::hash::value(this->operator[](0));
       for (std::size_t c = 1; c != size(); ++c) {
         sequant::hash::combine(val, sequant::hash::value(this->operator[](c)));
       }
       return val;
     }
   
    private:
     auto& base() { return static_cast<base_type&>(*this); }
   };
   
   template <std::size_t N, typename Tag, typename QNV>
   inline bool operator==(const QuantumNumberChange<Tag, QNV>& a,
                          const QuantumNumberChange<Tag, QNV>& b) {
     return a.operator==(b);
   }
   
   template <std::size_t N, typename Tag, typename QNV>
   inline bool operator!=(const QuantumNumberChange<Tag, QNV>& a,
                          const QuantumNumberChange<Tag, QNV>& b) {
     return !(a == b);
   }
   
   template <std::size_t N, typename Tag, typename QNV, typename I,
             typename = std::enable_if_t<N == 1 && std::is_integral_v<I>>>
   inline bool operator==(const QuantumNumberChange<Tag, QNV>& a, I b) {
     return a.operator==(b);
   }
   
   template <std::size_t N, typename Tag, typename QNV>
   inline bool equal(const QuantumNumberChange<Tag, QNV>& a,
                     const QuantumNumberChange<Tag, QNV>& b) {
     return operator==(a, b);
   }
   
   template <std::size_t N, typename Tag, typename QNV, typename I,
             typename = std::enable_if_t<N == 1 && std::is_integral_v<I>>>
   inline bool operator!=(const QuantumNumberChange<Tag, QNV>& a, I b) {
     return a.operator!=(b);
   }
   
   // clang-format off
   // clang-format on
   using qns_t = mbpt::QuantumNumberChange<>;
   using qninterval_t = typename qns_t::interval_t;
   using qnc_t = qns_t;
   using op_t = mbpt::Operator<qnc_t>;
   
   qns_t combine(qns_t, qns_t);
   
   qns_t excitation_type_qns(std::size_t k,
                             IndexSpace::QuantumNumbers SQN = Spin::any);
   
   qns_t interval_excitation_type_qns(std::size_t k,
                                      IndexSpace::QuantumNumbers SQN = Spin::any);
   
   qns_t deexcitation_type_qns(std::size_t k,
                               IndexSpace::QuantumNumbers SQN = Spin::any);
   
   qns_t interval_deexcitation_type_qns(
       std::size_t k, IndexSpace::QuantumNumbers SQN = Spin::any);
   
   qns_t general_type_qns(std::size_t k);
   
   qns_t generic_excitation_qns(std::size_t particle_rank, std::size_t hole_rank,
                                IndexSpace particle_space, IndexSpace hole_space,
                                IndexSpace::QuantumNumbers SQN = Spin::any);
   
   qns_t generic_deexcitation_qns(std::size_t particle_rank, std::size_t hole_rank,
                                  IndexSpace particle_space, IndexSpace hole_space,
                                  IndexSpace::QuantumNumbers SQN = Spin::any);
   
   inline namespace op {
   namespace tensor {
   ExprPtr vac_av(ExprPtr expr,
                  std::vector<std::pair<int, int>> nop_connections = {},
                  bool use_top = true);
   }  // namespace tensor
   }  // namespace op
   
   }  // namespace mbpt
   
   mbpt::qns_t adjoint(mbpt::qns_t qns);
   
   namespace mbpt {
   
   // clang-format off
   
   // clang-format on
   template <Statistics S>
   class OpMaker {
    public:
     template <typename IndexSpaceTypeRange1, typename IndexSpaceTypeRange2>
     OpMaker(OpType op, const cre<IndexSpaceTypeRange1>& cre_list,
             const ann<IndexSpaceTypeRange2>& ann_list)
         : op_(op),
           cre_spaces_(cre_list.begin(), cre_list.end()),
           ann_spaces_(ann_list.begin(), ann_list.end()) {
       assert(ncreators() > 0 || nannihilators() > 0);
     }
   
     OpMaker(OpType op, ncre nc, nann na);
   
     OpMaker(OpType op, std::size_t rank);
   
     OpMaker(OpType op, ncre nc, nann na, const cre<IndexSpace>& cre_space,
             const ann<IndexSpace>& ann_space);
   
     enum class UseDepIdx {
       Bra,
       Ket,
       None
     };
   
     // clang-format off
     // clang-format on
     ExprPtr operator()(std::optional<UseDepIdx> dep_opt = {},
                        std::optional<Symmetry> opsymm_opt = {}) const;
   
     template <typename TensorGenerator>
     static ExprPtr make(const container::svector<IndexSpace>& creators,
                         const container::svector<IndexSpace>& annihilators,
                         TensorGenerator&& tensor_generator,
                         UseDepIdx dep = UseDepIdx::None) {
       const bool symm =
           get_default_context().spbasis() ==
           SPBasis::spinorbital;  // antisymmetrize if spin-orbital basis
       const auto dep_bra = dep == UseDepIdx::Bra;
       const auto dep_ket = dep == UseDepIdx::Ket;
   
       // not sure what it means to use nonsymmetric operator if nbra != nket
       if (!symm) assert(ranges::size(creators) == ranges::size(annihilators));
   
       auto make_idx_vector = [](const auto& spacetypes) {
         container::svector<Index> result;
         const auto n = spacetypes.size();
         result.reserve(n);
         for (size_t i = 0; i != n; ++i) {
           auto space = spacetypes[i];
           result.push_back(Index::make_tmp_index(space));
         }
         return result;
       };
   
       auto make_depidx_vector = [](const auto& spacetypes, auto&& protoidxs) {
         const auto n = spacetypes.size();
         container::svector<Index> result;
         result.reserve(n);
         for (size_t i = 0; i != n; ++i) {
           auto space = spacetypes[i];
           result.push_back(Index::make_tmp_index(space, protoidxs, true));
         }
         return result;
       };
   
       container::svector<Index> creidxs, annidxs;
       if (dep_bra) {
         annidxs = make_idx_vector(annihilators);
         creidxs = make_depidx_vector(creators, annidxs);
       } else if (dep_ket) {
         creidxs = make_idx_vector(creators);
         annidxs = make_depidx_vector(annihilators, creidxs);
       } else {
         creidxs = make_idx_vector(creators);
         annidxs = make_idx_vector(annihilators);
       }
   
       using ranges::size;
       const auto mult =
           symm ? factorial(size(creators)) * factorial(size(annihilators))
                : factorial(size(creators));
       const auto opsymm = symm ? (S == Statistics::FermiDirac ? Symmetry::antisymm
                                                               : Symmetry::symm)
                                : Symmetry::nonsymm;
       return ex<Constant>(rational{1, mult}) *
              tensor_generator(creidxs, annidxs, opsymm) *
              ex<NormalOperator<S>>(cre(creidxs), ann(annidxs),
                                    get_default_context().vacuum());
     }
   
     template <typename TensorGenerator>
     static ExprPtr make(std::initializer_list<IndexSpace::Type> creators,
                         std::initializer_list<IndexSpace::Type> annihilators,
                         TensorGenerator&& tensor_generator,
                         UseDepIdx csv = UseDepIdx::None) {
       container::svector<IndexSpace> cre_vec(creators.begin(), creators.end());
       container::svector<IndexSpace> ann_vec(annihilators.begin(),
                                              annihilators.end());
       return OpMaker::make(cre_vec, ann_vec,
                            std::forward<TensorGenerator>(tensor_generator), csv);
     }
   
    protected:
     OpType op_;
     container::svector<IndexSpace> cre_spaces_;
     container::svector<IndexSpace> ann_spaces_;
   
     OpMaker(OpType op);
   
     [[nodiscard]] const auto ncreators() const { return cre_spaces_.size(); };
     [[nodiscard]] const auto nannihilators() const { return ann_spaces_.size(); };
   };
   
   extern template class OpMaker<Statistics::FermiDirac>;
   extern template class OpMaker<Statistics::BoseEinstein>;
   
   template <typename QuantumNumbers, Statistics S>
   class Operator : public Operator<void, S> {
     using this_type = Operator<QuantumNumbers, S>;
     using base_type = Operator<void, S>;
   
    protected:
     Operator();
   
    public:
     Operator(std::function<std::wstring_view()> label_generator,
              std::function<ExprPtr()> tensor_form_generator,
              std::function<void(QuantumNumbers&)> qn_action);
   
     virtual ~Operator();
   
     QuantumNumbers operator()(const QuantumNumbers& qns = {}) const;
   
     template <typename I, typename = std::enable_if_t<std::is_integral_v<I>>>
     QuantumNumbers operator()(std::initializer_list<I> qns) const {
       QuantumNumbers result(qns);
       this->apply_to(result);
       return result;
     }
   
     virtual QuantumNumbers& apply_to(QuantumNumbers& qns) const;
   
     bool static_less_than(const Expr& that) const override;
   
     bool commutes_with_atom(const Expr& that) const override;
   
     void adjoint() override;
   
     bool is_adjoint_ = false;
   
    private:
     std::function<void(QuantumNumbers&)> qn_action_;
   
     bool less_than_rank_of(const this_type& that) const;
   
     Expr::type_id_type type_id() const override;
   
     ExprPtr clone() const override;
   
     std::wstring to_latex() const override;
   
     Expr::hash_type memoizing_hash() const override;
   
   };  // class Operator
   
   extern template class Operator<qns_t, Statistics::FermiDirac>;
   extern template class Operator<qns_t, Statistics::BoseEinstein>;
   
   inline namespace op {
   namespace tensor {
   
   using mbpt::nₕ;
   using mbpt::nₚ;
   
   // clang-format off
   // clang-format on
   ExprPtr H_(std::size_t k);
   
   ExprPtr H(std::size_t k = 2);
   
   ExprPtr F(bool use_tensor = true, IndexSpace reference_occupied = {L"", 0});
   
   ExprPtr θ(std::size_t K);
   
   ExprPtr T_(std::size_t K);
   
   ExprPtr T(std::size_t K, bool skip1 = false);
   
   ExprPtr Λ_(std::size_t K);
   
   ExprPtr Λ(std::size_t K);
   
   ExprPtr R_(
       nann na, ncre nc,
       const cre<IndexSpace>& cre_space = cre(get_particle_space(Spin::any)),
       const ann<IndexSpace>& ann_space = ann(get_hole_space(Spin::any)));
   
   ExprPtr R_(nₚ np, nₕ nh);
   DEFINE_SINGLE_SIGNED_ARGUMENT_OP_VARIANT(R_);
   
   ExprPtr L_(
       nann na, ncre nc,
       const cre<IndexSpace>& cre_space = cre(get_hole_space(Spin::any)),
       const ann<IndexSpace>& ann_space = ann(get_particle_space(Spin::any)));
   
   ExprPtr L_(nₚ np, nₕ nh);
   DEFINE_SINGLE_SIGNED_ARGUMENT_OP_VARIANT(L_);
   
   // clang-format off
   // clang-format on
   ExprPtr P(nₚ np, nₕ nh);
   DEFINE_SINGLE_SIGNED_ARGUMENT_OP_VARIANT(P);
   
   // clang-format off
   // clang-format on
   ExprPtr A(nₚ np, nₕ nh);
   DEFINE_SINGLE_SIGNED_ARGUMENT_OP_VARIANT(A);
   
   ExprPtr S(std::int64_t K);
   
   ExprPtr H_pt(std::size_t order, std::size_t R);
   
   ExprPtr T_pt_(std::size_t order, std::size_t K);
   
   ExprPtr T_pt(std::size_t order, std::size_t K, bool skip1 = false);
   
   ExprPtr Λ_pt_(std::size_t order, std::size_t K);
   
   ExprPtr Λ_pt(std::size_t order, std::size_t K, bool skip1 = false);
   }  // namespace tensor
   }  // namespace op
   
   inline namespace op {
   // clang-format off
   // clang-format on
   ExprPtr H_(std::size_t k);
   
   ExprPtr H(std::size_t k = 2);
   
   ExprPtr F(bool use_tensor = true, IndexSpace reference_occupied = {L"", 0});
   
   ExprPtr θ(std::size_t K);
   
   ExprPtr T_(std::size_t K);
   
   ExprPtr T(std::size_t K, bool skip1 = false);
   
   ExprPtr Λ_(std::size_t K);
   
   ExprPtr Λ(std::size_t K);
   
   ExprPtr R_(
       nann na, ncre nc,
       const cre<IndexSpace>& cre_space = cre(get_particle_space(Spin::any)),
       const ann<IndexSpace>& ann_space = ann(get_hole_space(Spin::any)));
   
   ExprPtr R_(nₚ np, nₕ nh);
   DEFINE_SINGLE_SIGNED_ARGUMENT_OP_VARIANT(R_);
   
   ExprPtr L_(
       nann na, ncre nc,
       const cre<IndexSpace>& cre_space = cre(get_hole_space(Spin::any)),
       const ann<IndexSpace>& ann_space = ann(get_particle_space(Spin::any)));
   
   ExprPtr L_(nₚ np, nₕ nh);
   DEFINE_SINGLE_SIGNED_ARGUMENT_OP_VARIANT(L_);
   
   ExprPtr R(nann na, ncre nc,
             const cre<IndexSpace>& cre_space = cre(get_particle_space(Spin::any)),
             const ann<IndexSpace>& ann_space = ann(get_hole_space(Spin::any)));
   
   ExprPtr R(nₚ np, nₕ nh);
   DEFINE_SINGLE_SIGNED_ARGUMENT_OP_VARIANT(R);
   
   ExprPtr L(
       nann na, ncre nc,
       const cre<IndexSpace>& cre_space = cre(get_hole_space(Spin::any)),
       const ann<IndexSpace>& ann_space = ann(get_particle_space(Spin::any)));
   
   ExprPtr L(nₚ np, nₕ nh);
   DEFINE_SINGLE_SIGNED_ARGUMENT_OP_VARIANT(L);
   
   // clang-format off
   // clang-format on
   ExprPtr P(nₚ np, nₕ nh);
   DEFINE_SINGLE_SIGNED_ARGUMENT_OP_VARIANT(P);
   
   // clang-format off
   // clang-format on
   ExprPtr A(nₚ np, nₕ nh);
   DEFINE_SINGLE_SIGNED_ARGUMENT_OP_VARIANT(A);
   
   ExprPtr S(std::int64_t K);
   
   ExprPtr H_pt(std::size_t order, std::size_t R);
   
   ExprPtr T_pt_(std::size_t order, std::size_t K);
   
   ExprPtr T_pt(std::size_t order, std::size_t K, bool skip1 = false);
   
   ExprPtr Λ_pt_(std::size_t order, std::size_t K);
   
   ExprPtr Λ_pt(std::size_t order, std::size_t K, bool skip1 = false);
   
   bool raises_vacuum_up_to_rank(const ExprPtr& op_or_op_product,
                                 const unsigned long k);
   
   bool lowers_rank_or_lower_to_vacuum(const ExprPtr& op_or_op_product,
                                       const unsigned long k);
   
   bool raises_vacuum_to_rank(const ExprPtr& op_or_op_product,
                              const unsigned long k);
   
   bool lowers_rank_to_vacuum(const ExprPtr& op_or_op_product,
                              const unsigned long k);
   #include <SeQuant/domain/mbpt/vac_av.hpp>
   }  // namespace op
   }  // namespace mbpt
   }  // namespace sequant
   
   #endif  // SEQUANT_DOMAIN_MBPT_OP_HPP