Skip to content

Commit 20e6fbc

Browse files
author
Avirup Sircar
committed
Minor changes in AtomSphDataJSON
1 parent 65dbdb1 commit 20e6fbc

10 files changed

Lines changed: 234 additions & 118 deletions

File tree

‎analysis/classicalEnrichmentComparison/PSPnonlocal/KSDFTClassicalUniformQuad/CMakeLists.txt‎

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -6,6 +6,8 @@ set(CMAKE_CXX_STANDARD_REQUIRED ON)
66
set(MAIN_PROJECT_DIR $ENV{DFTEFE_PATH})
77
message(${MAIN_PROJECT_DIR})
88

9+
list(APPEND CMAKE_MODULE_PATH ${MAIN_PROJECT_DIR}/cmake)
10+
911
include_directories(${MAIN_PROJECT_DIR}/src)
1012
add_subdirectory(${MAIN_PROJECT_DIR}/src/utils ${MAIN_PROJECT_DIR}/test/basis/lib/buildUtils)
1113
add_subdirectory(${MAIN_PROJECT_DIR}/src/basis ${MAIN_PROJECT_DIR}/test/basis/lib/buildBasis)

‎analysis/classicalEnrichmentComparison/PSPnonlocal/KSDFTClassicalUniformQuad/TestKohnShamDft.cpp‎

Lines changed: 113 additions & 54 deletions
Original file line numberDiff line numberDiff line change
@@ -43,39 +43,62 @@ using namespace dftefe;
4343
const utils::MemorySpace Host = utils::MemorySpace::HOST;
4444

