Skip to content

Commit

Permalink
Modified ksdft to take individual hamiltonian quadrature; modified po…
Browse files Browse the repository at this point in the history
…issonsolver to remove gradNigradNj basisdatastorage dependence
  • Loading branch information
Avirup Sircar committed Jun 17, 2024
1 parent 7e2a21f commit a235350
Show file tree
Hide file tree
Showing 12 changed files with 805 additions and 200 deletions.
19 changes: 17 additions & 2 deletions analysis/classicalEnrichmentComparison/H.xml

Large diffs are not rendered by default.

Original file line number Diff line number Diff line change
Expand Up @@ -185,6 +185,9 @@ int main(int argc, char** argv)
bool evaluateEnergyEverySCF = readParameter<bool>(parameterInputFileName, "evaluateEnergyEverySCF", rootCout);
const size_type dim = 3;

unsigned int num1DGaussSizeVCorrecPlusPhi = readParameter<unsigned int>(parameterInputFileName, "num1DGaussSizeVCorrecPlusPhi", rootCout);
bool isNumericalNuclearSolve = readParameter<bool>(parameterInputFileName, "isNumericalNuclearSolve", rootCout);

// Set up Triangulation
std::shared_ptr<basis::TriangulationBase> triangulationBase =
std::make_shared<basis::TriangulationDealiiParallel<dim>>(comm);
Expand Down Expand Up @@ -276,7 +279,7 @@ int main(int argc, char** argv)
basisAttrMap[basis::BasisStorageAttributes::StoreGradient] = true;
basisAttrMap[basis::BasisStorageAttributes::StoreHessian] = false;
basisAttrMap[basis::BasisStorageAttributes::StoreOverlap] = false;
basisAttrMap[basis::BasisStorageAttributes::StoreGradNiGradNj] = true;
basisAttrMap[basis::BasisStorageAttributes::StoreGradNiGradNj] = false;
basisAttrMap[basis::BasisStorageAttributes::StoreJxW] = true;

// Set up the FE Basis Data Storage
Expand Down Expand Up @@ -397,52 +400,150 @@ std::shared_ptr<linearAlgebra::OperatorContext<double,
*feBasisDataGLL,
linAlgOpContext);
rootCout << "Entering KohnSham DFT Class....\n\n";
std::shared_ptr<ksdft::KohnShamDFT<double,
double,
double,
double,
Host,
dim>> dftefeSolve =
std::make_shared<ksdft::KohnShamDFT<double,

std::shared_ptr<const basis::FEBasisDataStorage<double, Host>> feBDTotalChargeStiffnessMatrix = feBasisData;
std::shared_ptr<const basis::FEBasisDataStorage<double, Host>> feBDTotalChargeRhs = feBasisData;

std::shared_ptr<const basis::FEBasisDataStorage<double,Host>> feBDKineticHamiltonian = feBasisData;
std::shared_ptr<const basis::FEBasisDataStorage<double, Host>> feBDEXCHamiltonian = feBasisData;

quadrature::QuadratureRuleAttributes quadAttrVCorrecPlusPhi(quadrature::QuadratureFamily::GAUSS,true,num1DGaussSizeVCorrecPlusPhi);

basisAttrMap[basis::BasisStorageAttributes::StoreValues] = true;
basisAttrMap[basis::BasisStorageAttributes::StoreGradient] = false;
basisAttrMap[basis::BasisStorageAttributes::StoreHessian] = false;
basisAttrMap[basis::BasisStorageAttributes::StoreOverlap] = false;
basisAttrMap[basis::BasisStorageAttributes::StoreGradNiGradNj] = false;
basisAttrMap[basis::BasisStorageAttributes::StoreJxW] = true;

std::shared_ptr<basis::FEBasisDataStorage<double, Host>> feBDElectrostaticsHamiltonian =
std::make_shared<basis::CFEBasisDataStorageDealii<double, double, Host,dim>>
(basisDofHandler, quadAttrVCorrecPlusPhi, basisAttrMap);

// evaluate basis data
feBDElectrostaticsHamiltonian->evaluateBasisData(quadAttrVCorrecPlusPhi, basisAttrMap);

