Skip to content

Commit

Permalink
Modified BasisDofHandler to account for no-dof processors from dealii
Browse files Browse the repository at this point in the history
  • Loading branch information
Avirup Sircar committed Jul 3, 2024
1 parent 6c88150 commit efd3aa4
Show file tree
Hide file tree
Showing 14 changed files with 566 additions and 114 deletions.
54 changes: 27 additions & 27 deletions analysis/classicalEnrichmentComparison/KSDFTClassical/TestKohnShamDft.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -353,13 +353,6 @@ int main(int argc, char** argv)
// evaluate basis data
feBasisData->evaluateBasisData(quadAttr, basisAttrMap);

std::shared_ptr<const quadrature::QuadratureRuleContainer> quadRuleContainer =
feBasisData->getQuadratureRuleContainer();

// scale the electronic charges
quadrature::QuadratureValuesContainer<double, Host>
electronChargeDensity(quadRuleContainer, 1, 0.0);

std::map<std::string, std::string> atomSymbolToFilename;
for (auto i:atomSymbolVec )
{
Expand All @@ -375,21 +368,7 @@ int main(int argc, char** argv)

std::shared_ptr<const utils::ScalarSpatialFunctionReal> rho = std::make_shared
<RhoFunction>(atomSphericalDataContainer, atomSymbolVec, atomChargesVec, atomCoordinatesVec);

for (size_type iCell = 0; iCell < electronChargeDensity.nCells(); iCell++)
{
for (size_type iComp = 0; iComp < 1; iComp++)
{
size_type quadId = 0;
std::vector<double> a(
electronChargeDensity.nCellQuadraturePoints(iCell));
a = (*rho)(quadRuleContainer->getCellRealPoints(iCell));
double *b = a.data();
electronChargeDensity.template
setCellQuadValues<Host>(iCell, iComp, b);
}
}


std::shared_ptr<const utils::ScalarSpatialFunctionReal>
zeroFunction = std::make_shared
<utils::ScalarZeroFunctionReal>();
Expand Down Expand Up @@ -475,15 +454,11 @@ std::shared_ptr<linearAlgebra::OperatorContext<double,
rootCout << "Entering KohnSham DFT Class....\n\n";

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::StoreGradient] = true;
basisAttrMap[basis::BasisStorageAttributes::StoreHessian] = false;
basisAttrMap[basis::BasisStorageAttributes::StoreOverlap] = false;
basisAttrMap[basis::BasisStorageAttributes::StoreGradNiGradNj] = false;
Expand All @@ -496,6 +471,31 @@ std::shared_ptr<linearAlgebra::OperatorContext<double,
// evaluate basis data
feBDElectrostaticsHamiltonian->evaluateBasisData(quadAttrVCorrecPlusPhi, basisAttrMap);

std::shared_ptr<const quadrature::QuadratureRuleContainer> quadRuleContainer =
feBDElectrostaticsHamiltonian->getQuadratureRuleContainer();

// scale the electronic charges
quadrature::QuadratureValuesContainer<double, Host>
electronChargeDensity(quadRuleContainer, 1, 0.0);

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

for (size_type iCell = 0; iCell < electronChargeDensity.nCells(); iCell++)
{
for (size_type iComp = 0; iComp < 1; iComp++)
{
size_type quadId = 0;
std::vector<double> a(
electronChargeDensity.nCellQuadraturePoints(iCell));
a = (*rho)(quadRuleContainer->getCellRealPoints(iCell));
double *b = a.data();
electronChargeDensity.template
setCellQuadValues<Host>(iCell, iComp, b);
}
}

