abacus-mc/source/source_lcao/module_dftu/dftu_force.h at 46d8278b85dcee03c1bf07b25958e50a18f18ec3 · mohanchen/abacus-mc · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
/// @file dftu_force.h
/// @brief Free-function helpers for DFT+U force and stress, extracted from
///        Plus_U. The top-level force_stress takes a Plus_U& because it needs
///        to call Plus_U::pot_onsite_real/complex; the four inner
///        functions are fully decoupled and take their dependencies as
///        explicit parameters (mirroring the folding helpers in the same
///        DFTU_LCAO namespace).
#ifndef DFTU_FORCE_H
#define DFTU_FORCE_H

#include "source_basis/module_ao/parallel_orbitals.h"
#include "source_base/matrix.h"
#include "source_base/vector3.h"
#include "source_cell/klist.h"
#include "source_cell/module_neighbor/sltk_grid_driver.h"
#include "source_cell/unitcell.h"
#include "source_lcao/force_stress_arrays.h"

#include <complex>
#include <string>
#include <vector>

#ifdef __LCAO

class Plus_U;

namespace DFTU_LCAO {

/// @brief Top-level entry: drives force/stress from DFT+U.
/// Takes Plus_U& because it calls dftu.pot_onsite_real/complex,
/// which are still members of Plus_U (defined in dftu_tools.cpp).
void force_stress(Plus_U& dftu,
                  const bool cal_force,
                  const bool cal_stress,
                  const UnitCell& ucell,
                  const Grid_Driver& gd,
                  std::vector<std::vector<double>>* dmk_d,
                  std::vector<std::vector<std::complex<double>>>* dmk_c,
                  const Parallel_Orbitals& pv,
                  ForceStressArrays& fsr,
                  ModuleBase::matrix& force_dftu,
                  ModuleBase::matrix& stress_dftu,
                  const K_Vectors& kv,
                  const int npol,
                  const bool gamma_only_local);

/// @brief Force contribution at a k-point (multik path).
void cal_force_k(int nlocal,
                 int npol,
                 const std::string& ks_solver,
                 const std::vector<double>& orb_cutoff,
                 const std::vector<int>& orbital_corr,
                 const std::vector<std::vector<std::vector<std::vector<std::vector<int>>>>>& iatlnmipol2iwt,
                 const UnitCell& ucell,
                 const Grid_Driver& gd,
                 ForceStressArrays& fsr,
                 const Parallel_Orbitals& pv,
                 const int ik,
                 const std::complex<double>* rho_VU,
                 ModuleBase::matrix& force_dftu,
                 const ModuleBase::Vector3<double>& kvec_d);

/// @brief Stress contribution at a k-point (multik path).
void cal_stress_k(int nlocal,
                  int npol,
                  const std::string& ks_solver,
                  const std::vector<double>& orb_cutoff,
                  const UnitCell& ucell,
                  const Grid_Driver& gd,
                  ForceStressArrays& fsr,
                  const Parallel_Orbitals& pv,
                  const int ik,
                  const std::complex<double>* rho_VU,
                  ModuleBase::matrix& stress_dftu,
                  const ModuleBase::Vector3<double>& kvec_d);

/// @brief Force contribution at gamma point.
void cal_force_gamma(int nlocal,
                     int npol,
                     const std::vector<int>& orbital_corr,
                     const std::vector<std::vector<std::vector<std::vector<std::vector<int>>>>>& iatlnmipol2iwt,
                     const UnitCell& ucell,
                     const double* rho_VU,
                     const Parallel_Orbitals& pv,
                     double* dsloc_x,
                     double* dsloc_y,
                     double* dsloc_z,
                     ModuleBase::matrix& force_dftu);

/// @brief Stress contribution at gamma point.
void cal_stress_gamma(int nlocal,
                      int npol,
                      const std::string& ks_solver,
                      const std::vector<double>& orb_cutoff,
                      const UnitCell& ucell,
                      const Parallel_Orbitals& pv,
                      const Grid_Driver* gd,
                      double* dsloc_x,
                      double* dsloc_y,
                      double* dsloc_z,
                      double* dh_r,
                      const double* rho_VU,
                      ModuleBase::matrix& stress_dftu);

} // namespace DFTU_LCAO

#endif // __LCAO

#endif // DFTU_FORCE_H