Skip to content

Commit cfb949b

Browse files
knikhil1995dsambit
knikhil1995
authored and
dsambit
committed
Merged in BFGSPreconOpt (pull request #553)
Fixed preconditioner construction for geometry optimization Approved-by: Sambit Das
2 parents 54f0a77 + aafcaed commit cfb949b

File tree

6 files changed

+131
-47
lines changed

6 files changed

+131
-47
lines changed

‎include/dft.h

+22-8
Original file line numberDiff line numberDiff line change
@@ -262,10 +262,10 @@ namespace dftfe
262262
double
263263
getFreeEnergy() const;
264264

265-
distributedCPUVec<double>
265+
const distributedCPUVec<double> &
266266
getRhoNodalOut() const;
267267

268-
distributedCPUVec<double>
268+
const distributedCPUVec<double> &
269269
getRhoNodalSplitOut() const;
270270

271271
double
@@ -389,14 +389,28 @@ namespace dftfe
389389
* @brief Gets the current atom Locations in cartesian form
390390
* (origin at center of domain) from dftClass
391391
*/
392-
std::vector<std::vector<double>>
392+
const std::vector<std::vector<double>> &
393393
getAtomLocationsCart() const;
394394

395+
396+
/**
397+
* @brief Gets the current image atom Locations in cartesian form
398+
* (origin at center of domain) from dftClass
399+
*/
400+
const std::vector<std::vector<double>> &
401+
getImageAtomLocationsCart() const;
402+
403+
/**
404+
* @brief Gets the current image atom ids from dftClass
405+
*/
406+
const std::vector<int> &
407+
getImageAtomIDs() const;
408+
395409
/**
396410
* @brief Gets the current atom Locations in fractional form
397411
* from dftClass (only applicable for periodic and semi-periodic BCs)
398412
*/
399-
std::vector<std::vector<double>>
413+
const std::vector<std::vector<double>> &
400414
getAtomLocationsFrac() const;
401415

402416

@@ -407,7 +421,7 @@ namespace dftfe
407421
* @return std::vector<std::vector<double>> 3 \times 3 matrix,lattice[i][j]
408422
* corresponds to jth component of ith lattice vector
409423
*/
410-
std::vector<std::vector<double>>
424+
const std::vector<std::vector<double>> &
411425
getCell() const;
412426

413427
/**
@@ -420,19 +434,19 @@ namespace dftfe
420434
/**
421435
* @brief Gets the current atom types from dftClass
422436
*/
423-
std::set<unsigned int>
437+
const std::set<unsigned int> &
424438
getAtomTypes() const;
425439

426440
/**
427441
* @brief Gets the current atomic forces from dftClass
428442
*/
429-
std::vector<double>
443+
const std::vector<double> &
430444
getForceonAtoms() const;
431445

432446
/**
433447
* @brief Gets the current cell stress from dftClass
434448
*/
435-
dealii::Tensor<2, 3, double>
449+
const dealii::Tensor<2, 3, double> &
436450
getCellStress() const;
437451

438452
/**

‎include/dftBase.h

+22-8
Original file line numberDiff line numberDiff line change
@@ -88,10 +88,10 @@ namespace dftfe
8888
virtual double
8989
getFreeEnergy() const = 0;
9090

91-
virtual distributedCPUVec<double>
91+
virtual const distributedCPUVec<double> &
9292
getRhoNodalOut() const = 0;
9393

94-
virtual distributedCPUVec<double>
94+
virtual const distributedCPUVec<double> &
9595
getRhoNodalSplitOut() const = 0;
9696

9797
virtual double
@@ -107,14 +107,28 @@ namespace dftfe
107107
* @brief Gets the current atom Locations in cartesian form
108108
* (origin at center of domain) from dftClass
109109
*/
110-
virtual std::vector<std::vector<double>>
110+
virtual const std::vector<std::vector<double>> &
111111
getAtomLocationsCart() const = 0;
112112

113+
/**
114+
* @brief Gets the current image atom Locations in cartesian form
115+
* (origin at center of domain) from dftClass
116+
*/
117+
virtual const std::vector<std::vector<double>> &
118+
getImageAtomLocationsCart() const = 0;
119+
120+
/**
121+
* @brief Gets the current image atom ids from dftClass
122+
*/
123+
virtual const std::vector<int> &
124+
getImageAtomIDs() const = 0;
125+
126+
113127
/**
114128
* @brief Gets the current atom Locations in fractional form
115129
* from dftClass (only applicable for periodic and semi-periodic BCs)
116130
*/
117-
virtual std::vector<std::vector<double>>
131+
virtual const std::vector<std::vector<double>> &
118132
getAtomLocationsFrac() const = 0;
119133

120134

@@ -125,7 +139,7 @@ namespace dftfe
125139
* @return std::vector<std::vector<double>> 3 \times 3 matrix,lattice[i][j]
126140
* corresponds to jth component of ith lattice vector
127141
*/
128-
virtual std::vector<std::vector<double>>
142+
virtual const std::vector<std::vector<double>> &
129143
getCell() const = 0;
130144

131145
/**
@@ -138,20 +152,20 @@ namespace dftfe
138152
/**
139153
* @brief Gets the current atom types from dftClass
140154
*/
141-
virtual std::set<unsigned int>
155+
virtual const std::set<unsigned int> &
142156
getAtomTypes() const = 0;
143157

144158
/**
145159
* @brief Gets the current atomic forces (configurational forces) from dftClass
146160
*/
147-
virtual std::vector<double>
161+
virtual const std::vector<double> &
148162
getForceonAtoms() const = 0;
149163

150164

151165
/**
152166
* @brief Gets the current cell stress from dftClass
153167
*/
154-
virtual dealii::Tensor<2, 3, double>
168+
virtual const dealii::Tensor<2, 3, double> &
155169
getCellStress() const = 0;
156170

157171
/**

‎include/force.h

+2-2
Original file line numberDiff line numberDiff line change
@@ -171,7 +171,7 @@ namespace dftfe
171171
* the three force components on each atom being the leading dimension.
172172
* Units- Hartree/Bohr
173173
*/
174-
std::vector<double>
174+
std::vector<double> &
175175
getAtomsForces();
176176

177177
/** @brief prints the currently stored configurational forces on atoms and the Gaussian generator constant
@@ -244,7 +244,7 @@ namespace dftfe
244244
*
245245
* @return dealii::Tensor<2,3,double> second order stress Tensor in Hartree/Bohr^3
246246
*/
247-
dealii::Tensor<2, 3, double>
247+
dealii::Tensor<2, 3, double> &
248248
getStress();
249249

250250
/** @brief get the value of Gaussian generator parameter (d_gaussianConstant).

‎src/dft/dft.cc

+26-8
Original file line numberDiff line numberDiff line change
@@ -4419,7 +4419,7 @@ namespace dftfe
44194419
template <unsigned int FEOrder,
44204420
unsigned int FEOrderElectro,
44214421
dftfe::utils::MemorySpace memorySpace>
4422-
std::vector<std::vector<double>>
4422+
const std::vector<std::vector<double>> &
44234423
dftClass<FEOrder, FEOrderElectro, memorySpace>::getAtomLocationsCart() const
44244424
{
44254425
return atomLocations;
@@ -4428,7 +4428,25 @@ namespace dftfe
44284428
template <unsigned int FEOrder,
44294429
unsigned int FEOrderElectro,
44304430
dftfe::utils::MemorySpace memorySpace>
4431-
std::vector<std::vector<double>>
4431+
const std::vector<std::vector<double>> &
4432+
dftClass<FEOrder, FEOrderElectro, memorySpace>::getImageAtomLocationsCart() const
4433+
{
4434+
return d_imagePositionsTrunc;
4435+
}
4436+
4437+
template <unsigned int FEOrder,
4438+
unsigned int FEOrderElectro,
4439+
dftfe::utils::MemorySpace memorySpace>
4440+
const std::vector<int> &
4441+
dftClass<FEOrder, FEOrderElectro, memorySpace>::getImageAtomIDs() const
4442+
{
4443+
return d_imageIdsTrunc;
4444+
}
4445+
4446+
template <unsigned int FEOrder,
4447+
unsigned int FEOrderElectro,
4448+
dftfe::utils::MemorySpace memorySpace>
4449+
const std::vector<std::vector<double>> &
44324450
dftClass<FEOrder, FEOrderElectro, memorySpace>::getAtomLocationsFrac() const
44334451
{
44344452
return atomLocationsFractional;
@@ -4437,7 +4455,7 @@ namespace dftfe
44374455
template <unsigned int FEOrder,
44384456
unsigned int FEOrderElectro,
44394457
dftfe::utils::MemorySpace memorySpace>
4440-
std::vector<std::vector<double>>
4458+
const std::vector<std::vector<double>> &
44414459
dftClass<FEOrder, FEOrderElectro, memorySpace>::getCell() const
44424460
{
44434461
return d_domainBoundingVectors;
@@ -4456,7 +4474,7 @@ namespace dftfe
44564474
template <unsigned int FEOrder,
44574475
unsigned int FEOrderElectro,
44584476
dftfe::utils::MemorySpace memorySpace>
4459-
std::set<unsigned int>
4477+
const std::set<unsigned int> &
44604478
dftClass<FEOrder, FEOrderElectro, memorySpace>::getAtomTypes() const
44614479
{
44624480
return atomTypes;
@@ -4465,7 +4483,7 @@ namespace dftfe
44654483
template <unsigned int FEOrder,
44664484
unsigned int FEOrderElectro,
44674485
dftfe::utils::MemorySpace memorySpace>
4468-
std::vector<double>
4486+
const std::vector<double> &
44694487
dftClass<FEOrder, FEOrderElectro, memorySpace>::getForceonAtoms() const
44704488
{
44714489
return (forcePtr->getAtomsForces());
@@ -4474,7 +4492,7 @@ namespace dftfe
44744492
template <unsigned int FEOrder,
44754493
unsigned int FEOrderElectro,
44764494
dftfe::utils::MemorySpace memorySpace>
4477-
dealii::Tensor<2, 3, double>
4495+
const dealii::Tensor<2, 3, double> &
44784496
dftClass<FEOrder, FEOrderElectro, memorySpace>::getCellStress() const
44794497
{
44804498
return (forcePtr->getStress());
@@ -4519,7 +4537,7 @@ namespace dftfe
45194537
template <unsigned int FEOrder,
45204538
unsigned int FEOrderElectro,
45214539
dftfe::utils::MemorySpace memorySpace>
4522-
distributedCPUVec<double>
4540+
const distributedCPUVec<double> &
45234541
dftClass<FEOrder, FEOrderElectro, memorySpace>::getRhoNodalOut() const
45244542
{
45254543
return d_densityOutNodalValues[0];
@@ -4528,7 +4546,7 @@ namespace dftfe
45284546
template <unsigned int FEOrder,
45294547
unsigned int FEOrderElectro,
45304548
dftfe::utils::MemorySpace memorySpace>
4531-
distributedCPUVec<double>
4549+
const distributedCPUVec<double> &
45324550
dftClass<FEOrder, FEOrderElectro, memorySpace>::getRhoNodalSplitOut() const
45334551
{
45344552
return d_rhoOutNodalValuesSplit;

‎src/force/force.cc

+2-2
Original file line numberDiff line numberDiff line change
@@ -569,7 +569,7 @@ namespace dftfe
569569
template <unsigned int FEOrder,
570570
unsigned int FEOrderElectro,
571571
dftfe::utils::MemorySpace memorySpace>
572-
std::vector<double>
572+
std::vector<double> &
573573
forceClass<FEOrder, FEOrderElectro, memorySpace>::getAtomsForces()
574574
{
575575
return d_globalAtomsForces;
@@ -578,7 +578,7 @@ namespace dftfe
578578
template <unsigned int FEOrder,
579579
unsigned int FEOrderElectro,
580580
dftfe::utils::MemorySpace memorySpace>
581-
dealii::Tensor<2, 3, double>
581+
dealii::Tensor<2, 3, double> &
582582
forceClass<FEOrder, FEOrderElectro, memorySpace>::getStress()
583583
{
584584
return d_stress;

‎src/geoOpt/geoOptIon.cc

+57-19
Original file line numberDiff line numberDiff line change
@@ -522,21 +522,37 @@ namespace dftfe
522522
else
523523
{
524524
const int numberGlobalAtoms = d_dftPtr->getAtomLocationsCart().size();
525+
const int numberImageAtoms =
526+
d_dftPtr->getImageAtomLocationsCart().size();
525527
std::vector<std::vector<double>> NNdistances(numberGlobalAtoms);
526528
double rNN = 0;
527529
for (int i = 0; i < numberGlobalAtoms; ++i)
528530
{
529531
double riMin = 0;
530-
for (int j = 0; j < numberGlobalAtoms; ++j)
532+
for (int j = 0; j < numberGlobalAtoms + numberImageAtoms; ++j)
531533
{
532534
double rij = 0;
533-
for (int k = 2; k < 5; ++k)
534-
{
535-
rij += (d_dftPtr->getAtomLocationsCart()[i][k] -
536-
d_dftPtr->getAtomLocationsCart()[j][k]) *
537-
(d_dftPtr->getAtomLocationsCart()[i][k] -
538-
d_dftPtr->getAtomLocationsCart()[j][k]);
539-
}
535+
if (j < numberGlobalAtoms)
536+
for (int k = 2; k < 5; ++k)
537+
{
538+
rij += (d_dftPtr->getAtomLocationsCart()[i][k] -
539+
d_dftPtr->getAtomLocationsCart()[j][k]) *
540+
(d_dftPtr->getAtomLocationsCart()[i][k] -
541+
d_dftPtr->getAtomLocationsCart()[j][k]);
542+
}
543+
else
544+
for (int k = 2; k < 5; ++k)
545+
{
546+
rij +=
547+
(d_dftPtr->getAtomLocationsCart()[i][k] -
548+
d_dftPtr
549+
->getImageAtomLocationsCart()[j - numberGlobalAtoms]
550+
[k - 2]) *
551+
(d_dftPtr->getAtomLocationsCart()[i][k] -
552+
d_dftPtr
553+
->getImageAtomLocationsCart()[j - numberGlobalAtoms]
554+
[k - 2]);
555+
}
540556
rij = std::sqrt(rij);
541557
if ((riMin > rij && i != j) || j == 0)
542558
{
@@ -554,22 +570,44 @@ namespace dftfe
554570
std::vector<double> L(numberGlobalAtoms * numberGlobalAtoms, 0.0);
555571
for (int i = 0; i < numberGlobalAtoms; ++i)
556572
{
557-
for (int j = i + 1; j < numberGlobalAtoms; ++j)
573+
for (int j = i + 1; j < numberGlobalAtoms + numberImageAtoms; ++j)
558574
{
559575
double rij = 0;
560-
for (int k = 2; k < 5; ++k)
561-
{
562-
rij += (d_dftPtr->getAtomLocationsCart()[i][k] -
563-
d_dftPtr->getAtomLocationsCart()[j][k]) *
564-
(d_dftPtr->getAtomLocationsCart()[i][k] -
565-
d_dftPtr->getAtomLocationsCart()[j][k]);
566-
}
576+
int jatomId =
577+
j < numberGlobalAtoms ?
578+
j :
579+
d_dftPtr->getImageAtomIDs()[j - numberGlobalAtoms];
580+
if (j < numberGlobalAtoms)
581+
for (int k = 2; k < 5; ++k)
582+
{
583+
rij += (d_dftPtr->getAtomLocationsCart()[i][k] -
584+
d_dftPtr->getAtomLocationsCart()[j][k]) *
585+
(d_dftPtr->getAtomLocationsCart()[i][k] -
586+
d_dftPtr->getAtomLocationsCart()[j][k]);
587+
}
588+
else
589+
for (int k = 2; k < 5; ++k)
590+
{
591+
rij +=
592+
(d_dftPtr->getAtomLocationsCart()[i][k] -
593+
d_dftPtr
594+
->getImageAtomLocationsCart()[j - numberGlobalAtoms]
595+
[k - 2]) *
596+
(d_dftPtr->getAtomLocationsCart()[i][k] -
597+
d_dftPtr
598+
->getImageAtomLocationsCart()[j - numberGlobalAtoms]
599+
[k - 2]);
600+
}
567601
rij = std::sqrt(rij);
568602
if (rij < rCut)
569603
{
570-
L[i * numberGlobalAtoms + j] =
571-
-std::exp(-3.0 * (rij / rNN - 1));
572-
L[j * numberGlobalAtoms + i] = L[i * numberGlobalAtoms + j];
604+
double preconVal = -std::exp(-3.0 * (rij / rNN - 1));
605+
if (preconVal < L[i * numberGlobalAtoms + jatomId])
606+
{
607+
L[i * numberGlobalAtoms + jatomId] = preconVal;
608+
L[jatomId * numberGlobalAtoms + i] =
609+
L[i * numberGlobalAtoms + jatomId];
610+
}
573611
}
574612
}
575613
}

0 commit comments

Comments
 (0)