4545
template<typename T>
46-
T readParameter(std::string ParamFile, std::string param, utils::ConditionalOStream &rootCout)
46+
T readParameter(const std::string &ParamFile,
47+
const std::string &param,
48+
utils::ConditionalOStream &rootCout,
49+
bool throwIfEmpty = true,
50+
bool throwIfMissing = true,
51+
T defaultValue = T{})
4752
{
48-
T t(0);
53+
T t = defaultValue;
4954
std::string line;
50-
std::fstream fstream;
51-
fstream.open(ParamFile, std::fstream::in);
52-
int count = 0;
55+
std::fstream fstream(ParamFile, std::fstream::in);
56+
bool found = false;
57+
5358
while (std::getline(fstream, line))
5459
{
55-
for (int i = 0; i < line.length(); i++)
60+
auto pos = line.find('=');
61+
if (pos == std::string::npos) continue;
62+
63+
std::string key = line.substr(0, pos);
64+
// trim spaces
65+
key.erase(std::remove_if(key.begin(), key.end(), ::isspace), key.end());
66+
67+
if (key == param)
5668
{
57-
if (line[i] == ' ')
58-
{
59-
line.erase(line.begin() + i);
60-
i--;
69+
found = true;
70+
std::string value = line.substr(pos + 1);
71+
// trim leading spaces
72+
value.erase(value.begin(),
73+
std::find_if(value.begin(), value.end(),
74+
[](unsigned char ch){ return !std::isspace(ch); }));
75+
76+
if (value.empty()) {
77+
if (throwIfEmpty) {
78+
utils::throwException(false, "Parameter found but empty: " + param);
6179
}
62-
}
63-
std::istringstream iss(line);
64-
std::string type;
65-
std::getline(iss, type, '=');
66-
if (type.compare(param) == 0)
67-
{
68-
iss >> t;
69-
count = 1;
80+
t = defaultValue;
81+
} else {
82+
if constexpr (std::is_same<T, std::string>::value) {
83+
t = value;
84+
} else {
85+
std::istringstream iss(value);
86+
iss >> t;
87+
}
88+
}
7089
break;
7190
}
7291
}
73-
if(count == 0)
74-
{
75-
utils::throwException(false, "The parameter is not found: "+ param);
92+
93+
if (!found) {
94+
if (throwIfMissing) {
95+
utils::throwException(false, "The parameter is not found: " + param);
96+
}
97+
t = defaultValue;
7698
}
99+
77100
fstream.close();
78-
rootCout << "Reading parameter -- " << param << " = "<<t<<std::endl;
101+
rootCout << "Reading parameter -- " << param << " = " << t << std::endl;
79102
return t;
80103
}
81104

@@ -317,12 +340,9 @@ int main(int argc, char** argv)
317340
utils::throwException(false,
318341
"dftefe_path does not exist!");
319342
}
320-
std::string atomDataFile = argv[1];
321-
std::string inputFileName = sourceDir + atomDataFile;
322-
std::string paramDataFile = argv[2];
343+
std::string paramDataFile = argv[1];
323344
std::string parameterInputFileName = sourceDir + paramDataFile;
324345

325-
rootCout << "Reading input file: "<<inputFileName<<std::endl;
326346
rootCout << "Reading parameter file: "<<parameterInputFileName<<std::endl;
327347

328348
// Read parameters
@@ -353,13 +373,22 @@ int main(int argc, char** argv)
353373
double atomPartitionTolerance = readParameter<double>(parameterInputFileName, "atomPartitionTolerance", rootCout);
354374
unsigned int num1DGaussSubdividedSizeElec = readParameter<unsigned int>(parameterInputFileName, "num1DGaussSubdividedSizeElec", rootCout);
355375
unsigned int gaussSubdividedCopiesElec = readParameter<unsigned int>(parameterInputFileName, "gaussSubdividedCopiesElec", rootCout);
376+
356377
unsigned int num1DGaussSubdividedSizeEigen = readParameter<unsigned int>(parameterInputFileName, "num1DGaussSubdividedSizeEigen", rootCout);
357378
unsigned int gaussSubdividedCopiesEigen = readParameter<unsigned int>(parameterInputFileName, "gaussSubdividedCopiesEigen", rootCout);
358379
bool isNumericalNuclearSolve = readParameter<bool>(parameterInputFileName, "isNumericalNuclearSolve", rootCout);
359380
bool isDeltaRhoPoissonSolve = readParameter<bool>(parameterInputFileName, "isDeltaRhoPoissonSolve", rootCout);
360381

361382
unsigned int num1DGaussSubdividedSizeGrad = readParameter<unsigned int>(parameterInputFileName, "num1DGaussSubdividedSizeGrad", rootCout);
362383
unsigned int gaussSubdividedCopiesGrad = readParameter<unsigned int>(parameterInputFileName, "gaussSubdividedCopiesGrad", rootCout);
384+
std::string coordinatesDataFile = readParameter<std::string>(parameterInputFileName, "coordinatesDataFile", rootCout , false , false);
385+
std::string basisDataFile = readParameter<std::string>(parameterInputFileName, "basisDataFile", rootCout , false , false);
386+
std::string PSPDataFile = readParameter<std::string>(parameterInputFileName, "PSPDataFile", rootCout , false , false);
387+
388+
std::string tciaFolder = readParameter<std::string>(parameterInputFileName, "tciaFolder", rootCout , false , false);
389+
std::string tciaOutFilePrefix = readParameter<std::string>(parameterInputFileName, "tciaOutFilePrefix", rootCout , false , false);
390+
391+
const atoms::TCIADataParams tciaparams{tciaFolder , tciaOutFilePrefix};
363392

364393
unsigned int num1DGaussSubdividedSizeNonLocOperator = 14;
365394
unsigned int gaussSubdividedCopiesNonLocOperator = 1;
@@ -375,47 +404,82 @@ int main(int argc, char** argv)
375404
domainVectors[2][2] = zmax;
376405

377406
std::fstream fstream;
378-
fstream.open(inputFileName, std::fstream::in);
379407

380408
// read the input file and create atomsymbol vector and atom coordinates vector.
381409
std::vector<utils::Point> atomCoordinatesVec(0,utils::Point(dim, 0.0));
382410
std::vector<double> coordinates;
383-
std::vector<std::string> pspFilePathVec(0);
384411
coordinates.resize(dim,0.);
385412
std::vector<std::string> atomSymbolVec(0);
386413
std::vector<double> atomChargesVec(0);
387-
std::string symbol;
388-
std::string pspFilePath;
414+
std::string symbol , basisFilePath, pspFilePath;
415+
std::map<std::string, double> atomSymbolToChargeMap;
389416
double valanceNumber;
390417
atomSymbolVec.resize(0);
391418
std::string line;
419+
420+
std::map<std::string, std::string> atomSymbolToBasisFileName;
421+
std::vector<std::string> matchString(0);
422+
fstream.open(basisDataFile, std::fstream::in);
392423
while (std::getline(fstream, line)){
393424
std::stringstream ss(line);
394425
ss >> symbol;
395-
ss >> valanceNumber;
426+
ss >> basisFilePath;
427+
atomSymbolToBasisFileName[symbol] = basisFilePath;
428+
if(std::find(matchString.begin(), matchString.end(), symbol) == matchString.end())
429+
matchString.push_back(symbol);
430+
else
431+
utils::throwException(false, "The atom Symbols were repeated for PSP filenames. ");
432+
}
433+
utils::mpi::MPIBarrier(comm);
434+
fstream.close();
435+
436+
std::map<std::string, std::string> atomSymbolToPSPFileName;
437+
matchString.clear();
438+
fstream.open(PSPDataFile, std::fstream::in);
439+
while (std::getline(fstream, line)){
440+
std::stringstream ss(line);
441+
ss >> symbol;
396442
ss >> pspFilePath;
443+
atomSymbolToPSPFileName[symbol] = pspFilePath;
444+
if(std::find(matchString.begin(), matchString.end(), symbol) == matchString.end())
445+
matchString.push_back(symbol);
446+
else
447+
utils::throwException(false, "The atom Symbols were repeated for PSP filenames. ");
448+
}
449+
utils::mpi::MPIBarrier(comm);
450+
fstream.close();
451+
452+
fstream.open(coordinatesDataFile, std::fstream::in);
453+
while (std::getline(fstream, line)){
454+
std::stringstream ss(line);
455+
ss >> symbol;
456+
ss >> valanceNumber;
397457
for(unsigned int i=0 ; i<dim ; i++){
398458
ss >> coordinates[i];
399459
}
400-
pspFilePathVec.push_back(pspFilePath);
401460
atomCoordinatesVec.push_back(coordinates);
402461
atomSymbolVec.push_back(symbol);
462+
if(atomSymbolToPSPFileName.find(symbol) == atomSymbolToPSPFileName.end())
463+
{
464+
utils::throwException(false, "PSP filename does not have the same atom symbol as Coordinate filename.");
465+
}
466+
if(atomSymbolToBasisFileName.find(symbol) == atomSymbolToBasisFileName.end())
467+
{
468+
utils::throwException(false, "Basis filename does not have the same atom symbol as Coordinate filename.");
469+
}
403470
atomChargesVec.push_back((-1.0)*valanceNumber);
471+
if(atomSymbolToChargeMap.find(symbol) == atomSymbolToChargeMap.end())
472+
atomSymbolToChargeMap[symbol] = valanceNumber;
404473
}
405474
utils::mpi::MPIBarrier(comm);
406475
fstream.close();
407476

408-
std::map<std::string, std::string> atomSymbolToPSPFilename;
409-
for (int i = 0 ; i < atomSymbolVec.size() ; i++)
410-
{
411-
atomSymbolToPSPFilename[atomSymbolVec[i]] = sourceDir + pspFilePathVec[i];
412-
}
413-
414477
size_type numElectrons = 0;
415478
for(auto &i : atomChargesVec)
416479
{
417480
numElectrons += (size_type)(std::abs(i));
418481
}
482+
419483
// Generate mesh
420484
std::shared_ptr<basis::CellMappingBase> cellMapping = std::make_shared<basis::LinearCellMappingDealii<dim>>();
421485

@@ -433,8 +497,6 @@ int main(int argc, char** argv)
433497

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

436-
std::vector<double> smearedChargeRadiusVec(atomCoordinatesVec.size(),rc);
437-
438500
// initialize the basis DofHandler
439501

440502
std::shared_ptr<const basis::FEBasisDofHandler<double, Host,dim>> basisDofHandlerTotalPot =
@@ -711,7 +773,7 @@ std::shared_ptr<linearAlgebra::OperatorContext<double,
711773
atomCoordinatesVec,
712774
atomChargesVec,
713775
atomSymbolVec,
714-
smearedChargeRadiusVec,
776+
rc,
715777
numElectrons,
716778
numWantedEigenvalues,
717779
smearingTemperature,
@@ -734,7 +796,7 @@ std::shared_ptr<linearAlgebra::OperatorContext<double,
734796
feBDElectrostaticsHamiltonian,
735797
feBDEXCHamiltonian,
736798
feBDAtomCenterNonLocalOperator,
737-
atomSymbolToPSPFilename,
799+
atomSymbolToPSPFileName,
738800
linAlgOpContext,
739801
*MContextForInv,
740802
/**MContextForInv,*/
@@ -743,20 +805,15 @@ std::shared_ptr<linearAlgebra::OperatorContext<double,
743805
}
744806
else if (!isNumericalNuclearSolve && isDeltaRhoPoissonSolve)
745807
{
746-
std::map<std::string, std::string> atomSymbolToFilename;
747-
for (auto i:atomSymbolVec )
748-
{
749-
atomSymbolToFilename[i] = sourceDir + i + ".xml";
750-
}
751-
752808
std::vector<std::string> fieldNames{"density", "vtotal"};
753-
std::vector<std::string> metadataNames{ "symbol", "Z", "charge", "NR", "r" };
809+
std::vector<std::string> metadataNames{ "symbol", "Z", "charge", "NR"};
754810
std::shared_ptr<atoms::AtomSphericalDataContainer> atomSphericalDataContainer =
755811
std::make_shared<atoms::AtomSphericalDataContainer>(
756812
atoms::AtomSphericalDataType::ENRICHMENT,
757-
atomSymbolToFilename,
813+
atomSymbolToBasisFileName,
758814
fieldNames,
759-
metadataNames);
815+
metadataNames,
816+
std::map<std::string, std::string>({{"rcsmear", std::to_string(rc)}, {"PSP/AE", "PSP"}}));
760817

761818
std::shared_ptr<utils::ScalarSpatialFunctionReal> smfuncAtTotPot =
762819
std::make_shared<AtomicTotalElectrostaticPotentialFunction>(atomSphericalDataContainer,
@@ -779,7 +836,7 @@ std::shared_ptr<linearAlgebra::OperatorContext<double,
779836
atomCoordinatesVec,
780837
atomChargesVec,
781838
atomSymbolVec,
782-
smearedChargeRadiusVec,
839+
rc,
783840
numElectrons,
784841
numWantedEigenvalues,
785842
smearingTemperature,
@@ -804,12 +861,14 @@ std::shared_ptr<linearAlgebra::OperatorContext<double,
804861
feBDElectrostaticsHamiltonian,
805862
feBDEXCHamiltonian,
806863
feBDAtomCenterNonLocalOperator,
807-
atomSymbolToPSPFilename,
864+
atomSymbolToPSPFileName,
808865
linAlgOpContext,
809866
*MContextForInv,
810867
/**MContextForInv,*/
811868
*MContext,
812-
*MInvContext);
869+
*MInvContext,
870+
true,
871+
tciaparams);
813872
}
814873
else
815874
{

0 commit comments

Comments
 (0)