if(isNumericalNuclearSolve)
{

Expand Down
30 changes: 15 additions & 15 deletions analysis/classicalEnrichmentComparison/KSDFTClassical/param.in
View file
Original file line number Diff line number Diff line change
@@ -1,27 +1,27 @@
xmax = 28.
ymax = 28.
zmax = 28.
xmax = 24.
ymax = 24.
zmax = 24.
radiusAtAtom = 0
meshSizeAtAtom = 0.5
radiusAroundAtom = 5.0
meshSizeAroundAtom = 0.5
meshSizeAtAtom = 0.01
radiusAroundAtom = 2.0
meshSizeAroundAtom = 0.1
rc = 0.6
feOrder = 4
num1DGaussSize = 5
num1DGLLSize = 5
feOrder = 5
num1DGaussSize = 6
num1DGLLSize = 6
smearingTemperature = 500.
fermiEnergyTolerance = 1e-8
fracOccupancyTolerance = 1e-2
eigenSolveResidualTolerance = 1e-2
chebyshevPolynomialDegree = 500
maxChebyshevFilterPass = 100
fracOccupancyTolerance = 1e-3
eigenSolveResidualTolerance = 1e-3
chebyshevPolynomialDegree = 1250
maxChebyshevFilterPass = 40
numWantedEigenvalues = 20
scfDensityResidualNormTolerance = 1e-5
maxSCFIter = 80
mixingHistory = 10
mixingParameter = 0.2
isAdaptiveAndersonMixingParameter = 0
evaluateEnergyEverySCF = 1
num1DGaussSizeVCorrecPlusPhi = 20
num1DGaussSizeVCorrecPlusPhi = 6
isNumericalNuclearSolve = 1
num1DGaussSizeSmearNucl = 20
num1DGaussSizeSmearNucl = 6
137 changes: 132 additions & 5 deletions analysis/classicalEnrichmentComparison/KSDFTOrthoEFE/TestKohnShamDft.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -39,6 +39,9 @@

#include <iostream>

#include <chrono>
using namespace std::chrono;

using namespace dftefe;
const utils::MemorySpace Host = utils::MemorySpace::HOST;

Expand Down Expand Up @@ -312,6 +315,9 @@ int main(int argc, char** argv)

utils::mpi::MPIComm comm = utils::mpi::MPICommWorld;

utils::mpi::MPIBarrier(comm);
auto startTotal = std::chrono::high_resolution_clock::now();

// Get the rank of the process
int rank;
utils::mpi::MPICommRank(comm, &rank);
Expand Down Expand Up @@ -567,6 +573,10 @@ int main(int argc, char** argv)

std::shared_ptr<basis::ParentToChildCellsManagerBase> parentToChildCellsManager = std::make_shared<basis::ParentToChildCellsManagerDealii<dim>>();

// add device synchronize for gpu
utils::mpi::MPIBarrier(comm);
auto start = std::chrono::high_resolution_clock::now();

std::shared_ptr<quadrature::QuadratureRuleContainer> quadRuleContainerAdaptiveElec =
std::make_shared<quadrature::QuadratureRuleContainer>
(quadAttrAdaptive,
Expand All @@ -581,6 +591,12 @@ int main(int argc, char** argv)
smallestCellVolume,
maxRecursion);

// add device synchronize for gpu
utils::mpi::MPIBarrier(comm);
auto stop = std::chrono::high_resolution_clock::now();

auto duration = std::chrono::duration_cast<std::chrono::microseconds>(stop - start);

unsigned int nQuad = quadRuleContainerAdaptiveElec->nQuadraturePoints();
int mpierr = utils::mpi::MPIAllreduce<Host>(
utils::mpi::MPIInPlace,
Expand All @@ -592,6 +608,8 @@ int main(int argc, char** argv)

rootCout << "Number of quadrature points in electrostatics adaptive quadrature: "<< nQuad<<"\n";

rootCout << "Time for adaptive quadrature creation is(in secs) : " << duration.count()/1e6 << std::endl;

// Compute Adaptive QuadratureRuleContainer for wavefn
/*
// Set up the vector of scalarSpatialRealFunctions for adaptive quadrature
Expand Down Expand Up @@ -636,6 +654,10 @@ int main(int argc, char** argv)

rootCout << "Number of quadrature points in wave function adaptive quadrature: "<<nQuad<<"\n";

// add device synchronize for gpu
utils::mpi::MPIBarrier(comm);
start = std::chrono::high_resolution_clock::now();

// Make orthogonalized EFE basis for all the fields

// 1. Make CFEBasisDataStorageDealii object for Rhs (ADAPTIVE with GAUSS and fns are N_i^2 - make quadrulecontainer), overlapmatrix (GLL)
Expand Down Expand Up @@ -725,6 +747,14 @@ int main(int argc, char** argv)
std::make_shared<basis::EFEBasisDofHandlerDealii<double, double,Host,dim>>(
enrichClassIntfceOrbital, feOrder, comm);

// add device synchronize for gpu
utils::mpi::MPIBarrier(comm);
stop = std::chrono::high_resolution_clock::now();

duration = std::chrono::duration_cast<std::chrono::microseconds>(stop - start);

rootCout << "Time for Ortho EFE basis manager creation is(in secs) : " << duration.count()/1e6 << std::endl;

std::map<global_size_type, utils::Point> dofCoords;
basisDofHandlerTotalPot->getBasisCenters(dofCoords);

Expand All @@ -745,8 +775,22 @@ int main(int argc, char** argv)
std::make_shared<basis::EFEBasisDataStorageDealii<double, double, Host,dim>>
(basisDofHandlerTotalPot, quadAttrAdaptive, basisAttrMap);

rootCout << "Number of quadrature points in wave function adaptive quadrature: "<<nQuad<<"\n";

// add device synchronize for gpu
utils::mpi::MPIBarrier(comm);
start = std::chrono::high_resolution_clock::now();

efeBasisDataAdaptiveTotPot->evaluateBasisData(quadAttrAdaptive, quadRuleContainerAdaptiveOrbital, basisAttrMap);

// add device synchronize for gpu
utils::mpi::MPIBarrier(comm);
stop = std::chrono::high_resolution_clock::now();

duration = std::chrono::duration_cast<std::chrono::microseconds>(stop - start);

rootCout << "Time for electrostatics basis datastorage evaluation is(in secs) : " << duration.count()/1e6 << std::endl;

// std::ofstream fout("filename.txt");
// utils::ConditionalOStream allCout(fout);
// auto basisDataAdapElec = efeBasisDataAdaptiveTotPot->getBasisDataInAllCells();
Expand Down Expand Up @@ -781,8 +825,20 @@ int main(int argc, char** argv)
std::make_shared<basis::EFEBasisDataStorageDealii<double, double, Host,dim>>
(basisDofHandlerWaveFn, quadAttrAdaptive, basisAttrMap);

// add device synchronize for gpu
utils::mpi::MPIBarrier(comm);
start = std::chrono::high_resolution_clock::now();

efeBasisDataAdaptiveOrbital->evaluateBasisData(quadAttrAdaptive, quadRuleContainerAdaptiveOrbital, basisAttrMap);

// add device synchronize for gpu
utils::mpi::MPIBarrier(comm);
stop = std::chrono::high_resolution_clock::now();

duration = std::chrono::duration_cast<std::chrono::microseconds>(stop - start);

rootCout << "Time for orbital basis datastorage evaluation is(in secs) : " << duration.count()/1e6 << std::endl;

std::shared_ptr<const quadrature::QuadratureRuleContainer> quadRuleContainerRho =
efeBasisDataAdaptiveOrbital->getQuadratureRuleContainer();

Expand Down Expand Up @@ -831,6 +887,12 @@ int main(int argc, char** argv)
const utils::ScalarSpatialFunctionReal *externalPotentialFunction = new
utils::PointChargePotentialFunction(atomCoordinatesVec, atomChargesVec);

rootCout << "Number of quadrature points in wave function adaptive quadrature: "<<nQuad<<"\n";

// add device synchronize for gpu
utils::mpi::MPIBarrier(comm);
start = std::chrono::high_resolution_clock::now();

// Create OperatorContext for Basisoverlap

std::shared_ptr<const basis::OrthoEFEOverlapOperatorContext<double,
Expand All @@ -848,6 +910,18 @@ int main(int argc, char** argv)
*cfeBasisDataStorageGLL,
50);

// add device synchronize for gpu
utils::mpi::MPIBarrier(comm);
stop = std::chrono::high_resolution_clock::now();

duration = std::chrono::duration_cast<std::chrono::microseconds>(stop - start);

rootCout << "Time for creation of MContext is(in secs) : " << duration.count()/1e6 << std::endl;

// add device synchronize for gpu
utils::mpi::MPIBarrier(comm);
start = std::chrono::high_resolution_clock::now();

std::shared_ptr<const basis::OrthoEFEOverlapOperatorContext<double,
double,
Host,
Expand All @@ -861,7 +935,19 @@ int main(int argc, char** argv)
*cfeBasisDataStorageGLL,
*efeBasisDataAdaptiveOrbital,
*cfeBasisDataStorageGLL,
50);
50);

// add device synchronize for gpu
utils::mpi::MPIBarrier(comm);
stop = std::chrono::high_resolution_clock::now();

duration = std::chrono::duration_cast<std::chrono::microseconds>(stop - start);

rootCout << "Time for creation of MContextForInv is(in secs) : " << duration.count()/1e6 << std::endl;

// add device synchronize for gpu
utils::mpi::MPIBarrier(comm);
start = std::chrono::high_resolution_clock::now();

std::shared_ptr<linearAlgebra::OperatorContext<double,
double,
Expand All @@ -873,10 +959,17 @@ int main(int argc, char** argv)
(*basisManagerWaveFn,
*cfeBasisDataStorageGLL,
*efeBasisDataAdaptiveOrbital,
linAlgOpContext);
linAlgOpContext);

rootCout << "Entering KohnSham DFT Class....\n\n";
// add device synchronize for gpu
utils::mpi::MPIBarrier(comm);
stop = std::chrono::high_resolution_clock::now();

duration = std::chrono::duration_cast<std::chrono::microseconds>(stop - start);

rootCout << "Time for creation of MInv is(in secs) : " << duration.count()/1e6 << std::endl;

rootCout << "Entering KohnSham DFT Class....\n\n";

std::shared_ptr<const basis::FEBasisDataStorage<double, Host>> feBDTotalChargeStiffnessMatrix = efeBasisDataAdaptiveTotPot;
std::shared_ptr<const basis::FEBasisDataStorage<double, Host>> feBDTotalChargeRhs = efeBasisDataAdaptiveTotPot;
Expand Down Expand Up @@ -971,7 +1064,20 @@ int main(int argc, char** argv)
*MContext,
*MInvContext);

dftefeSolve->solve();
// add device synchronize for gpu
utils::mpi::MPIBarrier(comm);
start = std::chrono::high_resolution_clock::now();

dftefeSolve->solve();

// add device synchronize for gpu
utils::mpi::MPIBarrier(comm);
stop = std::chrono::high_resolution_clock::now();

duration = std::chrono::duration_cast<std::chrono::microseconds>(stop - start);

rootCout << "Time for scf iterations is(in secs) : " << duration.count()/1e6 << std::endl;

}
else
{
Expand Down Expand Up @@ -1020,9 +1126,30 @@ int main(int argc, char** argv)
*MContext,
*MInvContext);

dftefeSolve->solve();
// add device synchronize for gpu
utils::mpi::MPIBarrier(comm);
start = std::chrono::high_resolution_clock::now();

dftefeSolve->solve();

// add device synchronize for gpu
utils::mpi::MPIBarrier(comm);
stop = std::chrono::high_resolution_clock::now();

duration = std::chrono::duration_cast<std::chrono::microseconds>(stop - start);

rootCout << "Time for scf iterations is(in secs) : " << duration.count()/1e6 << std::endl;

}

// add device synchronize for gpu
utils::mpi::MPIBarrier(comm);
auto stopTotal = std::chrono::high_resolution_clock::now();

auto durationTotal = std::chrono::duration_cast<std::chrono::microseconds>(stopTotal - startTotal);

rootCout << "Total wall time(in secs) : " << durationTotal.count()/1e6 << std::endl;

//gracefully end MPI

int mpiFinalFlag = 0;
Expand Down
Loading

0 comments on commit efd3aa4

Please sign in to comment.