density_matrix.h

Go to the documentation of this file.

#ifndef DENSITY_MATRIX_H
#define DENSITY_MATRIX_H
 
#include <string>
 
#include "source_cell/module_neighbor/sltk_grid_driver.h"
#include "source_lcao/record_adj.h"
#include "source_lcao/module_hcontainer/hcontainer.h"
 
namespace elecstate
{
template<typename T> struct ShiftRealComplex
{
    using type = void;
};
 
template<>
struct ShiftRealComplex<double> 
{
    using type = std::complex<double>;
};
 
template<>
struct ShiftRealComplex<std::complex<double>> 
{
    using type = double;
};
 
template <typename TK, typename TR>
class DensityMatrix
{
    using TRShift = typename ShiftRealComplex<TR>::type;
 
    public:
    ~DensityMatrix();
 
    DensityMatrix(const Parallel_Orbitals* _paraV, const int nspin, const std::vector<ModuleBase::Vector3<double>>& kvec_d, const int nk);
 
    DensityMatrix(const Parallel_Orbitals* _paraV, const int nspin);
 
    void init_DMR(const Grid_Driver* GridD_in, const UnitCell* ucell);
 
    void init_DMR(Record_adj& ra, const UnitCell* ucell);
 
    void init_DMR(const hamilt::HContainer<TR>& _DMR_in);
 
    // the old input type ``:HContainer<complex<double>` causes redefination error if TR = complex<double>
    void init_DMR(const hamilt::HContainer<TRShift>& _DMR_in);
 
    void set_DMK(const int ispin, const int ik, const int i, const int j, const TK value);
 
    void set_DMK_zero();
    
    TK get_DMK(const int ispin, const int ik, const int i, const int j) const;
 
    int get_DMK_nks() const;
    int get_DMK_size() const;
 
    int get_DMK_nrow() const;
 
    int get_DMK_ncol() const;
 
    hamilt::HContainer<TR>* get_DMR_pointer(const int ispin) const;
 
    const std::vector<hamilt::HContainer<TR>*>& get_DMR_vector() const {return this->_DMR;}
    std::vector<hamilt::HContainer<TR>*>& get_DMR_vector() {return this->_DMR;}
 
    const std::vector<std::vector<TR>>& get_DMR_save() const {return this->_DMR_save;}
    std::vector<std::vector<TR>>& get_DMR_save() {return this->_DMR_save;}
 
    TK* get_DMK_pointer(const int ik) const;
 
    const std::vector<std::vector<TK>>& get_DMK_vector() const {return this->_DMK;}
    std::vector<std::vector<TK>>& get_DMK_vector() {return this->_DMK;}
 
    void set_DMK_pointer(const int ik, TK* DMK_in);
 
    const Parallel_Orbitals* get_paraV_pointer() const {return this->_paraV;}
 
    const std::vector<ModuleBase::Vector3<double>>& get_kvec_d() const { return this->_kvec_d; }
 
    void cal_DMR(const int ik_in = -1);
 
    void cal_DMR_td(const UnitCell& ucell, const ModuleBase::Vector3<double> At, const int ik_in = -1);
 
    void cal_DMR_full(hamilt::HContainer<std::complex<double>>* dmR_out) const;
 
    void switch_dmr(const int mode);
 
    void write_DMK(const std::string directory, const int ispin, const int ik);
 
    void read_DMK(const std::string directory, const int ispin, const int ik);
 
    void save_DMR();
    
    std::vector<ModuleBase::ComplexMatrix> EDMK; // for TD-DFT
 
#ifdef __PEXSI
    std::vector<TK*> pexsi_EDM;
#endif
 
  private:
    std::vector<hamilt::HContainer<TR>*> _DMR;
    std::vector<std::vector<TR>> _DMR_save;
 
    std::vector<hamilt::HContainer<TR>*> _DMR_grid;
 
    // std::vector<ModuleBase::ComplexMatrix> _DMK;
    std::vector<std::vector<TK>> _DMK;
 
    const std::vector<ModuleBase::Vector3<double>> _kvec_d;
 
    const Parallel_Orbitals* _paraV = nullptr;
 
    int _nspin = 1;
 
    int _nk = 0;
 
    std::vector<TR> dmr_origin_;
    TR* dmr_tmp_ = nullptr;
 
};
 
} // namespace elecstate
 
#endif