diff --git a/.github/workflows/cmake.yml b/.github/workflows/cmake.yml
index 2459b5937f..d1a2c69203 100644
--- a/.github/workflows/cmake.yml
+++ b/.github/workflows/cmake.yml
@@ -125,7 +125,7 @@ jobs:
       working-directory: ${{github.workspace}}/build
       shell: bash
       # Execute the build.  You can specify a specific target with "--target <NAME>"
-      run: ccache -p && ccache -z && cmake --build . && ccache -s
+      run: ccache -p && ccache -z && cmake --build . -j 2 && ccache -s
 
     - name: Test
       if: ${{ !matrix.valgrind }}
diff --git a/SeQuant/domain/mbpt/models/cc.cpp b/SeQuant/domain/mbpt/models/cc.cpp
index def47f38f4..1eae80797b 100644
--- a/SeQuant/domain/mbpt/models/cc.cpp
+++ b/SeQuant/domain/mbpt/models/cc.cpp
@@ -89,8 +89,7 @@ std::vector<ExprPtr> CC::λ(size_t commutator_rank) {
   auto hbar = mbpt::lst(H(), T(N, skip_singles), commutator_rank - 1,
                         {.unitary = unitary()});
 
-  const auto One = ex<Constant>(1);
-  auto lhbar = simplify((One + Λ(N)) * hbar);
+  auto lhbar = simplify((1 + Λ(N)) * hbar);
 
   const auto op_connect = concat(default_op_connections(),
                                  OpConnections<OpType>{{OpType::h, OpType::A},
@@ -139,13 +138,13 @@ std::vector<ExprPtr> CC::λ(size_t commutator_rank) {
   return result;
 }
 
-std::vector<ExprPtr> CC::t_pt(size_t rank, size_t order,
-                              std::optional<size_t> nbatch) {
+std::vector<ExprPtr> CC::tʼ(size_t rank, size_t order,
+                            std::optional<size_t> nbatch) {
   SEQUANT_ASSERT(order == 1 &&
-                 "sequant::mbpt::CC::t_pt(): only first-order perturbation is "
+                 "sequant::mbpt::CC::tʼ(): only first-order perturbation is "
                  "supported now");
   SEQUANT_ASSERT(rank == 1 &&
-                 "sequant::mbpt::CC::t_pt(): only one-body perturbation "
+                 "sequant::mbpt::CC::tʼ(): only one-body perturbation "
                  "operator is supported now");
   SEQUANT_ASSERT(ansatz_ == Ansatz::T && "unitary ansatz is not yet supported");
 
@@ -153,12 +152,12 @@ std::vector<ExprPtr> CC::t_pt(size_t rank, size_t order,
   // truncate h1_bar at rank 2 for one-body perturbation
   // operator and at rank 4 for two-body perturbation operator
   const auto h1_truncate_at = rank == 1 ? 2 : 4;
-  const auto h1_bar = mbpt::lst(H_pt(rank, {.order = order, .nbatch = nbatch}),
+  const auto h1_bar = mbpt::lst(Hʼ(rank, {.order = order, .nbatch = nbatch}),
                                 T(N), h1_truncate_at);
 
   // construct [hbar, T(1)]
   const auto hbar_pert =
-      mbpt::lst(H(), T(N), 3) * T_pt(N, {.order = order, .nbatch = nbatch});
+      mbpt::lst(H(), T(N), 3) * Tʼ(N, {.order = order, .nbatch = nbatch});
 
   // [Eq. 34, WIREs Comput Mol Sci. 2019; 9:e1406]
   const auto expr = simplify(h1_bar + hbar_pert);
@@ -175,21 +174,21 @@ std::vector<ExprPtr> CC::t_pt(size_t rank, size_t order,
 
   std::vector<ExprPtr> result(N + 1);
   for (auto p = N; p >= 1; --p) {
-    const auto freq_term = ex<Variable>(L"ω") * P(nₚ(p)) *
-                           T_pt_(p, {.order = order, .nbatch = nbatch});
+    const auto freq_term =
+        L"ω" * P(nₚ(p)) * op::tʼ(p, {.order = order, .nbatch = nbatch});
     result.at(p) =
         this->ref_av(P(nₚ(p)) * expr, op_connect) - this->ref_av(freq_term);
   }
   return result;
 }
 
-std::vector<ExprPtr> CC::λ_pt(size_t rank, size_t order,
-                              std::optional<size_t> nbatch) {
+std::vector<ExprPtr> CC::λʼ(size_t rank, size_t order,
+                            std::optional<size_t> nbatch) {
   SEQUANT_ASSERT(order == 1 &&
-                 "sequant::mbpt::CC::λ_pt(): only first-order perturbation is "
+                 "sequant::mbpt::CC::λʼ(): only first-order perturbation is "
                  "supported now");
   SEQUANT_ASSERT(rank == 1 &&
-                 "sequant::mbpt::CC::λ_pt(): only one-body perturbation "
+                 "sequant::mbpt::CC::λʼ(): only one-body perturbation "
                  "operator is supported now");
   SEQUANT_ASSERT(ansatz_ == Ansatz::T && "unitary ansatz is not yet supported");
 
@@ -200,18 +199,16 @@ std::vector<ExprPtr> CC::λ_pt(size_t rank, size_t order,
   // truncate h1_bar at rank 2 for one-body perturbation
   // operator and at rank 4 for two-body perturbation operator
   const auto h1_truncate_at = rank == 1 ? 2 : 4;
-  const auto h1_bar = mbpt::lst(H_pt(rank, {.order = order, .nbatch = nbatch}),
+  const auto h1_bar = mbpt::lst(Hʼ(rank, {.order = order, .nbatch = nbatch}),
                                 T(N), h1_truncate_at);
 
   // construct [hbar, T(1)]
   const auto hbar_pert =
-      mbpt::lst(H(), T(N), 3) * T_pt(N, {.order = order, .nbatch = nbatch});
+      mbpt::lst(H(), T(N), 3) * Tʼ(N, {.order = order, .nbatch = nbatch});
 
   // [Eq. 35, WIREs Comput Mol Sci. 2019; 9:e1406]
-  const auto One = ex<Constant>(1);
-  const auto expr =
-      simplify((One + Λ(N)) * (h1_bar + hbar_pert) +
-               Λ_pt(N, {.order = order, .nbatch = nbatch}) * hbar);
+  const auto expr = simplify((1 + Λ(N)) * (h1_bar + hbar_pert) +
+                             Λʼ(N, {.order = order, .nbatch = nbatch}) * hbar);
 
   // connectivity:
   // t and t1 with {h,f,g}
@@ -235,9 +232,8 @@ std::vector<ExprPtr> CC::λ_pt(size_t rank, size_t order,
 
   std::vector<ExprPtr> result(N + 1);
   for (auto p = N; p >= 1; --p) {
-    const auto freq_term = ex<Variable>(L"ω") *
-                           Λ_pt_(p, {.order = order, .nbatch = nbatch}) *
-                           P(nₚ(-p));
+    const auto freq_term =
+        L"ω" * op::λʼ(p, {.order = order, .nbatch = nbatch}) * P(nₚ(-p));
     result.at(p) =
         this->ref_av(expr * P(nₚ(-p)), op_connect) + this->ref_av(freq_term);
   }
diff --git a/SeQuant/domain/mbpt/models/cc.hpp b/SeQuant/domain/mbpt/models/cc.hpp
index 0759d60f2b..45bc16b698 100644
--- a/SeQuant/domain/mbpt/models/cc.hpp
+++ b/SeQuant/domain/mbpt/models/cc.hpp
@@ -88,7 +88,7 @@ class CC {
   /// @pre `rank==1 && order==1`, only first order perturbation and one-body perturbation operator is supported now
   /// @return std::vector of perturbed t amplitude equations
   // clang-format on
-  [[nodiscard]] std::vector<ExprPtr> t_pt(
+  [[nodiscard]] std::vector<ExprPtr> tʼ(
       size_t rank = 1, size_t order = 1,
       std::optional<size_t> nbatch = std::nullopt);
 
@@ -100,7 +100,7 @@ class CC {
   /// @pre `rank==1 && order==1`, only first order perturbation and one-body perturbation operator is supported now
   /// @return std::vector of perturbed λ amplitude equations
   // clang-format on
-  [[nodiscard]] std::vector<ExprPtr> λ_pt(
+  [[nodiscard]] std::vector<ExprPtr> λʼ(
       size_t rank = 1, size_t order = 1,
       std::optional<size_t> nbatch = std::nullopt);
 
diff --git a/SeQuant/domain/mbpt/op.cpp b/SeQuant/domain/mbpt/op.cpp
index b2f87195e9..63ec23335a 100644
--- a/SeQuant/domain/mbpt/op.cpp
+++ b/SeQuant/domain/mbpt/op.cpp
@@ -613,7 +613,7 @@ template class Operator<qns_t, Statistics::BoseEinstein>;
 inline namespace op {
 
 namespace tensor {
-ExprPtr H_(std::size_t k) {
+ExprPtr h(std::size_t k) {
   SEQUANT_ASSERT(k > 0 && k <= 2);
   switch (k) {
     case 1:
@@ -636,7 +636,7 @@ ExprPtr H_(std::size_t k) {
 
 ExprPtr H(std::size_t k) {
   SEQUANT_ASSERT(k > 0 && k <= 2);
-  return k == 1 ? tensor::H_(1) : tensor::H_(1) + tensor::H_(2);
+  return k == 1 ? tensor::h(1) : tensor::h(1) + tensor::h(2);
 }
 
 ExprPtr F(bool use_tensor, IndexSpace reference_occupied) {
@@ -695,7 +695,7 @@ ExprPtr θ(std::size_t K) {
   return OpMaker<Statistics::FermiDirac>(OpType::θ, K)();
 }
 
-ExprPtr T_(std::size_t K) {
+ExprPtr t(std::size_t K) {
   return OpMaker<Statistics::FermiDirac>(OpType::t, K)();
 }
 
@@ -703,12 +703,12 @@ ExprPtr T(std::size_t K, bool skip1) {
   SEQUANT_ASSERT(K > (skip1 ? 1 : 0));
   ExprPtr result;
   for (auto k = skip1 ? 2ul : 1ul; k <= K; ++k) {
-    result += tensor::T_(k);
+    result += tensor::t(k);
   }
   return result;
 }
 
-ExprPtr Λ_(std::size_t K) {
+ExprPtr λ(std::size_t K) {
   return OpMaker<Statistics::FermiDirac>(OpType::λ, K)();
 }
 
@@ -717,29 +717,29 @@ ExprPtr Λ(std::size_t K, bool skip1) {
 
   ExprPtr result;
   for (auto k = (skip1 ? 2ul : 1ul); k <= K; ++k) {
-    result = k > 1 ? result + tensor::Λ_(k) : tensor::Λ_(k);
+    result = k > 1 ? result + tensor::λ(k) : tensor::λ(k);
   }
   return result;
 }
 
-ExprPtr R_(nann na, ncre nc, const cre<IndexSpace>& cre_space,
-           const ann<IndexSpace>& ann_space) {
+ExprPtr r(nann na, ncre nc, const cre<IndexSpace>& cre_space,
+          const ann<IndexSpace>& ann_space) {
   return OpMaker<Statistics::FermiDirac>(OpType::R, nc, na, cre_space,
                                          ann_space)();
 }
-ExprPtr R_(nₚ np, nₕ nh) {
+ExprPtr r(nₚ np, nₕ nh) {
   SEQUANT_ASSERT(np >= 0 && nh >= 0);
   return OpMaker<Statistics::FermiDirac>(OpType::R, ncre(np.value()),
                                          nann(nh.value()))();
 }
 
-ExprPtr L_(nann na, ncre nc, const cre<IndexSpace>& cre_space,
-           const ann<IndexSpace>& ann_space) {
+ExprPtr l(nann na, ncre nc, const cre<IndexSpace>& cre_space,
+          const ann<IndexSpace>& ann_space) {
   return OpMaker<Statistics::FermiDirac>(OpType::L, nc, na, cre_space,
                                          ann_space)();
 }
 
-ExprPtr L_(nₚ np, nₕ nh) {
+ExprPtr l(nₚ np, nₕ nh) {
   SEQUANT_ASSERT(np >= 0 && nh >= 0);
   return OpMaker<Statistics::FermiDirac>(OpType::L, ncre(nh.value()),
                                          nann(np.value()))();
@@ -811,65 +811,65 @@ ExprPtr S(std::int64_t K) {
       OpType::S, cre(creators), ann(annihilators))(dep, {Symmetry::Nonsymm});
 }
 
-ExprPtr H_pt(std::size_t R, const OpParams& params) {
+ExprPtr Hʼ(std::size_t R, const OpParams& params) {
   params.validate();
   SEQUANT_ASSERT(params.order == 1 &&
-                 "sequant::mbpt::tensor::H_pt: only supports first "
+                 "sequant::mbpt::tensor::Hʼ: only supports first "
                  "order perturbation");
   SEQUANT_ASSERT(R > 0 && "Operator rank must be > 0");
   return OpMaker<Statistics::FermiDirac>(OpType::h_1, ncre(R), nann(R),
                                          params)();
 }
 
-ExprPtr T_pt_(std::size_t K, const OpParams& params) {
+ExprPtr tʼ(std::size_t K, const OpParams& params) {
   params.validate();
   SEQUANT_ASSERT(params.order == 1 &&
-                 "sequant::mbpt::tensor::T_pt_: only supports first "
+                 "sequant::mbpt::tensor::tʼ: only supports first "
                  "order perturbation");
   SEQUANT_ASSERT(K > 0 && "Operator rank must be > 0");
   return OpMaker<Statistics::FermiDirac>(OpType::t_1, ncre(K), nann(K),
                                          params)();
 }
 
-ExprPtr T_pt(std::size_t K, const OpParams& params) {
+ExprPtr Tʼ(std::size_t K, const OpParams& params) {
   params.validate();
   if (params.skip1) SEQUANT_ASSERT(K > 1);
   ExprPtr result;
   for (auto k = (params.skip1 ? 2ul : 1ul); k <= K; ++k) {
-    result += tensor::T_pt_(k, {.order = params.order,
-                                .nbatch = params.nbatch,
-                                .batch_ordinals = params.batch_ordinals,
-                                .skip1 = false});
+    result += tensor::tʼ(k, {.order = params.order,
+                             .nbatch = params.nbatch,
+                             .batch_ordinals = params.batch_ordinals,
+                             .skip1 = false});
   }
   return result;
 }
 
-ExprPtr Λ_pt_(std::size_t K, const OpParams& params) {
+ExprPtr λʼ(std::size_t K, const OpParams& params) {
   params.validate();
   SEQUANT_ASSERT(params.order == 1 &&
-                 "sequant::mbpt::tensor::Λ_pt_: only supports first "
+                 "sequant::mbpt::tensor::Λʼ: only supports first "
                  "order perturbation");
   SEQUANT_ASSERT(K > 0 && "Operator rank must be > 0");
   return OpMaker<Statistics::FermiDirac>(OpType::λ_1, ncre(K), nann(K),
                                          params)();
 }
 
-ExprPtr Λ_pt(std::size_t K, const OpParams& params) {
+ExprPtr Λʼ(std::size_t K, const OpParams& params) {
   params.validate();
   if (params.skip1) SEQUANT_ASSERT(K > 1);
   ExprPtr result;
   for (auto k = (params.skip1 ? 2ul : 1ul); k <= K; ++k) {
-    result += tensor::Λ_pt_(k, {.order = params.order,
-                                .nbatch = params.nbatch,
-                                .batch_ordinals = params.batch_ordinals,
-                                .skip1 = false});
+    result += tensor::λʼ(k, {.order = params.order,
+                             .nbatch = params.nbatch,
+                             .batch_ordinals = params.batch_ordinals,
+                             .skip1 = false});
   }
   return result;
 }
 
 }  // namespace tensor
 
-ExprPtr H_(std::size_t k) {
+ExprPtr h(std::size_t k) {
   SEQUANT_ASSERT(k > 0 && k <= 2);
   switch (k) {
     case 1:
@@ -886,7 +886,7 @@ ExprPtr H_(std::size_t k) {
 
             SEQUANT_UNREACHABLE;
           },
-          [=]() -> ExprPtr { return tensor::H_(1); },
+          [=]() -> ExprPtr { return tensor::h(1); },
           [=](qnc_t& qns) {
             qnc_t op_qnc_t = general_type_qns(1);
             qns = combine(op_qnc_t, qns);
@@ -894,7 +894,7 @@ ExprPtr H_(std::size_t k) {
 
     case 2:
       return ex<op_t>([]() -> std::wstring_view { return L"g"; },
-                      [=]() -> ExprPtr { return tensor::H_(2); },
+                      [=]() -> ExprPtr { return tensor::h(2); },
                       [=](qnc_t& qns) {
                         qnc_t op_qnc_t = general_type_qns(2);
                         qns = combine(op_qnc_t, qns);
@@ -906,7 +906,7 @@ ExprPtr H_(std::size_t k) {
 
 ExprPtr H(std::size_t k) {
   SEQUANT_ASSERT(k > 0 && k <= 2);
-  return k == 1 ? H_(1) : H_(1) + H_(2);
+  return k == 1 ? h(1) : h(1) + h(2);
 }
 
 ExprPtr θ(std::size_t K) {
@@ -919,10 +919,10 @@ ExprPtr θ(std::size_t K) {
                   });
 }
 
-ExprPtr T_(std::size_t K) {
+ExprPtr t(std::size_t K) {
   SEQUANT_ASSERT(K > 0);
   return ex<op_t>([]() -> std::wstring_view { return L"t"; },
-                  [=]() -> ExprPtr { return tensor::T_(K); },
+                  [=]() -> ExprPtr { return tensor::t(K); },
                   [=](qnc_t& qns) {
                     qnc_t op_qnc_t = excitation_type_qns(K);
                     qns = combine(op_qnc_t, qns);
@@ -933,15 +933,15 @@ ExprPtr T(std::size_t K, bool skip1) {
   SEQUANT_ASSERT(K > (skip1 ? 1 : 0));
   ExprPtr result;
   for (auto k = skip1 ? 2ul : 1ul; k <= K; ++k) {
-    result += T_(k);
+    result += t(k);
   }
   return result;
 }
 
-ExprPtr Λ_(std::size_t K) {
+ExprPtr λ(std::size_t K) {
   SEQUANT_ASSERT(K > 0);
   return ex<op_t>([]() -> std::wstring_view { return L"λ"; },
-                  [=]() -> ExprPtr { return tensor::Λ_(K); },
+                  [=]() -> ExprPtr { return tensor::λ(K); },
                   [=](qnc_t& qns) {
                     qnc_t op_qnc_t = deexcitation_type_qns(K);
                     qns = combine(op_qnc_t, qns);
@@ -952,7 +952,7 @@ ExprPtr Λ(std::size_t K, bool skip1) {
   SEQUANT_ASSERT(K > (skip1 ? 1 : 0));
   ExprPtr result;
   for (auto k = (skip1 ? 2ul : 1ul); k <= K; ++k) {
-    result = k > 1 ? result + Λ_(k) : Λ_(k);
+    result = k > 1 ? result + λ(k) : λ(k);
   }
   return result;
 }
@@ -1053,71 +1053,71 @@ ExprPtr make_op_from_params(std::function<std::wstring_view()> label_gen,
 }
 }  // anonymous namespace
 
-ExprPtr H_pt(std::size_t R, const OpParams& params) {
+ExprPtr Hʼ(std::size_t R, const OpParams& params) {
   SEQUANT_ASSERT(R > 0);
   SEQUANT_ASSERT(params.order == 1 &&
                  "only first order perturbation is supported now");
 
   return make_op_from_params(
       []() -> std::wstring_view { return optype2label.at(OpType::h_1); },
-      [R, params]() -> ExprPtr { return tensor::H_pt(R, params); },
+      [R, params]() -> ExprPtr { return tensor::Hʼ(R, params); },
       [R](qnc_t& qns) { qns = combine(general_type_qns(R), qns); }, params);
 }
 
-ExprPtr T_pt_(std::size_t K, const OpParams& params) {
+ExprPtr tʼ(std::size_t K, const OpParams& params) {
   SEQUANT_ASSERT(K > 0);
   SEQUANT_ASSERT(params.order == 1 &&
                  "only first order perturbation is supported now");
 
   return make_op_from_params(
       []() -> std::wstring_view { return optype2label.at(OpType::t_1); },
-      [K, params]() -> ExprPtr { return tensor::T_pt_(K, params); },
+      [K, params]() -> ExprPtr { return tensor::tʼ(K, params); },
       [K](qnc_t& qns) { qns = combine(excitation_type_qns(K), qns); }, params);
 }
 
-ExprPtr T_pt(std::size_t K, const OpParams& params) {
+ExprPtr Tʼ(std::size_t K, const OpParams& params) {
   params.validate();
   SEQUANT_ASSERT(K > (params.skip1 ? 1 : 0));
   ExprPtr result;
   for (auto k = (params.skip1 ? 2ul : 1ul); k <= K; ++k) {
-    result += T_pt_(k, {.order = params.order,
-                        .nbatch = params.nbatch,
-                        .batch_ordinals = params.batch_ordinals,
-                        .skip1 = false});
+    result += tʼ(k, {.order = params.order,
+                     .nbatch = params.nbatch,
+                     .batch_ordinals = params.batch_ordinals,
+                     .skip1 = false});
   }
   return result;
 }
 
-ExprPtr Λ_pt_(std::size_t K, const OpParams& params) {
+ExprPtr λʼ(std::size_t K, const OpParams& params) {
   SEQUANT_ASSERT(K > 0);
   SEQUANT_ASSERT(params.order == 1 &&
                  "only first order perturbation is supported now");
 
   return make_op_from_params(
       []() -> std::wstring_view { return optype2label.at(OpType::λ_1); },
-      [K, params]() -> ExprPtr { return tensor::Λ_pt_(K, params); },
+      [K, params]() -> ExprPtr { return tensor::λʼ(K, params); },
       [K](qnc_t& qns) { qns = combine(deexcitation_type_qns(K), qns); },
       params);
 }
 
-ExprPtr Λ_pt(std::size_t K, const OpParams& params) {
+ExprPtr Λʼ(std::size_t K, const OpParams& params) {
   params.validate();
   SEQUANT_ASSERT(K > (params.skip1 ? 1 : 0));
   ExprPtr result;
   for (auto k = (params.skip1 ? 2ul : 1ul); k <= K; ++k) {
-    result += Λ_pt_(k, {.order = params.order,
-                        .nbatch = params.nbatch,
-                        .batch_ordinals = params.batch_ordinals,
-                        .skip1 = false});
+    result += λʼ(k, {.order = params.order,
+                     .nbatch = params.nbatch,
+                     .batch_ordinals = params.batch_ordinals,
+                     .skip1 = false});
   }
   return result;
 }
 
-ExprPtr R_(nann na, ncre nc, const cre<IndexSpace>& cre_space,
-           const ann<IndexSpace>& ann_space) {
+ExprPtr r(nann na, ncre nc, const cre<IndexSpace>& cre_space,
+          const ann<IndexSpace>& ann_space) {
   return ex<op_t>(
       []() -> std::wstring_view { return optype2label.at(OpType::R); },
-      [=]() -> ExprPtr { return tensor::R_(na, nc, cre_space, ann_space); },
+      [=]() -> ExprPtr { return tensor::r(na, nc, cre_space, ann_space); },
       [=](qnc_t& qns) {
         // ex -> creators in particle_space, annihilators in hole_space
         qns = combine(
@@ -1127,13 +1127,13 @@ ExprPtr R_(nann na, ncre nc, const cre<IndexSpace>& cre_space,
       });
 }
 
-ExprPtr R_(nₚ np, nₕ nh) { return R_(nann(nh), ncre(np)); }
+ExprPtr r(nₚ np, nₕ nh) { return r(nann(nh), ncre(np)); }
 
-ExprPtr L_(nann na, ncre nc, const cre<IndexSpace>& cre_space,
-           const ann<IndexSpace>& ann_space) {
+ExprPtr l(nann na, ncre nc, const cre<IndexSpace>& cre_space,
+          const ann<IndexSpace>& ann_space) {
   return ex<op_t>(
       []() -> std::wstring_view { return optype2label.at(OpType::L); },
-      [=]() -> ExprPtr { return tensor::L_(na, nc, cre_space, ann_space); },
+      [=]() -> ExprPtr { return tensor::l(na, nc, cre_space, ann_space); },
       [=](qnc_t& qns) {
         // deex -> creators in hole_space, annihilators in particle_space
         qns = combine(
@@ -1143,7 +1143,7 @@ ExprPtr L_(nann na, ncre nc, const cre<IndexSpace>& cre_space,
       });
 }
 
-ExprPtr L_(nₚ np, nₕ nh) { return L_(nann(np), ncre(nh)); }
+ExprPtr l(nₚ np, nₕ nh) { return l(nann(np), ncre(nh)); }
 
 ExprPtr R(nann na, ncre nc, const cre<IndexSpace>& cre_space,
           const ann<IndexSpace>& ann_space) {
@@ -1153,7 +1153,7 @@ ExprPtr R(nann na, ncre nc, const cre<IndexSpace>& cre_space,
   std::int64_t ra = na, rc = nc;
   while (ra >= 0 && rc >= 0) {
     if (ra == 0 && rc == 0) break;
-    result += R_(nann(ra), ncre(rc), cre_space, ann_space);
+    result += r(nann(ra), ncre(rc), cre_space, ann_space);
     if (ra == 0 || rc == 0) break;
     --ra;
     --rc;
@@ -1171,7 +1171,7 @@ ExprPtr L(nann na, ncre nc, const cre<IndexSpace>& cre_space,
   std::int64_t ra = na, rc = nc;
   while (ra >= 0 && rc >= 0) {
     if (ra == 0 && rc == 0) break;
-    result += L_(nann(ra), ncre(rc), cre_space, ann_space);
+    result += l(nann(ra), ncre(rc), cre_space, ann_space);
     if (ra == 0 || rc == 0) break;
     --ra;
     --rc;
diff --git a/SeQuant/domain/mbpt/op.hpp b/SeQuant/domain/mbpt/op.hpp
index 64f597ddc2..0c876ce147 100644
--- a/SeQuant/domain/mbpt/op.hpp
+++ b/SeQuant/domain/mbpt/op.hpp
@@ -976,7 +976,7 @@ using mbpt::nₚ;
 /// @param[in] k the rank of the particle interactions; only `k<=2` is
 /// supported
 // clang-format on
-ExprPtr H_(std::size_t k);
+ExprPtr h(std::size_t k);
 
 /// @brief Total Hamiltonian including up to `k`-body interactions
 /// @param[in] k the maximum rank of the particle interactions; only `k<=2` is
@@ -993,7 +993,7 @@ ExprPtr θ(std::size_t K);
 
 /// Makes particle-conserving excitation operator of rank \p K based on the
 /// defined context
-ExprPtr T_(std::size_t K);
+ExprPtr t(std::size_t K);
 
 /// Makes sum of particle-conserving excitation operators of all ranks up to \p
 /// K based on the defined context
@@ -1001,7 +1001,7 @@ ExprPtr T(std::size_t K, bool skip1 = false);
 
 /// Makes particle-conserving deexcitation operator of rank \p K based on the
 /// defined context
-ExprPtr Λ_(std::size_t K);
+ExprPtr λ(std::size_t K);
 
 /// Makes sum of particle-conserving deexcitation operators of all ranks up to
 /// \p K based on the defined context
@@ -1012,23 +1012,22 @@ ExprPtr Λ(std::size_t K, bool skip1 = false);
 /// @param nc number of creators
 /// @param cre_space IndexSpace on which creators act
 /// @param ann_space IndexSpace on which annihilators act
-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(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)));
 
 /// @brief Makes generic excitation operator
 /// @param np number of particle creators
 /// @param nh number of hole creators
-ExprPtr R_(nₚ np, nₕ nh);
-DEFINE_SINGLE_SIGNED_ARGUMENT_OP_VARIANT(R_);
+ExprPtr r(nₚ np, nₕ nh);
+DEFINE_SINGLE_SIGNED_ARGUMENT_OP_VARIANT(r);
 
 /// @brief Makes generic left-hand replacement operator
 /// @param na number of annihilators
 /// @param nc number of creators
 /// @param cre_space IndexSpace on which creators act
 /// @param ann_space IndexSpace on which annihilators act
-ExprPtr L_(
+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)));
@@ -1036,8 +1035,8 @@ ExprPtr L_(
 /// @brief Makes generic deexcitation operator
 /// @param np number of particle annihilators
 /// @param nh number of hole annihilators
-ExprPtr L_(nₚ np, nₕ nh);
-DEFINE_SINGLE_SIGNED_ARGUMENT_OP_VARIANT(L_);
+ExprPtr l(nₚ np, nₕ nh);
+DEFINE_SINGLE_SIGNED_ARGUMENT_OP_VARIANT(l);
 
 // clang-format off
 /// makes projector onto excited bra (if \p np > 0 && \p nh > 0) or ket (if \p np < 0 && \p nh <0) manifold
@@ -1068,14 +1067,14 @@ ExprPtr S(std::int64_t K);
 /// @param params OpParams for operator construction. Default: order=1
 /// @pre `params.order==1`, only first order perturbation is supported now
 /// @pre If batching is used, ISR must contain batching space
-ExprPtr H_pt(std::size_t R, const OpParams& params = {.order = 1});
+ExprPtr Hʼ(std::size_t R, const OpParams& params = {.order = 1});
 
 /// @brief Makes perturbed particle-conserving excitation operator
 /// @param K rank of the excitation operator
 /// @param params OpParams for operator construction. Default: order=1
 /// @pre `params.order==1`, only first order perturbation is supported now
 /// @pre If batching is used, ISR must contain batching space
-ExprPtr T_pt_(std::size_t K, const OpParams& params = {.order = 1});
+ExprPtr tʼ(std::size_t K, const OpParams& params = {.order = 1});
 
 /// @brief Makes sum of perturbed particle-conserving excitation operators
 /// @param K rank up to which the sum is to be formed
@@ -1083,15 +1082,15 @@ ExprPtr T_pt_(std::size_t K, const OpParams& params = {.order = 1});
 /// skip1=false
 /// @pre `params.order==1`, only first order perturbation is supported now
 /// @pre If batching is used, ISR must contain batching space
-ExprPtr T_pt(std::size_t K,
-             const OpParams& params = {.order = 1, .skip1 = false});
+ExprPtr Tʼ(std::size_t K,
+           const OpParams& params = {.order = 1, .skip1 = false});
 
 /// @brief Makes perturbed particle-conserving deexcitation operator
 /// @param K rank of the deexcitation operator
 /// @param params OpParams for operator construction. Default: order=1
 /// @pre `params.order==1`, only first order perturbation is supported now
 /// @pre If batching is used, ISR must contain batching space
-ExprPtr Λ_pt_(std::size_t K, const OpParams& params = {.order = 1});
+ExprPtr λʼ(std::size_t K, const OpParams& params = {.order = 1});
 
 /// @brief Makes sum of perturbed particle-conserving deexcitation operators
 /// @param K rank up to which the sum is to be formed
@@ -1099,8 +1098,8 @@ ExprPtr Λ_pt_(std::size_t K, const OpParams& params = {.order = 1});
 /// skip1=false
 /// @pre `params.order==1`, only first order perturbation is supported now
 /// @pre If batching is used, ISR must contain batching space
-ExprPtr Λ_pt(std::size_t K,
-             const OpParams& params = {.order = 1, .skip1 = false});
+ExprPtr Λʼ(std::size_t K,
+           const OpParams& params = {.order = 1, .skip1 = false});
 }  // namespace tensor
 }  // namespace op
 
@@ -1110,7 +1109,7 @@ inline namespace op {
 /// @param[in] k the rank of the particle interactions; only `k<=2` is
 /// supported
 // clang-format on
-ExprPtr H_(std::size_t k);
+ExprPtr h(std::size_t k);
 
 /// @brief Total Hamiltonian including up to `k`-body interactions
 /// @param[in] k the maximum rank of the particle interactions; only `k<=2` is
@@ -1126,14 +1125,14 @@ ExprPtr F(bool use_tensor = true, IndexSpace reference_occupied = {L"", 0});
 ExprPtr θ(std::size_t K);
 
 /// Makes particle-conserving excitation operator of rank \p K
-ExprPtr T_(std::size_t K);
+ExprPtr t(std::size_t K);
 
 /// Makes sum of particle-conserving excitation operators of all ranks up to \p
 /// K
 ExprPtr T(std::size_t K, bool skip1 = false);
 
 /// Makes particle-conserving deexcitation operator of rank \p K
-ExprPtr Λ_(std::size_t K);
+ExprPtr λ(std::size_t K);
 
 /// Makes sum of particle-conserving deexcitation operators of all ranks up to
 /// \p K
@@ -1144,23 +1143,22 @@ ExprPtr Λ(std::size_t K, bool skip1 = false);
 /// @param nc number of creators
 /// @param cre_space IndexSpace on which creators act
 /// @param ann_space IndexSpace on which annihilators act
-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(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)));
 
 /// @brief Makes generic excitation operator
 /// @param np number of particle creators
 /// @param nh number of hole creators
-ExprPtr R_(nₚ np, nₕ nh);
-DEFINE_SINGLE_SIGNED_ARGUMENT_OP_VARIANT(R_);
+ExprPtr r(nₚ np, nₕ nh);
+DEFINE_SINGLE_SIGNED_ARGUMENT_OP_VARIANT(r);
 
 /// @brief Makes generic deexcitation operator
 /// @param na number of annihilators
 /// @param nc number of creators
 /// @param cre_space IndexSpace on which creators act
 /// @param ann_space IndexSpace on which annihilators act
-ExprPtr L_(
+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)));
@@ -1168,15 +1166,15 @@ ExprPtr L_(
 /// @brief Makes generic deexcitation operator
 /// @param np number of particle annihilators
 /// @param nh number of hole annihilators
-ExprPtr L_(nₚ np, nₕ nh);
-DEFINE_SINGLE_SIGNED_ARGUMENT_OP_VARIANT(L_);
+ExprPtr l(nₚ np, nₕ nh);
+DEFINE_SINGLE_SIGNED_ARGUMENT_OP_VARIANT(l);
 
 /// @brief Makes sum of generic right-hand replacement operators up to max rank
 /// @param na number of annihilators
 /// @param nc number of creators
 /// @param cre_space IndexSpace on which creators act
 /// @param ann_space IndexSpace on which annihilators act
-/// @return `R_(na,nc) + R_(na-1,nc-1) + ...`
+/// @return `r(na,nc) + r(na-1,nc-1) + ...`
 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)));
@@ -1184,7 +1182,7 @@ ExprPtr R(nann na, ncre nc,
 /// @brief Makes sum of generic excitation operators up to max rank
 /// @param np max number of particle creators
 /// @param nh max number of hole creators
-/// @return `R_(np,nh) + R_(np-1,nh-1) + ...`
+/// @return `r(np,nh) + r(np-1,nh-1) + ...`
 ExprPtr R(nₚ np, nₕ nh);
 DEFINE_SINGLE_SIGNED_ARGUMENT_OP_VARIANT(R);
 
@@ -1193,7 +1191,7 @@ DEFINE_SINGLE_SIGNED_ARGUMENT_OP_VARIANT(R);
 /// @param nc number of creators
 /// @param cre_space IndexSpace on which creators act
 /// @param ann_space IndexSpace on which annihilators act
-/// @return `L_(na,nc) + L_(na-1,nc-1) + ...`
+/// @return `l(na,nc) + l(na-1,nc-1) + ...`
 ExprPtr L(
     nann na, ncre nc,
     const cre<IndexSpace>& cre_space = cre(get_hole_space(Spin::any)),
@@ -1202,7 +1200,7 @@ ExprPtr L(
 /// @brief Makes sum of deexcitation operators up to max rank
 /// @param np max number of particle annihilators
 /// @param nh max number of hole annihilators
-/// @return `L_(np,nh) + L_(np-1,nh-1) + ...`
+/// @return `l(np,nh) + l(np-1,nh-1) + ...`
 ExprPtr L(nₚ np, nₕ nh);
 DEFINE_SINGLE_SIGNED_ARGUMENT_OP_VARIANT(L);
 
@@ -1234,14 +1232,14 @@ ExprPtr S(std::int64_t K);
 /// @param params OpParams for operator construction. Default: order=1
 /// @pre `params.order==1`, only first order perturbation is supported now
 /// @pre If batching is used, ISR must contain batching space
-ExprPtr H_pt(std::size_t R, const OpParams& params = {.order = 1});
+ExprPtr Hʼ(std::size_t R, const OpParams& params = {.order = 1});
 
 /// @brief Makes perturbed particle-conserving excitation operator from OpParams
 /// @param K rank of the excitation operator
 /// @param params OpParams for operator construction. Default: order=1
 /// @pre `params.order==1`, only first order perturbation is supported now
 /// @pre If batching is used, ISR must contain batching space
-ExprPtr T_pt_(std::size_t K, const OpParams& params = {.order = 1});
+ExprPtr tʼ(std::size_t K, const OpParams& params = {.order = 1});
 
 /// @brief Makes sum of perturbed particle-conserving excitation operators
 /// @param K rank up to which the sum is to be formed
@@ -1249,15 +1247,15 @@ ExprPtr T_pt_(std::size_t K, const OpParams& params = {.order = 1});
 /// skip1=false
 /// @pre `params.order==1`, only first order perturbation is supported now
 /// @pre If batching is used, ISR must contain batching space
-ExprPtr T_pt(std::size_t K,
-             const OpParams& params = {.order = 1, .skip1 = false});
+ExprPtr Tʼ(std::size_t K,
+           const OpParams& params = {.order = 1, .skip1 = false});
 
 /// @brief Makes perturbed particle-conserving deexcitation operator
 /// @param K rank of the deexcitation operator
 /// @param params OpParams for operator construction. Default: order=1
 /// @pre `params.order==1`, only first order perturbation is supported now
 /// @pre If batching is used, ISR must contain batching space
-ExprPtr Λ_pt_(std::size_t K, const OpParams& params = {.order = 1});
+ExprPtr λʼ(std::size_t K, const OpParams& params = {.order = 1});
 
 /// @brief Makes sum of perturbed particle-conserving deexcitation operators
 /// @param K rank up to which the sum is to be formed
@@ -1265,8 +1263,8 @@ ExprPtr Λ_pt_(std::size_t K, const OpParams& params = {.order = 1});
 /// skip1=false
 /// @pre `params.order==1`, only first order perturbation is supported now
 /// @pre If batching is used, ISR must contain batching space
-ExprPtr Λ_pt(std::size_t K,
-             const OpParams& params = {.order = 1, .skip1 = false});
+ExprPtr Λʼ(std::size_t K,
+           const OpParams& params = {.order = 1, .skip1 = false});
 
 /// @brief computes the quantum number change effected by a given Operator or
 /// Operator Product when applied to the vacuum state
diff --git a/benchmarks/wick.cpp b/benchmarks/wick.cpp
index 5141d2695a..c5e7f61878 100644
--- a/benchmarks/wick.cpp
+++ b/benchmarks/wick.cpp
@@ -135,14 +135,13 @@ VacAvPair get_mbpt_expr(std::size_t i) {
       return {t::A(nₚ(-2)) * t::H(2) * t ::T(2) * t::T(2),
               VacAvPair::Connections{{1, 2}, {1, 3}}};
     case 2:
-      return {t::A(nₚ(-4)) * t::H_(2) * t ::T_(3) * t::T_(3),
+      return {t::A(nₚ(-4)) * t::h(2) * t ::t(3) * t::t(3),
               VacAvPair::Connections{{1, 2}, {1, 3}}};
     case 3:
-      return {t::A(-5) * t::H_(2) * t::T_(2) * t::T_(2) * t::T_(3),
+      return {t::A(-5) * t::h(2) * t::t(2) * t::t(2) * t::t(3),
               VacAvPair::Connections{{1, 2}, {1, 3}, {1, 4}}};
     case 4:
-      return {t::L_(nₚ(2), nₕ(1)) * t::H_(2) * t::R_(nₚ(1), nₕ(0)),
-              std::nullopt};
+      return {t::l(nₚ(2), nₕ(1)) * t::h(2) * t::r(nₚ(1), nₕ(0)), std::nullopt};
   }
 
   throw "Invalid index";
diff --git a/doc/examples/synopsis/synopsis6.cpp b/doc/examples/synopsis/synopsis6.cpp
index 2c1817ca5d..1e78e0926c 100644
--- a/doc/examples/synopsis/synopsis6.cpp
+++ b/doc/examples/synopsis/synopsis6.cpp
@@ -13,8 +13,8 @@ int main() {
                        .vacuum = Vacuum::SingleProduct});
 
   auto hbar =
-      H(2) + commutator(H(2), T_(2)) +
-      ex<Constant>(rational(1, 2)) * commutator(commutator(H(2), T_(2)), T_(2));
+      H(2) + commutator(H(2), t(2)) +
+      ex<Constant>(rational(1, 2)) * commutator(commutator(H(2), t(2)), t(2));
   auto ccd_eq = op::vac_av(P(nₚ(2)) * hbar);
   std::wcout << "<" << to_latex(P(nₚ(2)) * hbar) << "> = " << to_latex(ccd_eq)
              << std::endl;
diff --git a/doc/examples/user/cc.cpp b/doc/examples/user/cc.cpp
index 91fcac7ce2..21bdcb28c2 100644
--- a/doc/examples/user/cc.cpp
+++ b/doc/examples/user/cc.cpp
@@ -62,12 +62,12 @@ int main() {
 
   // start-snippet-3
   // First-order perturbed amplitude equations
-  auto t_pt = CC{2}.t_pt(1, 1);
+  auto t_pt = CC{2}.tʼ(1, 1);
   std::wcout << "T1 perturbed: " << to_latex(t_pt[1]) << "\n"
              << "T2 perturbed: " << to_latex(t_pt[2]) << "\n";
 
   // First-order perturbed Lambda amplitude equations
-  auto l_pt = CC{2}.λ_pt(1, 1);
+  auto l_pt = CC{2}.λʼ(1, 1);
   std::wcout << "λ1 perturbed: " << to_latex(l_pt[1]) << "\n"
              << "λ2 perturbed: " << to_latex(l_pt[2]) << "\n";
   // end-snippet-3
diff --git a/doc/examples/user/getting_started/ccd.cpp b/doc/examples/user/getting_started/ccd.cpp
index e83bcb350b..2924155fc7 100644
--- a/doc/examples/user/getting_started/ccd.cpp
+++ b/doc/examples/user/getting_started/ccd.cpp
@@ -13,8 +13,8 @@ int main() {
                        .vacuum = Vacuum::SingleProduct});
 
   auto hbar =
-      H(2) + commutator(H(2), T_(2)) +
-      ex<Constant>(rational(1, 2)) * commutator(commutator(H(2), T_(2)), T_(2));
+      H(2) + commutator(H(2), t(2)) +
+      ex<Constant>(rational(1, 2)) * commutator(commutator(H(2), t(2)), t(2));
   auto ccd_eq = vac_av(P(2) * hbar);
   std::wcout << "<" << to_latex(P(2) * hbar) << "> = " << to_latex(ccd_eq)
              << std::endl;
diff --git a/doc/examples/user/operator.cpp b/doc/examples/user/operator.cpp
index 4c487d7945..c1b0ab3d75 100644
--- a/doc/examples/user/operator.cpp
+++ b/doc/examples/user/operator.cpp
@@ -45,13 +45,13 @@ int main() {
   using namespace sequant::mbpt::op;
 
   // Construct a double excitation operator (T2)
-  auto T2 = T_(2);
+  auto t2 = t(2);
 
   // Construct a two-body Hamiltonian
-  auto H2 = H_(2);
+  auto h2 = h(2);
 
-  // Product: H2 * T2 * T2
-  auto expr1 = H2 * T2 * T2;
+  // Product: h2 * t2 * t2
+  auto expr1 = h2 * t2 * t2;
 
   // Screening: can an expression raise the vacuum to double excitation?
   bool raises_to_double = raises_vacuum_to_rank(expr1, 2);  // true
@@ -60,9 +60,7 @@ int main() {
   bool raises_up_to_double = raises_vacuum_up_to_rank(expr1, 2);  // true
 
   // Screening: can an expression lower a double excitation to the vacuum?
-  auto lambda2 = Λ_(2);
-  auto expr2 = lambda2 * H_(1);
-  bool lowers_to_vacuum = lowers_rank_to_vacuum(expr2, 2);  // true
+  bool lowers_to_vacuum = lowers_rank_to_vacuum(λ(2) * h(1), 2);  // true
   // end-snippet-3
   (void)raises_to_double;
   (void)raises_up_to_double;
diff --git a/doc/user/guide/cc.rst b/doc/user/guide/cc.rst
index 29a95c4b25..e621b88ba0 100644
--- a/doc/user/guide/cc.rst
+++ b/doc/user/guide/cc.rst
@@ -77,8 +77,8 @@ Coupled-Cluster Response
 
 .. code-block:: cpp
 
-   std::vector<ExprPtr> t_pt(size_t order = 1, size_t rank = 1);
-   std::vector<ExprPtr> λ_pt(size_t order = 1, size_t rank = 1);
+   std::vector<ExprPtr> tʼ(size_t rank = 1, size_t order = 1);
+   std::vector<ExprPtr> λʼ(size_t rank = 1, size_t order = 1);
 
 Derives perturbed amplitude equations for response theory calculations.
 
diff --git a/python/src/sequant/mbpt.h b/python/src/sequant/mbpt.h
index a027c75e84..a2780d2066 100644
--- a/python/src/sequant/mbpt.h
+++ b/python/src/sequant/mbpt.h
@@ -35,8 +35,7 @@ inline void __init__(py::module m) {
   py::enum_<sequant::mbpt::OpType>(m, "OpType")
       .value("h", sequant::mbpt::OpType::h)
       .value("f", sequant::mbpt::OpType::f)
-      .value("t", sequant::mbpt::OpType::t)
-      .export_values();
+      .value("t", sequant::mbpt::OpType::t);
 
   m.def("F", &sequant::mbpt::F);
   m.def("H", &sequant::mbpt::H,
@@ -45,7 +44,7 @@ inline void __init__(py::module m) {
 
   m.def(SR_OP(A));
   m.def(SR_OP(T));
-  m.def(SR_OP(T_));
+  m.def(SR_OP(t));
 
   m.def("VacuumAverage", &VacuumAverage<>);
   m.def("VacuumAverage",
diff --git a/python/test_sequant.py b/python/test_sequant.py
index 5e04827e89..b146603af6 100644
--- a/python/test_sequant.py
+++ b/python/test_sequant.py
@@ -30,8 +30,8 @@ def _test_expressions(self):
     self.assertTrue((p+s).latex)
 
   def test_ccsd(self):
-    from _sequant.mbpt import A,H,T,T_,VacuumAverage,OpType
-    ccd = VacuumAverage( A(-2) * H() * T_(2) * T_(2), [(OpType.h, OpType.t)] );  # explicit list of operator connections ..
+    from _sequant.mbpt import A,H,T,t,VacuumAverage,OpType
+    ccd = VacuumAverage( A(-2) * H() * t(2) * t(2), [(OpType.h, OpType.t)] );  # explicit list of operator connections ..
     ccsd = VacuumAverage( A(-2) * H() * T(2) * T(2));                            # .. is not needed since H and T are connected by default
     print (ccsd.latex)
 
diff --git a/tests/unit/test_mbpt.cpp b/tests/unit/test_mbpt.cpp
index 828247c06b..e95e012a84 100644
--- a/tests/unit/test_mbpt.cpp
+++ b/tests/unit/test_mbpt.cpp
@@ -171,11 +171,11 @@ TEST_CASE("mbpt", "[mbpt]") {
 
       auto f = F();
       auto t1 = T(1);
-      auto t2 = T_(2);
-      auto lambda1 = Λ_(1);
-      auto lambda2 = Λ_(2);
-      auto r_2_1 = R_(nₚ(1), nₕ(2));
-      auto r_1_2 = R_(nₚ(2), nₕ(1));
+      auto t2 = t(2);
+      auto lambda1 = λ(1);
+      auto lambda2 = λ(2);
+      auto r_2_1 = r(nₚ(1), nₕ(2));
+      auto r_1_2 = r(nₚ(2), nₕ(1));
       auto theta2 = θ(2);
       REQUIRE(to_latex(theta2) == L"{\\hat{\\theta}_{2}}");
       REQUIRE(to_latex(f) == L"{\\hat{f}}");
@@ -201,11 +201,11 @@ TEST_CASE("mbpt", "[mbpt]") {
       using qns_t [[maybe_unused]] = mbpt::qns_t;
       using namespace sequant::mbpt;
       auto f = F();
-      auto t1 = T_(1);
-      auto l1 = Λ_(1);
-      auto t2 = T_(2);
-      auto l2 = Λ_(2);
-      auto h_pt = H_pt(1, {.order = 1});
+      auto t1 = t(1);
+      auto l1 = λ(1);
+      auto t2 = t(2);
+      auto l2 = λ(2);
+      auto h_pt = Hʼ(1, {.order = 1});
       REQUIRE(to_latex(f * t1 * t2) == to_latex(canonicalize(f * t2 * t1)));
       REQUIRE(to_latex(canonicalize(f * t1 * t2)) ==
               to_latex(canonicalize(f * t2 * t1)));
@@ -260,9 +260,9 @@ TEST_CASE("mbpt", "[mbpt]") {
       using op_t = mbpt::Operator<qns_t>;
       using namespace mbpt;
       op_t f = F()->as<op_t>();
-      op_t t1 = T_(1)->as<op_t>();
-      op_t lambda2 = Λ_(2)->as<op_t>();
-      op_t r_1_2 = R_(nₚ(2), nₕ(1))->as<op_t>();
+      op_t t1 = t(1)->as<op_t>();
+      op_t lambda2 = λ(2)->as<op_t>();
+      op_t r_1_2 = r(nₚ(2), nₕ(1))->as<op_t>();
 
       REQUIRE_NOTHROW(adjoint(f));
       REQUIRE_NOTHROW(adjoint(t1));
@@ -272,7 +272,7 @@ TEST_CASE("mbpt", "[mbpt]") {
       REQUIRE(adjoint(f)() == mbpt::general_type_qns(1));
       REQUIRE(adjoint(t1)() == mbpt::deexcitation_type_qns(1));
       REQUIRE(adjoint(lambda2)() == mbpt::excitation_type_qns(2));
-      REQUIRE(adjoint(r_1_2)() == L_(nₚ(2), nₕ(1))->as<op_t>()());
+      REQUIRE(adjoint(r_1_2)() == l(nₚ(2), nₕ(1))->as<op_t>()());
 
       // adjoint(adjoint(Op)) = Op
       REQUIRE(adjoint(adjoint(t1))() == t1());
@@ -306,14 +306,14 @@ TEST_CASE("mbpt", "[mbpt]") {
 
     SECTION("screen") {
       using namespace sequant::mbpt;
-      auto g_t2_t2 = H_(2) * T_(2) * T_(2);
+      auto g_t2_t2 = h(2) * t(2) * t(2);
       REQUIRE(raises_vacuum_to_rank(g_t2_t2, 2));
       REQUIRE(raises_vacuum_up_to_rank(g_t2_t2, 2));
 
-      auto g_t2 = H_(2) * T_(2);
+      auto g_t2 = h(2) * t(2);
       REQUIRE(raises_vacuum_to_rank(g_t2, 3));
 
-      auto lambda2_f = Λ_(2) * H_(1);
+      auto lambda2_f = λ(2) * h(1);
       REQUIRE(lowers_rank_to_vacuum(lambda2_f, 2));
 
       auto expr1 = P(nₚ(0), nₕ(1)) * H() * R(nₚ(0), nₕ(1));
@@ -329,21 +329,20 @@ TEST_CASE("mbpt", "[mbpt]") {
       REQUIRE(to_latex(vev2_op) == to_latex(vev2_t));
 
       // Test screen_vac_av
-      auto hbar = mbpt::lst(H(), T_(2), 4);  // CCD Hbar
+      auto hbar = mbpt::lst(H(), t(2), 4);  // CCD Hbar
       auto screened_hbar = screen_vac_av(hbar);
-      auto expected = H_(2) * T_(2);
+      auto expected = h(2) * t(2);
       REQUIRE(simplify(screened_hbar - expected) == ex<Constant>(0));
 
-      auto expr3 = P(2) * hbar * R_(nₚ(2), nₕ(2));
+      auto expr3 = P(2) * hbar * r(nₚ(2), nₕ(2));
       auto screened_expr3 = screen_vac_av(expr3);
-      auto expected3 =
-          op::P(2) * (H_(2) * T_(2) + H_(2) + H_(1)) * R_(nₚ(2), nₕ(2));
+      auto expected3 = op::P(2) * (h(2) * t(2) + h(2) + h(1)) * r(nₚ(2), nₕ(2));
       REQUIRE(simplify(screened_expr3 - expected3) == ex<Constant>(0));
 
       auto expr4 = P(nₚ(2), nₕ(1)) * hbar * R(nₚ(1), nₕ(0));
       auto screened_expr4 = screen_vac_av(expr4);
       auto expected4 = P(nₚ(2), nₕ(1)) *
-                       (H_(1) + H_(1) * T_(2) + H_(2) * T_(2) + H_(2)) *
+                       (h(1) + h(1) * t(2) + h(2) * t(2) + h(2)) *
                        R(nₚ(1), nₕ(0));
       REQUIRE(simplify(screened_expr4 - expected4) == ex<Constant>(0));
 
@@ -358,41 +357,40 @@ TEST_CASE("mbpt", "[mbpt]") {
 
       // non-unitary, rank 3
       auto expr1 =
-          lst(H(), T_(2), 3, {.unitary = false, .use_commutators = false});
+          lst(H(), t(2), 3, {.unitary = false, .use_commutators = false});
       auto expected1 =
-          H() * (ex<Constant>(1) + T_(2) +
-                 ex<Constant>(rational{1, 2}) * T_(2) * T_(2) +
-                 ex<Constant>(rational{1, 6}) * T_(2) * T_(2) * T_(2));
+          H() *
+          (ex<Constant>(1) + t(2) + ex<Constant>(rational{1, 2}) * t(2) * t(2) +
+           ex<Constant>(rational{1, 6}) * t(2) * t(2) * t(2));
       REQUIRE(simplify(expr1 - expected1) == ex<Constant>(0));
 
       auto expr2 =
-          lst(H(), T_(2), 3, {.unitary = false, .use_commutators = true});
+          lst(H(), t(2), 3, {.unitary = false, .use_commutators = true});
       auto expected2 =
-          H() + commutator(H(), T_(2)) +
+          H() + commutator(H(), t(2)) +
           ex<Constant>(rational{1, 2}) *
-              commutator(commutator(H(), T_(2)), T_(2)) +
+              commutator(commutator(H(), t(2)), t(2)) +
           ex<Constant>(rational{1, 6}) *
-              commutator(commutator(commutator(H(), T_(2)), T_(2)), T_(2));
+              commutator(commutator(commutator(H(), t(2)), t(2)), t(2));
       REQUIRE(simplify(expr2 - expected2) == ex<Constant>(0));
 
       // unitary, rank 2
       using sequant::adjoint;
       auto expr3 =
-          lst(H(), T_(2), 2, {.unitary = true, .use_commutators = false});
+          lst(H(), t(2), 2, {.unitary = true, .use_commutators = false});
       auto expected3 =
-          H() + H() * T_(2) + adjoint(T_(2)) * H() +
-          adjoint(T_(2)) * H() * T_(2) +
-          H() * ex<Constant>(rational{1, 2}) * T_(2) * T_(2) +
-          ex<Constant>(rational{1, 2}) * adjoint(T_(2)) * adjoint(T_(2)) * H();
+          H() + H() * t(2) + adjoint(t(2)) * H() + adjoint(t(2)) * H() * t(2) +
+          H() * ex<Constant>(rational{1, 2}) * t(2) * t(2) +
+          ex<Constant>(rational{1, 2}) * adjoint(t(2)) * adjoint(t(2)) * H();
       REQUIRE(simplify(expr3 - expected3) == ex<Constant>(0));
 
       auto expr4 =
-          lst(H(), T_(2), 2, {.unitary = true, .use_commutators = true});
-      auto generator = commutator(H(), T_(2)) - commutator(H(), adjoint(T_(2)));
-      auto expected4 = H() + generator +
-                       ex<Constant>(rational{1, 2}) *
-                           (commutator(generator, T_(2)) -
-                            commutator(generator, adjoint(T_(2))));
+          lst(H(), t(2), 2, {.unitary = true, .use_commutators = true});
+      auto generator = commutator(H(), t(2)) - commutator(H(), adjoint(t(2)));
+      auto expected4 =
+          H() + generator +
+          ex<Constant>(rational{1, 2}) * (commutator(generator, t(2)) -
+                                          commutator(generator, adjoint(t(2))));
       REQUIRE(simplify(expr4 - expected4) == ex<Constant>(0));
     }  // SECTION("lst")
 
@@ -408,7 +406,7 @@ TEST_CASE("mbpt", "[mbpt]") {
       REQUIRE(to_latex(simplify(theta1.tensor_form())) ==
               L"{{\\theta^{{p_2}}_{{p_1}}}{\\tilde{a}^{{p_1}}_{{p_2}}}}");
 
-      auto R_2 = R_(2)->as<op_t>();
+      auto R_2 = r(2)->as<op_t>();
       //    std::wcout << "R_2: " << to_latex(simplify(R_2.tensor_form())) <<
       //    std::endl;
       REQUIRE(
@@ -417,7 +415,7 @@ TEST_CASE("mbpt", "[mbpt]") {
           L"{{a_1}{"
           L"a_2}}_{{i_1}{i_2}}}}");
 
-      auto L_3 = L_(3)->as<op_t>();
+      auto L_3 = l(3)->as<op_t>();
       //    std::wcout << "L_3: " << to_latex(simplify(L_3.tensor_form())) <<
       //    std::endl;
       REQUIRE(
@@ -425,7 +423,7 @@ TEST_CASE("mbpt", "[mbpt]") {
           L"{{{\\frac{1}{36}}}{\\bar{L}^{{a_1}{a_2}{a_3}}_{{i_1}{i_2}{i_3}}}{"
           L"\\tilde{a}^{{i_1}{i_2}{i_3}}_{{a_1}{a_2}{a_3}}}}");
 
-      auto R_2_3 = R_(nₚ(3), nₕ(2))->as<op_t>();
+      auto R_2_3 = r(nₚ(3), nₕ(2))->as<op_t>();
       //    std::wcout << "R_2_3: " << to_latex(simplify(R_2_3.tensor_form()))
       //    << std::endl;
       REQUIRE(to_latex(simplify(R_2_3.tensor_form())) ==
@@ -433,8 +431,8 @@ TEST_CASE("mbpt", "[mbpt]") {
               L"\\tilde{a}^{"
               L"{a_1}{a_2}{a_3}}_{\\textvisiblespace\\,{i_1}{i_2}}}}");
 
-      auto L_1_2 = L_(nₚ(1), nₕ(2))->as<op_t>();
-      // std::wcout << "L_(1,2): " << to_latex(simplify(L_1_2.tensor_form())) <<
+      auto L_1_2 = l(nₚ(1), nₕ(2))->as<op_t>();
+      // std::wcout << "l(1,2): " << to_latex(simplify(L_1_2.tensor_form())) <<
       // std::endl;
       REQUIRE(
           to_latex(simplify(L_1_2.tensor_form())) ==
@@ -530,32 +528,31 @@ TEST_CASE("mbpt", "[mbpt]") {
           get_default_context().index_space_registry()->retrieve(L"z"));
 
       using namespace mbpt;
-      REQUIRE_NOTHROW(op::H_pt(1, {.nbatch = 1}));
-      REQUIRE_NOTHROW(op::H_pt(2, {.nbatch = 2}));
-      REQUIRE_NOTHROW(
-          op::Λ_pt(3, {.batch_ordinals = {3, 4, 5}, .skip1 = true}));
-      REQUIRE_NOTHROW(op::T_pt(1, {.nbatch = 20}));
+      REQUIRE_NOTHROW(op::Hʼ(1, {.nbatch = 1}));
+      REQUIRE_NOTHROW(op::Hʼ(2, {.nbatch = 2}));
+      REQUIRE_NOTHROW(op::Λʼ(3, {.batch_ordinals = {3, 4, 5}, .skip1 = true}));
+      REQUIRE_NOTHROW(op::Tʼ(1, {.nbatch = 20}));
 
       // invalid usages
 #if SEQUANT_ASSERT_BEHAVIOR == SEQUANT_ASSERT_THROW
       // cannot set both nbatch and batch_ordinals
-      REQUIRE_THROWS_AS(op::H_pt(2, {.nbatch = 2, .batch_ordinals = {1, 2}}),
+      REQUIRE_THROWS_AS(op::Hʼ(2, {.nbatch = 2, .batch_ordinals = {1, 2}}),
                         sequant::Exception);
       // all ordinals must be unique
-      REQUIRE_THROWS_AS(op::H_pt(2, {.batch_ordinals = {1, 2, 2}}),
+      REQUIRE_THROWS_AS(op::Hʼ(2, {.batch_ordinals = {1, 2, 2}}),
                         sequant::Exception);
       // ordinals must be sorted
-      REQUIRE_THROWS_AS(op::H_pt(1, {.batch_ordinals = {3, 2}}),
+      REQUIRE_THROWS_AS(op::Hʼ(1, {.batch_ordinals = {3, 2}}),
                         sequant::Exception);
 #endif
 
       // operations
-      auto h0 = op::H_pt(1);
+      auto h0 = op::Hʼ(1);
       REQUIRE(to_latex(h0) == L"{\\hat{h¹}}");
 
-      auto h1 = op::H_pt(1, {.nbatch = 1});
-      auto h1_2 = op::H_pt(1, {.batch_ordinals = {1, 2}});
-      auto pt1 = op::T_pt(2, {.batch_ordinals = {1}});
+      auto h1 = op::Hʼ(1, {.nbatch = 1});
+      auto h1_2 = op::Hʼ(1, {.batch_ordinals = {1, 2}});
+      auto pt1 = op::Tʼ(2, {.batch_ordinals = {1}});
 
       auto sum0 = h0 + h1;
       simplify(sum0);
@@ -646,8 +643,8 @@ TEST_CASE("mbpt", "[mbpt]") {
 
   // H2**T3**T3 -> R4
   SECTION("wick(H2**T3**T3 -> R4)") {
-    auto result = t::vac_av(t::A(nₚ(-4)) * t::H_(2) * t::T_(3) * t::T_(3),
-                            {{1, 2}, {1, 3}});
+    auto result =
+        t::vac_av(t::A(nₚ(-4)) * t::h(2) * t::t(3) * t::t(3), {{1, 2}, {1, 3}});
 
     // std::wcout << "H2**T3**T3 -> R4 = " << to_latex_align(result, 20)
     //            << std::endl;
@@ -660,9 +657,8 @@ TEST_CASE("mbpt", "[mbpt]") {
   {
     ExprPtr ref_result;
     SECTION("wick(H2**T2**T2**T3 -> R5)") {
-      ref_result =
-          t::vac_av(t::A(-5) * t::H_(2) * t::T_(2) * t::T_(2) * t::T_(3),
-                    {{1, 2}, {1, 3}, {1, 4}});
+      ref_result = t::vac_av(t::A(-5) * t::h(2) * t::t(2) * t::t(2) * t::t(3),
+                             {{1, 2}, {1, 3}, {1, 4}});
       REQUIRE(ref_result->size() == 7);
     }
   }
@@ -672,7 +668,7 @@ TEST_CASE("mbpt", "[mbpt]") {
 SECTION("SRSO Fock") {
   // <2p1h|H2|1p> ->
   SECTION("wick(<2p1h|H2|1p>)") {
-    auto input = t::L_(nₚ(2), nₕ(1)) * t::H_(2) * t::R_(nₚ(1), nₕ(0));
+    auto input = t::l(nₚ(2), nₕ(1)) * t::h(2) * t::r(nₚ(1), nₕ(0));
     auto result = t::vac_av(input);
 
     REQUIRE(result->is<Product>());  // product ...
@@ -681,7 +677,7 @@ SECTION("SRSO Fock") {
 
   // <2p1h|H2|2p1h(c)> ->
   SECTION("wick(<2p1h|H2|2p1h(c)>)") {
-    auto input = t::L_(nₚ(2), nₕ(1)) * t::H() * t::R_(nₚ(2), nₕ(1));
+    auto input = t::l(nₚ(2), nₕ(1)) * t::H() * t::r(nₚ(2), nₕ(1));
     auto result = t::vac_av(input);
 
     // std::wcout << "<2p1h|H|2p1h(c)> = " << to_latex(result)
@@ -698,8 +694,8 @@ SECTION("SRSO-PNO") {
 
   // H2**T2**T2 -> R2
   SECTION("wick(H2**T2**T2 -> R2)") {
-    auto result = t::vac_av(t::A(nₚ(-2)) * t::H_(2) * t::T_(2) * t::T_(2),
-                            {{1, 2}, {1, 3}});
+    auto result =
+        t::vac_av(t::A(nₚ(-2)) * t::h(2) * t::t(2) * t::t(2), {{1, 2}, {1, 3}});
 
     REQUIRE(result->size() == 4);
   }
@@ -712,7 +708,7 @@ SECTION("SRSF") {
 
   // H2 -> R2
   SECTION("wick(H2 -> R2)") {
-    auto result = t::vac_av(t::S(-2) * t::H_(2));
+    auto result = t::vac_av(t::S(-2) * t::h(2));
 
     {
       // std::wcout << "H2 -> R2 = " << to_latex_align(result, 0, 1)
@@ -724,7 +720,7 @@ SECTION("SRSF") {
 
   // H2**T2 -> R2
   SECTION("wick(H2**T2 -> R2)") {
-    auto result = t::vac_av(t::S(-2) * t::H_(2) * t::T_(2), {{1, 2}});
+    auto result = t::vac_av(t::S(-2) * t::h(2) * t::t(2), {{1, 2}});
 
     {
       // std::wcout << "H2**T2 -> R2 = " << to_latex_align(result, 0, 1)
@@ -742,9 +738,9 @@ SECTION("MRSO") {
 
   SECTION("wick(H2**T2 -> 0)") {
     {
-      auto result = t::ref_av(t::H_(2) * t::T_(2), {{0, 1}});
+      auto result = t::ref_av(t::h(2) * t::t(2), {{0, 1}});
 
-      auto result_wo_top = t::ref_av(t::H_(2) * t::T_(2), {{0, 1}},
+      auto result_wo_top = t::ref_av(t::h(2) * t::t(2), {{0, 1}},
                                      /* use_topology = */ false);
       REQUIRE(simplify(result - result_wo_top) == ex<Constant>(0));
     }
@@ -755,17 +751,17 @@ SECTION("MRSO") {
       ctx.set(mbpt::make_mr_spaces());
       ctx.set(Vacuum::Physical);
       auto ctx_resetter = set_scoped_default_context(ctx);
-      auto result_phys = t::ref_av(t::H_(2) * t::T_(2), {{0, 1}});
+      auto result_phys = t::ref_av(t::h(2) * t::t(2), {{0, 1}});
     }
   }
 
   // H2 ** T2 ** T2 -> 0
   SECTION("wick(H2**T2**T2 -> 0)") {
     // first without use of topology
-    auto result = t::ref_av(t::H_(2) * t::T_(2) * t::T_(2), {{0, 1}},
+    auto result = t::ref_av(t::h(2) * t::t(2) * t::t(2), {{0, 1}},
                             /* use_topology = */ false);
     // now with topology use
-    auto result_top = t::ref_av(t::H_(2) * t::T_(2) * t ::T_(2), {{0, 1}},
+    auto result_top = t::ref_av(t::h(2) * t::t(2) * t ::t(2), {{0, 1}},
                                 /* use_topology = */ true);
 
     REQUIRE(simplify(result - result_top) == ex<Constant>(0));
@@ -774,7 +770,7 @@ SECTION("MRSO") {
 #if 0
     // H**T12 -> R2
     SECTION("wick(H**T2 -> R2)") {
-      auto result = t::ref_av(t::A(-2) * t::H() * t::T_(2), {{1, 2}});
+      auto result = t::ref_av(t::A(-2) * t::H() * t::t(2), {{1, 2}});
 
       {
         std::wcout << "H*T2 -> R2 = " << to_latex_align(result, 0, 1)
@@ -792,11 +788,11 @@ SECTION("MRSF") {
   auto ctx_resetter = set_scoped_default_context(ctx);
 
   SECTION("wick(H2**T2 -> 0)") {
-    auto result = t::ref_av(t::H_(2) * t::T_(2), {{0, 1}});
+    auto result = t::ref_av(t::h(2) * t::t(2), {{0, 1}});
 
     {
       // make sure get same result without use of topology
-      auto result_wo_top = t::ref_av(t::H_(2) * t::T_(2), {{0, 1}},
+      auto result_wo_top = t::ref_av(t::h(2) * t::t(2), {{0, 1}},
                                      /* use_topology = */ false);
 
       REQUIRE(simplify(result - result_wo_top) == ex<Constant>(0));
@@ -804,7 +800,7 @@ SECTION("MRSF") {
 
     {
       // make sure get same result using operators
-      auto result_op = o::ref_av(o::H_(2) * o::T_(2));
+      auto result_op = o::ref_av(o::h(2) * o::t(2));
 
       REQUIRE(result_op->size() == result->size());
       REQUIRE(simplify(result - result_op) == ex<Constant>(0));
@@ -879,4 +875,115 @@ SECTION("rules") {
     }
   }
 }  // SECTION("rules")
+
+SECTION("manuscript-examples") {
+  using namespace sequant::mbpt;
+
+  /// When order == 0, returns sim. transformed of zeroth order Hamiltonian.
+  /// When order > 0, returns sim. transformed perturbation operator or given
+  /// order.
+  auto H̅ = [](size_t order = 0) {
+    auto hbar0 = lst(op::H(), T(2), 4);
+    if (order == 0) return hbar0;
+    // only one-body perturbation operator
+    auto hbar_pt = lst(op::Hʼ(/*rank*/ 1, {.order = order}), T(2), 2);
+    return hbar_pt;
+  };
+
+  SECTION("CCD Term") {
+    auto expr = ref_av(P(2) * H() * t(2) * t(2));
+    REQUIRE(expr.size() == 4);
+  }
+
+  SECTION("Ground State Amplitudes") {
+    // connectivity info for t and λ amplitude equations
+    const auto t_connect = default_op_connections();
+    const auto l_connect =
+        concat(default_op_connections(),
+               OpConnections<OpType>{{OpType::h, OpType::A},
+                                     {OpType::f, OpType::A},
+                                     {OpType::g, OpType::A}});
+
+    auto t = ref_av(P(2) * H̅(), t_connect);
+    auto λ = ref_av((1 + Λ(2)) * H̅() * P(-2), l_connect);
+
+    // numer of terms are verified against srcc results
+    REQUIRE(t.size() == 31);
+    REQUIRE(λ.size() == 32);
+  }
+
+  SECTION("CC LR Function") {
+    const int N = 2;  // CC rank
+
+    auto θ̅ = lst(θ(1), T(N), 2);
+    auto expr = (1 + Λ(N)) * θ̅ * Tʼ(N) + Λʼ(N) * θ̅;
+    auto result = ref_av(expr, {{L"θ", L"t"}, {L"θ", L"t¹"}});
+    // number of terms is verified against MPQC4 implementation
+    REQUIRE(result.size() == 21);
+  }
+
+  SECTION("EOM-CC Equations") {
+    // connectivity info for right and left amplitude equations
+    const auto r_connect =
+        concat(default_op_connections(),
+               OpConnections<OpType>{{OpType::h, OpType::R},
+                                     {OpType::f, OpType::R},
+                                     {OpType::g, OpType::R}});
+    const auto l_connect =
+        concat(default_op_connections(),
+               OpConnections<OpType>{{OpType::h, OpType::A},
+                                     {OpType::f, OpType::A},
+                                     {OpType::g, OpType::A}});
+
+    // EE
+    auto r_EE = ref_av(P(2) * H̅() * R(2), r_connect);
+    auto l_EE = ref_av(L(2) * H̅() * P(-2), l_connect);
+    // EA
+    auto r_EA = ref_av(P(nₚ(2), nₕ(1)) * H̅() * R(nₚ(2), nₕ(1)), r_connect);
+    auto l_EA = ref_av(L(nₚ(2), nₕ(1)) * H̅() * P(nₚ(-2), nₕ(-1)), l_connect);
+    // IP
+    auto r_IP = ref_av(P(nₚ(1), nₕ(2)) * H̅() * R(nₚ(1), nₕ(2)), r_connect);
+    auto l_IP = ref_av(L(nₚ(1), nₕ(2)) * H̅() * P(nₚ(-1), nₕ(-2)), l_connect);
+
+    // number of terms are verified against eomcc results
+    REQUIRE(r_EE.size() == 53);
+    REQUIRE(l_EE.size() == 31);
+    REQUIRE(r_EA.size() == 32);
+    REQUIRE(l_EA.size() == 24);
+    REQUIRE(r_IP.size() == 32);
+    REQUIRE(l_IP.size() == 24);
+  }
+
+  SECTION("CC Perturbed Amplitudes") {
+    // connectivity info for perturbed t and λ amplitude equations
+    const auto t_connect =
+        concat(default_op_connections(),
+               OpConnections<OpType>{{OpType::h, OpType::t_1},
+                                     {OpType::f, OpType::t_1},
+                                     {OpType::g, OpType::t_1},
+                                     {OpType::h_1, OpType::t}});
+
+    const auto l_connect =
+        concat(default_op_connections(),
+               OpConnections<OpType>{{OpType::h, OpType::t_1},
+                                     {OpType::f, OpType::t_1},
+                                     {OpType::g, OpType::t_1},
+                                     {OpType::h_1, OpType::t},
+                                     {OpType::h, OpType::A},
+                                     {OpType::f, OpType::A},
+                                     {OpType::g, OpType::A},
+                                     {OpType::h_1, OpType::A}});
+
+    // perturbed t amplitudes (Eq 18 in SQ Manuscript #2)
+    auto t = ref_av(P(2) * (H̅(1) + H̅() * Tʼ(2) - "ω" * tʼ(2)), t_connect);
+    // perturbed λ amplitudes (Eq 19 in SQ Manuscript #2)
+    auto λ = ref_av(
+        ((1 + Λ(2)) * (H̅(1) + H̅() * Tʼ(2)) + Λʼ(2) * H̅() + "ω" * λʼ(2)) * P(-2),
+        l_connect);
+
+    // number of terms are verified against MPQC4 implementation
+    REQUIRE(t.size() == 58);
+    REQUIRE(λ.size() == 63);
+  }
+}  // SECTION("manuscript-examples")
 }
diff --git a/tests/unit/test_tensor_network.cpp b/tests/unit/test_tensor_network.cpp
index 3450543ae0..4e092eb34b 100644
--- a/tests/unit/test_tensor_network.cpp
+++ b/tests/unit/test_tensor_network.cpp
@@ -389,7 +389,7 @@ TEST_CASE("tensor_network", "[elements]") {
     }
 
     {  // with Tensors and NormalOperators
-      auto tmp = t::A(nₚ(-2)) * t::H_(2) * t::T_(2) * t::T_(2);
+      auto tmp = t::A(nₚ(-2)) * t::h(2) * t::t(2) * t::t(2);
       REQUIRE_NOTHROW(TN(tmp->as<Product>().factors()));
     }
 
@@ -516,7 +516,7 @@ TEST_CASE("tensor_network", "[elements]") {
     // can't use operator expression (due to unspecified order of evaluation of
     // function arguments), must use initializer list
     auto tmp = ex<Product, std::initializer_list<ExprPtr>>(
-        {t::A(nₚ(-2)), t::H_(2), t::T_(2), t::T_(2), t::T_(2)});
+        {t::A(nₚ(-2)), t::h(2), t::t(2), t::t(2), t::t(2)});
     // canonicalize to avoid dependence on the implementation details of
     // mbpt::sr::make_op
     // N.B. graph structure before canonicalization will depend on the input
@@ -930,7 +930,7 @@ TEST_CASE("tensor_network_v2", "[elements]") {
     }
 
     {  // with Tensors and NormalOperators
-      auto tmp = t::A(nₚ(-2)) * t::H_(2) * t::T_(2) * t::T_(2);
+      auto tmp = t::A(nₚ(-2)) * t::h(2) * t::t(2) * t::t(2);
       REQUIRE_NOTHROW(TensorNetworkV2(tmp->as<Product>().factors()));
     }
 
@@ -1572,7 +1572,7 @@ TEST_CASE("tensor_network_v3", "[elements]") {
     }
 
     {  // with Tensors and NormalOperators
-      auto tmp = t::A(nₚ(-2)) * t::H_(2) * t::T_(2) * t::T_(2);
+      auto tmp = t::A(nₚ(-2)) * t::h(2) * t::t(2) * t::t(2);
       REQUIRE_NOTHROW(TN(tmp->as<Product>().factors()));
     }