Skip to content

Commit

Permalink
TestKohnShamDFTClassical passes
Browse files Browse the repository at this point in the history
  • Loading branch information
Avirup Sircar committed Jun 6, 2024
1 parent 114c3c8 commit b509805
Show file tree
Hide file tree
Showing 7 changed files with 1,171 additions and 16 deletions.
4 changes: 3 additions & 1 deletion src/ksdft/ElectrostaticAllElectronFE.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -253,7 +253,9 @@ namespace dftefe
dim>>
d_feBasisOpHamiltonian;
std::shared_ptr<linearAlgebra::LinAlgOpContext<memorySpace>>
d_linAlgOpContext;
d_linAlgOpContext;
std::vector<RealType> d_nuclearChargeQuad;

}; // end of class ElectrostaticAllElectronFE
} // end of namespace ksdft
} // end of namespace dftefe
Expand Down
39 changes: 33 additions & 6 deletions src/ksdft/ElectrostaticAllElectronFE.t.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -428,17 +428,20 @@ namespace dftefe
d_atomCharges,
d_smearedChargeRadius);

RealType totNuclearChargeQuad = 0;
for (size_type iCell = 0; iCell < quadRuleContainer->nCells(); iCell++)
{
for (size_type iComp = 0; iComp < d_numComponents; iComp++)
{
size_type quadId = 0;
std::vector<double> jxw = quadRuleContainer->getCellJxW(iCell);
std::vector<RealType> a(
quadRuleContainer->nCellQuadraturePoints(iCell));
for (auto j : quadRuleContainer->getCellRealPoints(iCell))
{
a[quadId] = (RealType)smfunc(j);
quadId = quadId + 1;
totNuclearChargeQuad += (RealType)smfunc(j) * jxw[quadId];
quadId = quadId + 1;
}
RealType *b = a.data();
d_nuclearChargesDensity
Expand All @@ -448,6 +451,22 @@ namespace dftefe
}
}

utils::mpi::MPIAllreduce<utils::MemorySpace::HOST>(
utils::mpi::MPIInPlace,
&totNuclearChargeQuad,
1,
utils::mpi::Types<RealType>::getMPIDatatype(),
utils::mpi::MPISum,
d_feBMTotalCharge->getMPIPatternP2P()->mpiCommunicator());

double totalAtomCharges =
std::accumulate(d_atomCharges.begin(), d_atomCharges.end(), (double)0);
// d_nuclearChargesDensity by totalAtomCharges/totNuclearChargeQuad
quadrature::scale((RealType)std::abs(totalAtomCharges /
totNuclearChargeQuad),
*d_nuclearChargesDensity,
*d_linAlgOpContext);

// create the correction quadValuesContainer for numerical solve

nuclearPotentialSolve(feBDNuclearChargeStiffnessMatrix,
Expand Down Expand Up @@ -805,6 +824,7 @@ namespace dftefe
quadrature::QuadratureValuesContainer<RealType, memorySpace>
nuclearChargeDensity(quadRuleContainer, d_numComponents);

d_nuclearChargeQuad.resize(d_numAtoms, 0);
for (unsigned int iAtom = 0; iAtom < d_numAtoms; iAtom++)
{
std::shared_ptr<const utils::ScalarSpatialFunctionReal> smfunc =
Expand All @@ -825,7 +845,7 @@ namespace dftefe
d_atomCharges[iAtom],
d_smearedChargeRadius[iAtom]);

RealType nuclearChargeQuad = 0;
d_nuclearChargeQuad[iAtom] = 0;
for (size_type iCell = 0; iCell < quadRuleContainer->nCells();
iCell++)
{
Expand All @@ -839,7 +859,8 @@ namespace dftefe
for (auto j : quadRuleContainer->getCellRealPoints(iCell))
{
a[quadId] = (RealType)(*smfunc)(j);
nuclearChargeQuad += (RealType)(*smfunc)(j)*jxw[quadId];
d_nuclearChargeQuad[iAtom] +=
(RealType)(*smfunc)(j)*jxw[quadId];
quadId = quadId + 1;
}
RealType *b = a.data();
Expand All @@ -852,16 +873,16 @@ namespace dftefe

utils::mpi::MPIAllreduce<utils::MemorySpace::HOST>(
utils::mpi::MPIInPlace,
&nuclearChargeQuad,
&d_nuclearChargeQuad[iAtom],
1,
utils::mpi::Types<RealType>::getMPIDatatype(),
utils::mpi::MPISum,
d_feBMTotalCharge->getMPIPatternP2P()->mpiCommunicator());

// Scale by 4\pi * d_atomCharges[iAtom]/nuclearChargeQuad
// Scale by 4\pi * d_atomCharges[iAtom]/d_nuclearChargeQuad[iAtom]
quadrature::scale((RealType)std::abs(4 * utils::mathConstants::pi *
d_atomCharges[iAtom] /
nuclearChargeQuad),
d_nuclearChargeQuad[iAtom]),
nuclearChargeDensity,
*d_linAlgOpContext);

Expand Down Expand Up @@ -997,6 +1018,12 @@ namespace dftefe
}
}

