FazBrowse GitHub Viewer | Trending |
URL:
| Home
Tools: [Download Repo ZIP]   [Original HTTPS Page]

Fix CPCM regression and Factory initialization · PCMSolver/pcmsolver@edcaf9d · GitHub

  1. PCMSolver
  2. pcmsolver

Repository navigation

Commit edcaf9d

Browse files
robertodr
committed
Fix CPCM regression and Factory initialization
* The Factory is no longer implemented as a Singleton, thus avoiding the need for static methods returning references to static, local data. These are, in fact, not cleaned up until program execution is fully completed. * Reverted use of LLt (Cholesky) to LDLt (robust Cholesky) in CPCMSolver. Addition tests have been introduced for CPCM.
1 parent 1864952 commit edcaf9d

21 files changed

Lines changed: 236 additions & 138 deletions

‎CHANGELOG.md‎

Lines changed: 10 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1,5 +1,13 @@
11
# Change Log
22

3+
## [Unreleased]
4+
5+
### Changed
6+
7+
- `Factory` is no longer implemented as a Singleton.
8+
- Revert to use [Robust Cholesky decomposition](https://eigen.tuxfamily.org/dox/classEigen_1_1LDLT.html)
9+
to compute the inverse of the S matrix in `CPCMSolver`.
10+
311
## [Version 1.1.8] - 2017-02-06
412

513
### Added
@@ -288,7 +296,8 @@
288296

289297
## v1.0.0 - 2014-09-30 [YANKED]
290298

291-
[Unreleased]: https://github.com/PCMSolver/pcmsolver/compare/v1.1.7...HEAD
299+
[Unreleased]: https://github.com/PCMSolver/pcmsolver/compare/v1.1.8...HEAD
300+
[Version 1.1.8]: https://github.com/PCMSolver/pcmsolver/compare/v1.1.7...v1.1.8
292301
[Version 1.1.7]: https://github.com/PCMSolver/pcmsolver/compare/v1.1.6...v1.1.7
293302
[Version 1.1.6]: https://github.com/PCMSolver/pcmsolver/compare/v1.1.5...v1.1.6
294303
[Version 1.1.5]: https://github.com/PCMSolver/pcmsolver/compare/v1.1.4...v1.1.5

‎src/bi_operators/BoundaryIntegralOperator.hpp‎

Lines changed: 8 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -60,24 +60,25 @@ inline IBoundaryIntegralOperator * createPurisima(const BIOperatorData & data) {
6060
return new Purisima(data.scaling);
6161
}
6262

63-
inline void bootstrapFactory() {
63+
inline Factory<IBoundaryIntegralOperator, BIOperatorData> bootstrapFactory() {
64+
Factory<IBoundaryIntegralOperator, BIOperatorData> factory_;
65+
6466
const bool registeredCollocation =
65-
Factory<IBoundaryIntegralOperator, BIOperatorData>::TheFactory()
66-
.registerObject("COLLOCATION", createCollocation);
67+
factory_.registerObject("COLLOCATION", createCollocation);
6768
if (!registeredCollocation)
6869
PCMSOLVER_ERROR("Subscription of collocation integrator to factory failed!");
6970

7071
const bool registeredNumerical =
71-
Factory<IBoundaryIntegralOperator, BIOperatorData>::TheFactory()
72-
.registerObject("NUMERICAL", createNumerical);
72+
factory_.registerObject("NUMERICAL", createNumerical);
7373
if (!registeredNumerical)
7474
PCMSOLVER_ERROR("Subscription of numerical integrator to factory failed!");
7575

7676
const bool registeredPurisima =
77-
Factory<IBoundaryIntegralOperator, BIOperatorData>::TheFactory()
78-
.registerObject("PURISIMA", createPurisima);
77+
factory_.registerObject("PURISIMA", createPurisima);
7978
if (!registeredPurisima)
8079
PCMSOLVER_ERROR("Subscription of Purisima integrator to factory failed!");
80+
81+
return factory_;
8182
}
8283
} // namespace bi_operators
8384
} // namespace pcm

