Skip to content

Commit 43706df

Browse files
anirudhsg03phanimotamarri
anirudhsg03
authored and
phanimotamarri
committed
Merged in ldosPR (pull request #714)
LDOS Implementation Approved-by: Sambit Das Approved-by: Kartick Ramakrishnan Approved-by: Phani Motamarri
2 parents 4ef0937 + b3d08d0 commit 43706df

38 files changed

Lines changed: 3227 additions & 100 deletions

‎CMakeLists.txt‎

Lines changed: 3 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -156,6 +156,7 @@ SET(TARGET_SRC
156156
./src/poisson/poissonSolverProblem.cc
157157
./src/poisson/MultiVectorPoissonLinearSolverProblem.cpp
158158
./src/helmholtz/kerkerSolverProblem.cc
159+
./src/helmholtz/ldosSolverProblem.cc
159160
./src/dftOperator/KohnShamDFTBaseOperator.cc
160161
./src/dftOperator/KohnShamDFTStandardEigenOperator.cc
161162
./src/dftOperator/veffPrimeForLRDM.cc
@@ -229,6 +230,7 @@ SET(TARGET_SRC
229230
./src/dft/moveAtoms.cc
230231
./src/dft/moveMeshToAtoms.cc
231232
./src/dft/kerker.cc
233+
./src/dft/ldos.cc
232234
./src/dft/pRefinedDoFHandler.cc
233235
./src/dft/psiInitialGuess.cc
234236
./src/dft/publicMethods.cc
@@ -354,6 +356,7 @@ IF (WITH_GPU)
354356
./src/linAlg/rayleighRitzDevice.cc
355357
./src/poisson/poissonSolverProblemDevice.cc
356358
./src/helmholtz/kerkerSolverProblemDevice.cc
359+
./src/helmholtz/ldosSolverProblemDevice.cc
357360
./src/solvers/eigenSolvers/chebyshevOrthogonalizedSubspaceIterationSolverDevice.cc
358361
./src/solvers/linearSolverDevice.cc
359362
./src/solvers/linearSolverProblemDevice.cc

‎include/dftfe/MatrixFree.h‎

Lines changed: 13 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -80,6 +80,16 @@ namespace dftfe
8080
inline void
8181
initOperatorCoeffs(T coeffHelmholtz);
8282

83+
/**
84+
* @brief Initialize Helmholtz operator coefficient as a per-(cell,
85+
* quadrature-point) field. nQuadTotal must equal
86+
* nCells * nQuadPointsPerDim^3 and the source pointer must live in
87+
* @p memorySpace.
88+
*
89+
*/
90+
inline void
91+
initOperatorCoeffs(const T *coeffHelmholtz, dftfe::uInt nQuadTotal);
92+
8393
/**
8494
* @brief Compute Laplace operator multipled by X
8595
*
@@ -123,7 +133,9 @@ namespace dftfe
123133
d_localBlockSize, d_localSize, d_ghostBlockSize, d_ghostSize,
124134
d_nOMPThreads;
125135

126-
T d_coeffHelmholtz;
136+
T d_coeffHelmholtz;
137+
bool d_coeffHelmholtzIsPointer = false;
138+
dftfe::utils::MemoryStorage<T, memorySpace> d_coeffHelmholtzVec;
127139

128140
static constexpr std::uint32_t d_quadODim = nQuadPointsPerDim / 2;
129141
static constexpr std::uint32_t d_quadEDim =

‎include/dftfe/MatrixFreeDevice.h‎

Lines changed: 21 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -79,6 +79,27 @@ namespace dftfe
7979
dftfe::uInt nCells,
8080
dftfe::uInt nBatch);
8181

82+
static inline void
83+
computeHelmholtzX(T *dst,
84+
T *src,
85+
T *jacobianFactor,
86+
dftfe::uInt *map,
87+
T *shapeBuffer,
88+
T *coeffHelmholtz,
89+
dftfe::uInt nCells,
90+
dftfe::uInt nBatch);
91+
92+
template <typename CoeffOp>
93+
static inline void
94+
computeHelmholtzXImpl(T *dst,
95+
T *src,
96+
T *jacobianFactor,
97+
dftfe::uInt *map,
98+
T *shapeBuffer,
99+
CoeffOp coeffHelmholtz,
100+
dftfe::uInt nCells,
101+
dftfe::uInt nBatch);
102+
82103
static inline void
83104
constraintsDistribute(T *src,
84105
const dftfe::uInt *constrainingNodeBuckets,

‎include/dftfe/MatrixFreeWrapper.h‎

Lines changed: 13 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -175,6 +175,19 @@ namespace dftfe
175175
d_MatrixFreeObject);
176176
}
177177

178+
/**
179+
* @brief Initialize Helmholtz operator coefficient as a per-(cell,
180+
* quadrature-point) field.
181+
*
182+
*/
183+
inline void
184+
initOperatorCoeffs(const T *coeffHelmholtz, dftfe::uInt nQuadTotal)
185+
{
186+
std::visit(
187+
[&](auto &t) { t->initOperatorCoeffs(coeffHelmholtz, nQuadTotal); },
188+
d_MatrixFreeObject);
189+
}
190+
178191
/**
179192
* @brief Compute Laplace operator multipled by X
180193
*

‎include/dftfe/densityCalculator.h‎

Lines changed: 10 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -48,13 +48,16 @@ namespace dftfe
4848
&gradDensityValues,
4949
std::vector<
5050
dftfe::utils::MemoryStorage<double, dftfe::utils::MemorySpace::HOST>>
51-
&tauValues,
52-
const bool isEvaluateGradRho,
53-
const bool isEvaluateTau,
54-
const MPI_Comm &mpiCommParent,
55-
const MPI_Comm &interpoolcomm,
56-
const MPI_Comm &interBandGroupComm,
57-
const dftParameters &dftParams);
51+
&tauValues,
52+
const bool isEvaluateGradRho,
53+
const bool isEvaluateTau,
54+
const MPI_Comm &mpiCommParent,
55+
const MPI_Comm &interpoolcomm,
56+
const MPI_Comm &interBandGroupComm,
57+
const dftParameters &dftParams,
58+
const std::vector<std::vector<double>> *ldosOccupancies = nullptr,
59+
dftfe::utils::MemoryStorage<double, dftfe::utils::MemorySpace::HOST>
60+
*ldosQuadValues = nullptr);
5861

5962

6063
template <typename NumberType>

‎include/dftfe/dft.h‎

Lines changed: 27 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -73,6 +73,7 @@
7373
#include <dftfe/atomCenteredPostProcessing.h>
7474
#include <dftfe/poissonSolverProblemWrapper.h>
7575
#include <dftfe/kerkerSolverProblemWrapper.h>
76+
#include <dftfe/ldosSolverProblemWrapper.h>
7677
#include <dftfe/groupSymmetry.h>
7778

7879
namespace dftfe
@@ -552,6 +553,9 @@ namespace dftfe
552553
const std::vector<std::vector<double>> &eigenValuesInput,
553554
const double numElectronsInput);
554555

556+
void
557+
compute_ldosOccupanciesAndTotalDOS();
558+
555559
/**
556560
*@brief find HOMO eigenvalue for pure state
557561
*/
@@ -1179,6 +1183,9 @@ namespace dftfe
11791183
void
11801184
compute_rhoOut(const bool isGroundState = false);
11811185

1186+
void
1187+
compute_ldosOut();
1188+
11821189
/**
11831190
*@brief Mixing schemes for mixing electron-density
11841191
*/
@@ -1196,6 +1203,18 @@ namespace dftfe
11961203
distributedCPUVec<double> &residualRho,
11971204
distributedCPUVec<double> &preCondTotalDensityResidualVector);
11981205

1206+
void
1207+
applyLdosPreconditionerToTotalDensityResidual(
1208+
#ifdef DFTFE_WITH_DEVICE
1209+
ldosSolverProblemDeviceWrapperClass
1210+
&ldosPreconditionedResidualSolverProblemDevice,
1211+
linearSolverCGDevice &CGSolverDevice,
1212+
#endif
1213+
ldosSolverProblemWrapperClass &ldosPreconditionedResidualSolverProblem,
1214+
dealiiLinearSolver &CGSolver,
1215+
distributedCPUVec<double> &residualRho,
1216+
distributedCPUVec<double> &preCondTotalDensityResidualVector);
1217+
11991218
double
12001219
lowrankApproxScfDielectricMatrixInv(const dftfe::uInt scfIter);
12011220

@@ -1663,6 +1682,14 @@ namespace dftfe
16631682
std::vector<std::vector<double>> eigenValues;
16641683
std::vector<std::vector<double>> d_partialOccupancies;
16651684

1685+
std::vector<std::vector<double>> d_ldosOccupancies;
1686+
double d_totalDOS;
1687+
dftfe::utils::MemoryStorage<double, dftfe::utils::MemorySpace::HOST>
1688+
d_ldosQuadValues;
1689+
distributedCPUVec<double> d_ldosNodalValues;
1690+
dftfe::utils::MemoryStorage<double, dftfe::utils::MemorySpace::HOST>
1691+
d_ldosAxQuadValuesElectro;
1692+
16661693
/**
16671694
* data storage for the occupancy of Kohn-Sham wavefunctions
16681695
*/

‎include/dftfe/dftParameters.h‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -58,7 +58,7 @@ namespace dftfe
5858
spinMixingEnhancementFactor;
5959
bool adaptAndersonMixingParameter;
6060
double absLinearSolverTolerance, selfConsistentSolverTolerance, TVal,
61-
selfConsistentSolverEnergyTolerance, tot_magnetization,
61+
LDOSTemperature, selfConsistentSolverEnergyTolerance, tot_magnetization,
6262
absLinearSolverToleranceHelmholtz, smearTval, intervalSize;
6363

6464
bool useAtomicMagnetizationGuessConstraintMag;

‎include/dftfe/ldosSolverProblem.h‎

Lines changed: 192 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,192 @@
1+
#include <dftfe/dealiiLinearSolverProblem.h>
2+
#include <dftfe/triangulationManager.h>
3+
#include <dftfe/FEBasisOperations.h>
4+
5+
#ifndef ldosSolverProblem_H_
6+
# define ldosSolverProblem_H_
7+
8+
namespace dftfe
9+
{
10+
template <dftfe::uInt FEOrderElectro>
11+
class ldosSolverProblem : public dealiiLinearSolverProblem
12+
{
13+
public:
14+
/// Constructor
15+
ldosSolverProblem(const MPI_Comm &mpi_comm_parent,
16+
const MPI_Comm &mpi_comm_domain);
17+
18+
/**
19+
* @brief initialize the matrix-free data structures
20+
*
21+
*/
22+
void
23+
init(std::shared_ptr<
24+
dftfe::basis::
25+
FEBasisOperations<double, double, dftfe::utils::MemorySpace::HOST>>
26+
&basisOperationsPtr,
27+
dealii::AffineConstraints<double> &constraintMatrix,
28+
distributedCPUVec<double> &x,
29+
const dftfe::uInt matrixFreeVectorComponent,
30+
const dftfe::uInt matrixFreeQuadratureComponent,
31+
const dftfe::uInt matrixFreeAxQuadratureComponent,
32+
const bool isComputeMeanValueConstraint);
33+
34+
/**
35+
* @brief reinitialize data structures .
36+
*
37+
*/
38+
void
39+
reinit(
40+
distributedCPUVec<double> &x,
41+
const dftfe::utils::MemoryStorage<double, dftfe::utils::MemorySpace::HOST>
42+
&residualQuadValues,
43+
const dftfe::utils::MemoryStorage<double, dftfe::utils::MemorySpace::HOST>
44+
&ldosAxQuadValues);
45+
46+
/**
47+
* @brief get the reference to x field
48+
*
49+
*/
50+
distributedCPUVec<double> &
51+
getX();
52+
53+
/**
54+
* @brief Compute A matrix multipled by x.
55+
*
56+
*/
57+
void
58+
vmult(distributedCPUVec<double> &Ax, distributedCPUVec<double> &x);
59+
60+
/**
61+
* @brief Compute right hand side vector for the problem Ax = rhs.
62+
*
63+
*/
64+
void
65+
computeRhs(distributedCPUVec<double> &rhs);
66+
67+
/**
68+
* @brief Jacobi preconditioning.
69+
*
70+
*/
71+
void
72+
precondition_Jacobi(distributedCPUVec<double> &dst,
73+
const distributedCPUVec<double> &src,
74+
const double omega) const;
75+
76+
/**
77+
* @brief distribute x to the constrained nodes.
78+
*
79+
*/
80+
void
81+
distributeX();
82+
83+
/// function needed by dealii to mimic SparseMatrix for Jacobi
84+
/// preconditioning
85+
void
86+
subscribe(std::atomic<bool> *const validity,
87+
const std::string &identifier = "") const {};
88+
89+
/// function needed by dealii to mimic SparseMatrix for Jacobi
90+
/// preconditioning
91+
void
92+
unsubscribe(std::atomic<bool> *const validity,
93+
const std::string &identifier = "") const {};
94+
95+
/// function needed by dealii to mimic SparseMatrix
96+
bool
97+
operator!=(double val) const
98+
{
99+
return true;
100+
};
101+
102+
const distributedCPUVec<double> &
103+
getDlocMassVector() const;
104+
105+
double
106+
computeDlocIntegral(const distributedCPUVec<double> &phi) const;
107+
108+
double
109+
getTotalDOS() const;
110+
111+
private:
112+
void
113+
computeMeanValueConstraint();
114+
115+
void
116+
meanValueConstraintDistribute(distributedCPUVec<double> &vec) const;
117+
118+
void
119+
meanValueConstraintDistributeSlaveToMaster(
120+
distributedCPUVec<double> &vec) const;
121+
122+
void
123+
meanValueConstraintSetZero(distributedCPUVec<double> &vec) const;
124+
125+
/**
126+
* @brief required for the cell_loop operation in dealii's MatrixFree class
127+
*
128+
*/
129+
void
130+
AX(const dealii::MatrixFree<3, double> &matrixFreeData,
131+
distributedCPUVec<double> &dst,
132+
const distributedCPUVec<double> &src,
133+
const std::pair<dftfe::uInt, dftfe::uInt> &cell_range) const;
134+
135+
/**
136+
* @brief Compute the diagonal of A.
137+
*
138+
*/
139+
void
140+
computeDiagonalA();
141+
142+
void
143+
computeProjectedQuadToNodalField(
144+
const dftfe::utils::MemoryStorage<double, dftfe::utils::MemorySpace::HOST>
145+
&quadValues,
146+
distributedCPUVec<double> &nodalField);
147+
148+
/// storage for diagonal of the A matrix
149+
distributedCPUVec<double> d_diagonalA;
150+
distributedCPUVec<double> d_dlocMassVector;
151+
distributedCPUVec<double> d_meanValueConstraintVec;
152+
153+
bool d_isMeanValueConstraintComputed;
154+
dealii::types::global_dof_index d_meanValueConstraintNodeId;
155+
dftfe::uInt d_meanValueConstraintProcId;
156+
157+
/// pointer to the x vector being solved for
158+
distributedCPUVec<double> *d_xPtr;
159+
160+
/// Integral of the discrete LDOS field used in the Helmholtz operator.
161+
double d_totalDOS;
162+
163+
/// matrix free index required to access the DofHandler and
164+
/// dealii::AffineConstraints<double> objects corresponding to the problem
165+
dftfe::uInt d_matrixFreeVectorComponent;
166+
167+
/// matrix free quadrature index
168+
dftfe::uInt d_matrixFreeQuadratureComponent;
169+
dftfe::uInt d_matrixFreeAxQuadratureComponent;
170+
171+
/// pointer to electron density cell and grad residual data
172+
const dftfe::utils::MemoryStorage<double, dftfe::utils::MemorySpace::HOST>
173+
*d_residualQuadValuesPtr;
174+
const dftfe::utils::MemoryStorage<double, dftfe::utils::MemorySpace::HOST>
175+
*d_ldosAxQuadValuesPtr;
176+
const dealii::DoFHandler<3> *d_dofHandlerPRefinedPtr;
177+
const dealii::AffineConstraints<double> *d_constraintMatrixPRefinedPtr;
178+
const dealii::MatrixFree<3, double> *d_matrixFreeDataPRefinedPtr;
179+
std::shared_ptr<
180+
dftfe::basis::
181+
FEBasisOperations<double, double, dftfe::utils::MemorySpace::HOST>>
182+
d_basisOperationsPtr;
183+
184+
const MPI_Comm d_mpiCommParent;
185+
const MPI_Comm mpi_communicator;
186+
const dftfe::uInt n_mpi_processes;
187+
const dftfe::uInt this_mpi_process;
188+
dealii::ConditionalOStream pcout;
189+
};
190+
191+
} // namespace dftfe
192+
#endif // ldosSolverProblem_H_

0 commit comments

Comments
 (0)