// Scale by d_atomCharges[iAtom]/d_nuclearChargeQuad[iAtom]
quadrature::scale((RealType)std::abs(d_atomCharges[iAtom] /
d_nuclearChargeQuad[iAtom]),
nuclearChargeDensity,
*d_linAlgOpContext);

basis::FEBasisOperations<ValueTypeBasisCoeff,
ValueTypeBasisData,
memorySpace,
Expand Down
69 changes: 69 additions & 0 deletions src/ksdft/KohnShamDFT.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -71,6 +71,75 @@ namespace dftefe
memorySpace>::OpContext;

public:
KohnShamDFT(
/* Atom related info */
const std::vector<utils::Point> &atomCoordinates,
const std::vector<double> & atomCharges,
const std::vector<double> & smearedChargeRadius,
const size_type numElectrons,
/* SCF related info */
const size_type numWantedEigenvalues,
const double smearingTemperature,
const double fermiEnergyTolerance,
const double fracOccupancyTolerance,
const double eigenSolveResidualTolerance,
const double scfDensityResidualNormTolerance,
const size_type chebyshevPolynomialDegree,
const size_type maxChebyshevFilterPass,
const size_type maxSCFIter,
const bool evaluateEnergyEverySCF,
/* Mixing related info */
const size_type mixingHistory,
const double mixingParameter,
const bool isAdaptiveAndersonMixingParameter,
/* Electron density related info */
const quadrature::QuadratureValuesContainer<RealType, memorySpace>
&electronChargeDensityInput,
/* Basis related info */
/* Field boundary */
std::shared_ptr<
const basis::FEBasisManager<ValueTypeElectrostaticsCoeff,
ValueTypeElectrostaticsBasis,
memorySpace,
dim>> feBMTotalCharge,
std::shared_ptr<const basis::FEBasisManager<ValueTypeWaveFunctionCoeff,
ValueTypeWaveFunctionBasis,
memorySpace,
dim>> feBMWaveFn,
/* Field data storages */
std::shared_ptr<
const basis::FEBasisDataStorage<ValueTypeElectrostaticsBasis,
memorySpace>>
feBDTotalChargeStiffnessMatrix,
std::shared_ptr<
const basis::FEBasisDataStorage<ValueTypeElectrostaticsBasis,
memorySpace>> feBDTotalChargeRhs,
std::shared_ptr<
const basis::FEBasisDataStorage<ValueTypeWaveFunctionBasis,
memorySpace>> feBDHamiltonian,
/* PSP/AE related info */
const utils::ScalarSpatialFunctionReal &externalPotentialFunction,
/* linAgOperations Context*/
std::shared_ptr<linearAlgebra::LinAlgOpContext<memorySpace>>
linAlgOpContext,
/* Batch size related info */
const size_type cellBlockSize,
const size_type waveFunctionBatchSize,
/* basis overlap related info */
const OpContext &MContextForInv =
linearAlgebra::IdentityOperatorContext<ValueTypeOperator,
ValueTypeOperand,
memorySpace>(),
const OpContext &MContext =
linearAlgebra::IdentityOperatorContext<ValueTypeOperator,
ValueTypeOperand,
memorySpace>(),
const OpContext &MInvContext =
linearAlgebra::IdentityOperatorContext<ValueTypeOperator,
ValueTypeOperand,
memorySpace>());


KohnShamDFT(
/* Atom related info */
const std::vector<utils::Point> &atomCoordinates,
Expand Down
Loading

0 comments on commit b509805

Please sign in to comment.