‎src/cavity/Cavity.hpp‎

Lines changed: 7 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -53,18 +53,19 @@ inline ICavity * createRestartCavity(const CavityData & data) {
5353
return new RestartCavity(data.filename);
5454
}
5555

56-
inline void bootstrapFactory() {
57-
const bool registeredGePol =
58-
Factory<ICavity, CavityData>::TheFactory().registerObject("GEPOL",
59-
createGePolCavity);
56+
inline Factory<ICavity, CavityData> bootstrapFactory() {
57+
Factory<ICavity, CavityData> factory_;
58+
59+
const bool registeredGePol = factory_.registerObject("GEPOL", createGePolCavity);
6060
if (!registeredGePol)
6161
PCMSOLVER_ERROR("Subscription of GePol cavity to factory failed!");
6262

6363
const bool registeredRestart =
64-
Factory<ICavity, CavityData>::TheFactory().registerObject("RESTART",
65-
createRestartCavity);
64+
factory_.registerObject("RESTART", createRestartCavity);
6665
if (!registeredRestart)
6766
PCMSOLVER_ERROR("Subscription of restart cavity to factory failed!");
67+
68+
return factory_;
6869
}
6970
} // namespace cavity
7071
} // namespace pcm

‎src/green/Green.hpp‎

Lines changed: 10 additions & 12 deletions
Original file line numberDiff line numberDiff line change
@@ -74,40 +74,38 @@ inline IGreensFunction * createSphericalDiffuse(const GreenData & data) {
7474
build, data, data.howProfile);
7575
}
7676

77-
inline void bootstrapFactory() {
78-
const bool registeredVacuum =
79-
Factory<IGreensFunction, GreenData>::TheFactory().registerObject("VACUUM",
80-
createVacuum);
77+
inline Factory<IGreensFunction, GreenData> bootstrapFactory() {
78+
Factory<IGreensFunction, GreenData> factory_;
79+
80+
const bool registeredVacuum = factory_.registerObject("VACUUM", createVacuum);
8181
if (!registeredVacuum)
8282
PCMSOLVER_ERROR("Subscription of vacuum Green's function to factory failed!");
8383

8484
const bool registeredUniformDielectric =
85-
Factory<IGreensFunction, GreenData>::TheFactory().registerObject(
86-
"UNIFORMDIELECTRIC", createUniformDielectric);
85+
factory_.registerObject("UNIFORMDIELECTRIC", createUniformDielectric);
8786
if (!registeredUniformDielectric)
8887
PCMSOLVER_ERROR(
8988
"Subscription of uniform dielectric Green's function to factory failed!");
9089

9190
const bool registeredSphericalDiffuse =
92-
Factory<IGreensFunction, GreenData>::TheFactory().registerObject(
93-
"SPHERICALDIFFUSE", createSphericalDiffuse);
91+
factory_.registerObject("SPHERICALDIFFUSE", createSphericalDiffuse);
9492
if (!registeredSphericalDiffuse)
9593
PCMSOLVER_ERROR(
9694
"Subscription of spherical diffuse Green's function to factory failed!");
9795

9896
const bool registeredIonicLiquid =
99-
Factory<IGreensFunction, GreenData>::TheFactory().registerObject(
100-
"IONICLIQUID", createIonicLiquid);
97+
factory_.registerObject("IONICLIQUID", createIonicLiquid);
10198
if (!registeredIonicLiquid)
10299
PCMSOLVER_ERROR(
103100
"Subscription of ionic liquid Green's function to factory failed!");
104101

105102
const bool registeredAnisotropicLiquid =
106-
Factory<IGreensFunction, GreenData>::TheFactory().registerObject(
107-
"ANISOTROPICLIQUID", createAnisotropicLiquid);
103+
factory_.registerObject("ANISOTROPICLIQUID", createAnisotropicLiquid);
108104
if (!registeredAnisotropicLiquid)
109105
PCMSOLVER_ERROR(
110106
"Subscription of anisotropic liquid Green's function to factory failed!");
107+
108+
return factory_;
111109
}
112110
} // namespace green
113111
} // namespace pcm

