.. _program_listing_file_SeQuant_domain_mbpt_rules_df.cpp:

Program Listing for File df.cpp
===============================

|exhale_lsh| :ref:`Return to documentation for file <file_SeQuant_domain_mbpt_rules_df.cpp>` (``SeQuant/domain/mbpt/rules/df.cpp``)

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

.. code-block:: cpp

   //
   // Created by Eduard Valeyev on 3/8/25.
   //
   
   #include <SeQuant/domain/mbpt/rules/df.hpp>
   
   #include <SeQuant/core/expr.hpp>
   #include <SeQuant/core/space.hpp>
   #include <SeQuant/core/utility/macros.hpp>
   
   #include <range/v3/view.hpp>
   
   #include <string_view>
   
   namespace sequant::mbpt {
   
   ExprPtr density_fit_impl(Tensor const& tnsr, Index const& aux_idx,
                            std::wstring_view factor_label) {
     SEQUANT_ASSERT(tnsr.bra_rank() == 2     //
                    && tnsr.ket_rank() == 2  //
                    && tnsr.aux_rank() == 0);
   
     auto t1 = ex<Tensor>(factor_label, bra({ranges::front(tnsr.bra())}),
                          ket({ranges::front(tnsr.ket())}), aux({aux_idx}));
   
     auto t2 = ex<Tensor>(factor_label, bra({ranges::back(tnsr.bra())}),
                          ket({ranges::back(tnsr.ket())}), aux({aux_idx}));
   
     if (tnsr.symmetry() == Symmetry::Antisymm) {
       auto t3 = ex<Tensor>(factor_label, bra({ranges::back(tnsr.bra())}),
                            ket({ranges::front(tnsr.ket())}), aux({aux_idx}));
   
       auto t4 = ex<Tensor>(factor_label, bra({ranges::front(tnsr.bra())}),
                            ket({ranges::back(tnsr.ket())}), aux({aux_idx}));
       return t1 * t2 - t3 * t4;
     }
   
     return t1 * t2;
   }
   
   ExprPtr density_fit(ExprPtr const& expr, IndexSpace aux_space,
                       std::wstring_view tensor_label,
                       std::wstring_view factor_label) {
     using ranges::views::transform;
     if (expr->is<Sum>())
       return ex<Sum>(*expr | transform([&](auto&& x) {
         return density_fit(x, aux_space, tensor_label, factor_label);
       }));
   
     else if (expr->is<Tensor>()) {
       auto const& tensor = expr->as<Tensor>();
       if (tensor.label() == tensor_label  //
           && tensor.bra_rank() == 2       //
           && tensor.ket_rank() == 2)
         return density_fit_impl(tensor, Index(aux_space, 1), factor_label);
       else
         return expr;
     } else if (expr->is<Product>()) {
       auto const& prod = expr->as<Product>();
   
       Product result;
       result.scale(prod.scalar());
       size_t aux_ix = 0;
       for (auto&& f : prod.factors())
         if (f.is<Tensor>() && f.as<Tensor>().label() == tensor_label) {
           auto const& g = f->as<Tensor>();
           auto g_df =
               density_fit_impl(g, Index(aux_space, ++aux_ix), factor_label);
           result.append(1, std::move(g_df), Product::Flatten::Yes);
         } else {
           result.append(1, f, Product::Flatten::No);
         }
       return ex<Product>(std::move(result));
     } else
       return expr;
   }
   
   }  // namespace sequant::mbpt