two_center_table.h

Go to the documentation of this file.

#ifndef TWO_CENTER_TABLE_H
#define TWO_CENTER_TABLE_H
 
#include <ATen/tensor.h>
#include "source_base/spherical_bessel_transformer.h"
#include "source_basis/module_nao/radial_collection.h"
 
class TwoCenterTable
{
  public:
    TwoCenterTable() = default;
    ~TwoCenterTable() { delete[] rgrid_; }
 
    TwoCenterTable(const TwoCenterTable&) = delete;
    TwoCenterTable& operator=(const TwoCenterTable&) = delete;
 
    void build(const RadialCollection& bra,  
               const RadialCollection& ket,  
               const char op,                
               const int nr,                 
               const double cutoff           
    );
 
    char op() const { return op_; }
 
    int nr() const { return nr_; }
 
    // returns the number of table entries
    int ntab() const { return ntab_; }
 
    double rmax() const { return rmax_; }
 
    const double* table(const int itype1,        
                        const int l1,            
                        const int izeta1,        
                        const int itype2,        
                        const int l2,            
                        const int izeta2,        
                        const int l,             
                        const bool deriv = false 
    ) const;
 
    void lookup(const int itype1,      
                const int l1,          
                const int izeta1,      
                const int itype2,      
                const int l2,          
                const int izeta2,      
                const int l,           
                const double R,        
                double* val,           
                double* dval = nullptr 
    ) const;
 
    // NOTE: "nchi_ket" and "lmax_ket" are added for the purpose of "snap" in TwoCenterIntegrator
    // which generates a batch of two-center integrals depending on those information.
    // This might not be the intended purpose of this class.
 
    int nchi_ket(const int itype, const int l) const
    {
        assert(itype >= 0 && itype < nchi_ket_.shape().dim_size(0));
        assert(l >= 0 && l < nchi_ket_.shape().dim_size(1));
        return nchi_ket_.get_value<int>(itype, l);
    }
 
    int lmax_ket() const { return nchi_ket_.shape().dim_size(1) - 1; }
 
    size_t memory() const {
        return (table_.NumElements() + dtable_.NumElements()
                + nchi_ket_.NumElements() + index_map_.NumElements() + nr_) * sizeof(double);
    }
 
  private:
    char op_ = '\0';   
    int ntab_ = 0;     
    int nr_ = 0;       
    double rmax_= 0.0; 
    double* rgrid_ = nullptr;
 
    container::Tensor nchi_ket_{container::DataType::DT_INT, container::TensorShape({0})};
 
    container::Tensor table_{container::DataType::DT_DOUBLE, container::TensorShape({0})};
 
    container::Tensor dtable_{container::DataType::DT_DOUBLE, container::TensorShape({0})};
 
    container::Tensor index_map_{container::DataType::DT_INT, container::TensorShape({0})};
 
    int& table_index(const NumericalRadial* ptr_rad1, const NumericalRadial* ptr_rad2, const int l);
 
    void cleanup();
 
    bool is_present(const int itype1,
                    const int l1,
                    const int izeta1,
                    const int itype2,
                    const int l2,
                    const int izeta2,
                    const int l) const;
 
    double dfact(int l) const;
 
    typedef void(TwoCenterTable::*looped_func)(const NumericalRadial*, const NumericalRadial*, const int l);
 
    void two_center_loop(const RadialCollection& bra, 
                         const RadialCollection& ket, 
                         looped_func f);
 
    void _indexing(const NumericalRadial* it1, const NumericalRadial* it2, const int l);
    void _tabulate(const NumericalRadial* it1, const NumericalRadial* it2, const int l);
};
 
#endif