if(isNumericalNuclearSolve)
{

unsigned int num1DGaussSizeSmearNucl = readParameter<unsigned int>(parameterInputFileName, "num1DGaussSizeSmearNucl", rootCout);
quadrature::QuadratureRuleAttributes quadAttrSmearNucl(quadrature::QuadratureFamily::GAUSS,true,num1DGaussSizeSmearNucl);

basisAttrMap[basis::BasisStorageAttributes::StoreValues] = true;
basisAttrMap[basis::BasisStorageAttributes::StoreGradient] = false;
basisAttrMap[basis::BasisStorageAttributes::StoreHessian] = false;
basisAttrMap[basis::BasisStorageAttributes::StoreOverlap] = false;
basisAttrMap[basis::BasisStorageAttributes::StoreGradNiGradNj] = false;
basisAttrMap[basis::BasisStorageAttributes::StoreJxW] = true;

// Set up the FE Basis Data Storage
std::shared_ptr<basis::FEBasisDataStorage<double, Host>> feBDNuclearChargeRhs =
std::make_shared<basis::CFEBasisDataStorageDealii<double, double, Host,dim>>
(basisDofHandler, quadAttrSmearNucl, basisAttrMap);

// evaluate basis data
feBDNuclearChargeRhs->evaluateBasisData(quadAttrSmearNucl, basisAttrMap);

std::shared_ptr<const basis::FEBasisDataStorage<double,Host>> feBDNuclearChargeStiffnessMatrix = feBasisData;

std::shared_ptr<ksdft::KohnShamDFT<double,
double,
double,
double,
Host,
dim>> dftefeSolve =
std::make_shared<ksdft::KohnShamDFT<double,
double,
double,
double,
Host,
dim>>(
atomCoordinatesVec,
atomChargesVec,
smearedChargeRadiusVec,
numElectrons,
numWantedEigenvalues,
smearingTemperature,
fermiEnergyTolerance,
fracOccupancyTolerance,
eigenSolveResidualTolerance,
scfDensityResidualNormTolerance,
chebyshevPolynomialDegree,
maxChebyshevFilterPass,
maxSCFIter,
evaluateEnergyEverySCF,
mixingHistory,
mixingParameter,
isAdaptiveAndersonMixingParameter,
electronChargeDensity,
basisManagerTotalPot,
basisManagerWaveFn,
feBDTotalChargeStiffnessMatrix,
feBDTotalChargeRhs,
feBDNuclearChargeStiffnessMatrix,
feBDNuclearChargeRhs,
feBDKineticHamiltonian,
feBDElectrostaticsHamiltonian,
feBDEXCHamiltonian,
*externalPotentialFunction,
linAlgOpContext,
50,
50,
*MContextForInv,
*MContext,
*MInvContext);

dftefeSolve->solve();
}
else
{
std::shared_ptr<ksdft::KohnShamDFT<double,
double,
double,
double,
Host,
dim>>(
atomCoordinatesVec,
atomChargesVec,
smearedChargeRadiusVec,
numElectrons,
numWantedEigenvalues,
smearingTemperature,
fermiEnergyTolerance,
fracOccupancyTolerance,
eigenSolveResidualTolerance,
scfDensityResidualNormTolerance,
chebyshevPolynomialDegree,
maxChebyshevFilterPass,
maxSCFIter,
evaluateEnergyEverySCF,
mixingHistory,
mixingParameter,
isAdaptiveAndersonMixingParameter,
electronChargeDensity,
basisManagerTotalPot,
basisManagerWaveFn,
feBasisData,
feBasisData,
feBasisData,
feBasisData,
feBasisData,
*externalPotentialFunction,
linAlgOpContext,
50,
50,
*MContextForInv,
*MContext,
*MInvContext);

dftefeSolve->solve();
dim>> dftefeSolve =
std::make_shared<ksdft::KohnShamDFT<double,
double,
double,
double,
Host,
dim>>(
atomCoordinatesVec,
atomChargesVec,
smearedChargeRadiusVec,
numElectrons,
numWantedEigenvalues,
smearingTemperature,
fermiEnergyTolerance,
fracOccupancyTolerance,
eigenSolveResidualTolerance,
scfDensityResidualNormTolerance,
chebyshevPolynomialDegree,
maxChebyshevFilterPass,
maxSCFIter,
evaluateEnergyEverySCF,
mixingHistory,
mixingParameter,
isAdaptiveAndersonMixingParameter,
electronChargeDensity,
basisManagerTotalPot,
basisManagerWaveFn,
feBDTotalChargeStiffnessMatrix,
feBDTotalChargeRhs,
feBDKineticHamiltonian,
feBDElectrostaticsHamiltonian,
feBDEXCHamiltonian,
*externalPotentialFunction,
linAlgOpContext,
50,
50,
*MContextForInv,
*MContext,
*MInvContext);

dftefeSolve->solve();
}

//gracefully end MPI

Expand Down
Original file line number Diff line number Diff line change
@@ -1,7 +1,7 @@
xmax = 24.
ymax = 24.
zmax = 24.
radiusAtAtom = 1.0
radiusAtAtom = 0
meshSizeAtAtom = 0.1
radiusAroundAtom = 2.0
meshSizeAroundAtom = 0.5
Expand All @@ -22,4 +22,7 @@ maxSCFIter = 80
mixingHistory = 10
mixingParameter = 0.2
isAdaptiveAndersonMixingParameter = 1
evaluateEnergyEverySCF = 1
evaluateEnergyEverySCF = 1
num1DGaussSizeVCorrecPlusPhi = 6
isNumericalNuclearSolve = 1
num1DGaussSizeSmearNucl = 6
Loading

0 comments on commit a235350

Please sign in to comment.