‎src/interface/Input.cpp‎

Lines changed: 2 additions & 10 deletions
Original file line numberDiff line numberDiff line change
@@ -270,16 +270,8 @@ void Input::initMolecule() {
270270
std::vector<Atom> radiiSet;
271271
std::vector<Atom> atoms;
272272
atoms.reserve(nuclei);
273-
// FIXME This was globally bootstrapped...
274-
using utils::Factory;
275-
radiiSet = Factory::TheFactory().create(radiiSet_);
276-
if (radiiSet_ == "UFF") {
277-
radiiSetName_ = "UFF";
278-
} else if (radiiSet_ == "BONDI") {
279-
radiiSetName_ = "Bondi-Mantina";
280-
} else {
281-
radiiSetName_ = "Allinger's MM3";
282-
}
273+
// FIXME Code duplication in function initMolecule in interface/Meddle.cpp
274+
tie(radiiSetName_, radiiSet) = utils::bootstrapRadiiSet().create(radiiSet_);
283275
for (int i = 0; i < charges.size(); ++i) {
284276
int index = int(charges(i)) - 1;
285277
atoms.push_back(radiiSet[index]);

‎src/interface/Meddle.cpp‎

Lines changed: 19 additions & 31 deletions
Original file line numberDiff line numberDiff line change
@@ -191,7 +191,6 @@ Printer hostWriter;
191191

192192
Meddle::Meddle(const std::string & inputFileName, const HostWriter & write)
193193
: input_(Input(inputFileName)), hasDynamic_(false) {
194-
bootstrap();
195194
input_.initMolecule();
196195
hostWriter.writer_ = write;
197196
infoStream_ << std::endl;
@@ -226,7 +225,6 @@ Meddle::Meddle(pcmsolver_reader_t input_reading,
226225
const PCMInput & host_input,
227226
const HostWriter & write)
228227
: hasDynamic_(false) {
229-
bootstrap();
230228
hostWriter.writer_ = write;
231229
TIMER_ON("Meddle::initInput");
232230
initInput(
@@ -454,8 +452,8 @@ void Meddle::initInput(pcmsolver_reader_t input_reading,
454452
}
455453

456454
void Meddle::initCavity() {
457-
cavity_ = Factory<ICavity, CavityData>::TheFactory().create(input_.cavityType(),
458-
input_.cavityParams());
455+
cavity_ =
456+
cavity::bootstrapFactory().create(input_.cavityType(), input_.cavityParams());
459457
cavity_->saveCavity();
460458

461459
infoStream_ << "========== Cavity " << std::endl;
@@ -464,16 +462,16 @@ void Meddle::initCavity() {
464462
}
465463

466464
void Meddle::initStaticSolver() {
467-
IGreensFunction * gf_i = Factory<IGreensFunction, GreenData>::TheFactory().create(
465+
IGreensFunction * gf_i = green::bootstrapFactory().create(
468466
input_.greenInsideType(), input_.insideGreenParams());
469-
IGreensFunction * gf_o = Factory<IGreensFunction, GreenData>::TheFactory().create(
467+
IGreensFunction * gf_o = green::bootstrapFactory().create(
470468
input_.greenOutsideType(), input_.outsideStaticGreenParams());
471-
std::string modelType = input_.solverType();
472-
K_0_ = Factory<ISolver, SolverData>::TheFactory().create(modelType,
473-
input_.solverParams());
474-
IBoundaryIntegralOperator * biop =
475-
Factory<IBoundaryIntegralOperator, BIOperatorData>::TheFactory().create(
476-
input_.integratorType(), input_.integratorParams());
469+
470+
K_0_ =
471+
solver::bootstrapFactory().create(input_.solverType(), input_.solverParams());
472+
473+
IBoundaryIntegralOperator * biop = bi_operators::bootstrapFactory().create(
474+
input_.integratorType(), input_.integratorParams());
477475
K_0_->buildSystemMatrix(*cavity_, *gf_i, *gf_o, *biop);
478476
delete biop;
479477

@@ -485,17 +483,16 @@ void Meddle::initStaticSolver() {
485483
}
486484

487485
void Meddle::initDynamicSolver() {
488-
IGreensFunction * gf_i = Factory<IGreensFunction, GreenData>::TheFactory().create(
486+
IGreensFunction * gf_i = green::bootstrapFactory().create(
489487
input_.greenInsideType(), input_.insideGreenParams());
490-
IGreensFunction * gf_o = Factory<IGreensFunction, GreenData>::TheFactory().create(
488+
IGreensFunction * gf_o = green::bootstrapFactory().create(
491489
input_.greenOutsideType(), input_.outsideDynamicGreenParams());
492-
std::string modelType = input_.solverType();
493-
K_d_ = Factory<ISolver, SolverData>::TheFactory().create(modelType,
494-
input_.solverParams());
495490

496-
IBoundaryIntegralOperator * biop =
497-
Factory<IBoundaryIntegralOperator, BIOperatorData>::TheFactory().create(
498-
input_.integratorType(), input_.integratorParams());
491+
K_d_ =
492+
solver::bootstrapFactory().create(input_.solverType(), input_.solverParams());
493+
494+
IBoundaryIntegralOperator * biop = bi_operators::bootstrapFactory().create(
495+
input_.integratorType(), input_.integratorParams());
499496
K_d_->buildSystemMatrix(*cavity_, *gf_i, *gf_o, *biop);
500497
hasDynamic_ = true;
501498
delete biop;
@@ -529,14 +526,6 @@ void Meddle::printInfo() const {
529526
hostWriter(infoStream_);
530527
}
531528

532-
void bootstrap() {
533-
bi_operators::bootstrapFactory();
534-
cavity::bootstrapFactory();
535-
green::bootstrapFactory();
536-
solver::bootstrapFactory();
537-
utils::bootstrapRadiiSet();
538-
}
539-
540529
void initMolecule(const Input & inp,
541530
const Symmetry & pg,
542531
int nuclei,
@@ -548,9 +537,8 @@ void initMolecule(const Input & inp,
548537
std::vector<Atom> radiiSet;
549538
std::vector<Atom> atoms;
550539
atoms.reserve(nuclei);
551-
// FIXME This is globally initialized...
552-
using utils::Factory;
553-
radiiSet = Factory::TheFactory().create(set);
540+
// FIXME Code duplication in function initMolecule in interface/Input.cpp
541+
tie(ignore, radiiSet) = utils::bootstrapRadiiSet().create(set);
554542
std::vector<Sphere> spheres;
555543
spheres.reserve(nuclei);
556544
for (int i = 0; i < charges.size(); ++i) {

‎src/solver/CPCMSolver.cpp‎

Lines changed: 5 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -29,6 +29,7 @@
2929

3030
#include "Config.hpp"
3131

32+
#include <Eigen/Cholesky>
3233
#include <Eigen/Core>
3334

3435
#include "bi_operators/IBoundaryIntegralOperator.hpp"
@@ -74,18 +75,19 @@ Eigen::VectorXd CPCMSolver::computeCharge_impl(const Eigen::VectorXd & potential
7475
int nrBlocks = blockS_.size();
7576
int irrDim = fullDim / nrBlocks;
7677
charge.segment(irrep * irrDim, irrDim) =
77-
-blockS_[irrep].llt().solve(potential.segment(irrep * irrDim, irrDim));
78+
-blockS_[irrep].ldlt().solve(potential.segment(irrep * irrDim, irrDim));
7879

7980
return charge;
8081
}
8182

8283
std::ostream & CPCMSolver::printSolver(std::ostream & os) {
8384
os << "Solver Type: C-PCM" << std::endl;
8485
if (hermitivitize_) {
85-
os << "PCM matrix hermitivitized";
86+
os << "PCM matrix hermitivitized" << std::endl;
8687
} else {
87-
os << "PCM matrix NOT hermitivitized (matches old DALTON)";
88+
os << "PCM matrix NOT hermitivitized (matches old DALTON)" << std::endl;
8889
}
90+
os << "Correction = " << correction_;
8991

9092
return os;
9193
}

‎src/solver/Solver.hpp‎

Lines changed: 7 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -52,18 +52,20 @@ inline ISolver * createIEFSolver(const SolverData & data) {
5252
return new IEFSolver(data.hermitivitize);
5353
}
5454

55-
inline void bootstrapFactory() {
55+
inline Factory<ISolver, SolverData> bootstrapFactory() {
56+
Factory<ISolver, SolverData> factory_;
57+
5658
const bool registeredCPCMSolver =
57-
Factory<ISolver, SolverData>::TheFactory().registerObject("CPCM",
58-
createCPCMSolver);
59+
factory_.registerObject("CPCM", createCPCMSolver);
5960
if (!registeredCPCMSolver)
6061
PCMSOLVER_ERROR("Subscription of solver for CPCM to factory failed!");
6162

6263
const bool registeredIEFSolver =
63-
Factory<ISolver, SolverData>::TheFactory().registerObject("IEFPCM",
64-
createIEFSolver);
64+
factory_.registerObject("IEFPCM", createIEFSolver);
6565
if (!registeredIEFSolver)
6666
PCMSOLVER_ERROR("Subscription of solver for IEF to factory failed!");
67+
68+
return factory_;
6769
}
6870
} // namespace solver
6971
} // namespace pcm

‎src/utils/Atom.cpp‎

Lines changed: 9 additions & 10 deletions
Original file line numberDiff line numberDiff line change
@@ -38,9 +38,8 @@ bool invalid(const Atom & atom) {
3838

3939
namespace utils {
4040
namespace detail {
41-
std::vector<Atom> & initBondi() {
42-
43-
static std::vector<Atom> Bondi(118);
41+
RadiiSet initBondi() {
42+
std::vector<Atom> Bondi(118);
4443
Eigen::Vector3d Origin = Eigen::Vector3d::Zero();
4544

4645
// ------------------------------------------------------------
@@ -167,12 +166,12 @@ std::vector<Atom> & initBondi() {
167166

168167
// ------------------------------------------------------------
169168

170-
return Bondi;
169+
return make_tuple("Bondi-Mantina", Bondi);
171170
}
172171

173-
std::vector<Atom> & initUFF() {
172+
RadiiSet initUFF() {
174173

175-
static std::vector<Atom> UFF(118);
174+
std::vector<Atom> UFF(118);
176175
Eigen::Vector3d Origin = Eigen::Vector3d::Zero();
177176

178177
// ------------------------------------------------------------
@@ -298,12 +297,12 @@ std::vector<Atom> & initUFF() {
298297
UFF[117] = Atom("Ununoctium", "Uuo", 118.0, 294.0000000, 0.0000, Origin, 1.0);
299298
// ------------------------------------------------------------
300299

301-
return UFF;
300+
return make_tuple("UFF", UFF);
302301
}
303302

304-
std::vector<Atom> & initAllinger() {
303+
RadiiSet initAllinger() {
305304

306-
static std::vector<Atom> Allinger(118);
305+
std::vector<Atom> Allinger(118);
307306
Eigen::Vector3d Origin = Eigen::Vector3d::Zero();
308307

309308
// ------------------------------------------------------------
@@ -447,7 +446,7 @@ std::vector<Atom> & initAllinger() {
447446
Atom("Ununoctium", "Uuo", 118.0, 294.0000000, 0.000000, Origin, 1.0);
448447
// ------------------------------------------------------------
449448

450-
return Allinger;
449+
return make_tuple("Allinger's MM3", Allinger);
451450
}
452451
} // namespace detail
453452
} // namespace utils

0 commit comments

Comments
 (0)
Back | FazBrowse Home | New Git URL