| 46 |  | #include "minimizers/Minimizer.hpp" | 
| 47 |  | #include "primitives/Molecule.hpp" | 
| 48 |  | namespace oopse { | 
| 49 | < | double dotProduct(const std::vector<double>& v1, const std::vector<double>& v2) { | 
| 49 | > | RealType dotProduct(const std::vector<RealType>& v1, const std::vector<RealType>& v2) { | 
| 50 |  | if (v1.size() != v2.size()) { | 
| 51 |  |  | 
| 52 |  | } | 
| 53 |  |  | 
| 54 |  |  | 
| 55 | < | double result = 0.0; | 
| 55 | > | RealType result = 0.0; | 
| 56 |  | for (unsigned int i = 0; i < v1.size(); ++i) { | 
| 57 |  | result += v1[i] * v2[i]; | 
| 58 |  | } | 
| 76 |  | delete paramSet; | 
| 77 |  | } | 
| 78 |  |  | 
| 79 | < | void Minimizer::calcEnergyGradient(std::vector<double> &x, | 
| 80 | < | std::vector<double> &grad, double&energy, int&status) { | 
| 79 | > | void Minimizer::calcEnergyGradient(std::vector<RealType> &x, | 
| 80 | > | std::vector<RealType> &grad, RealType&energy, int&status) { | 
| 81 |  |  | 
| 82 |  | SimInfo::MoleculeIterator i; | 
| 83 |  | Molecule::IntegrableObjectIterator  j; | 
| 84 |  | Molecule* mol; | 
| 85 |  | StuntDouble* integrableObject; | 
| 86 | < | std::vector<double> myGrad; | 
| 86 | > | std::vector<RealType> myGrad; | 
| 87 |  | int shakeStatus; | 
| 88 |  |  | 
| 89 |  | status = 1; | 
| 110 |  |  | 
| 111 |  | myGrad = integrableObject->getGrad(); | 
| 112 |  | for (unsigned int k = 0; k < myGrad.size(); ++k) { | 
| 113 | < | //gradient is equal to -f | 
| 114 | < | grad[index++] = -myGrad[k]; | 
| 113 | > |  | 
| 114 | > | grad[index++] = myGrad[k]; | 
| 115 |  | } | 
| 116 |  | } | 
| 117 |  | } | 
| 118 |  |  | 
| 119 |  | } | 
| 120 |  |  | 
| 121 | < | void Minimizer::setCoor(std::vector<double> &x) { | 
| 121 | > | void Minimizer::setCoor(std::vector<RealType> &x) { | 
| 122 |  | Vector3d position; | 
| 123 |  | Vector3d eulerAngle; | 
| 124 |  | SimInfo::MoleculeIterator i; | 
| 149 |  |  | 
| 150 |  | } | 
| 151 |  |  | 
| 152 | < | std::vector<double> Minimizer::getCoor() { | 
| 152 | > | std::vector<RealType> Minimizer::getCoor() { | 
| 153 |  | Vector3d position; | 
| 154 |  | Vector3d eulerAngle; | 
| 155 |  | SimInfo::MoleculeIterator i; | 
| 157 |  | Molecule* mol; | 
| 158 |  | StuntDouble* integrableObject; | 
| 159 |  | int index = 0; | 
| 160 | < | std::vector<double> x(getDim()); | 
| 160 | > | std::vector<RealType> x(getDim()); | 
| 161 |  |  | 
| 162 |  | for (mol = info->beginMolecule(i); mol != NULL; mol = info->nextMolecule(i)) { | 
| 163 |  | for (integrableObject = mol->beginIntegrableObject(j); integrableObject != NULL; | 
| 186 |  |  | 
| 187 |  | int    done; | 
| 188 |  |  | 
| 189 | < | double posA[3], posB[3]; | 
| 189 | > | RealType posA[3], posB[3]; | 
| 190 |  |  | 
| 191 | < | double velA[3], velB[3]; | 
| 191 | > | RealType velA[3], velB[3]; | 
| 192 |  |  | 
| 193 | < | double pab[3]; | 
| 193 | > | RealType pab[3]; | 
| 194 |  |  | 
| 195 | < | double rab[3]; | 
| 195 | > | RealType rab[3]; | 
| 196 |  |  | 
| 197 |  | int    a,       b, | 
| 198 |  | ax,      ay, | 
| 199 |  | az,      bx, | 
| 200 |  | by,      bz; | 
| 201 |  |  | 
| 202 | < | double rma,     rmb; | 
| 202 | > | RealType rma,     rmb; | 
| 203 |  |  | 
| 204 | < | double dx,      dy, | 
| 204 | > | RealType dx,      dy, | 
| 205 |  | dz; | 
| 206 |  |  | 
| 207 | < | double rpab; | 
| 207 | > | RealType rpab; | 
| 208 |  |  | 
| 209 | < | double rabsq,   pabsq, | 
| 209 | > | RealType rabsq,   pabsq, | 
| 210 |  | rpabsq; | 
| 211 |  |  | 
| 212 | < | double diffsq; | 
| 212 | > | RealType diffsq; | 
| 213 |  |  | 
| 214 | < | double gab; | 
| 214 | > | RealType gab; | 
| 215 |  |  | 
| 216 |  | int    iteration; | 
| 217 |  |  | 
| 377 |  |  | 
| 378 |  | int    done; | 
| 379 |  |  | 
| 380 | < | double posA[3], posB[3]; | 
| 380 | > | RealType posA[3], posB[3]; | 
| 381 |  |  | 
| 382 | < | double frcA[3], frcB[3]; | 
| 382 | > | RealType frcA[3], frcB[3]; | 
| 383 |  |  | 
| 384 | < | double rab[3],  fpab[3]; | 
| 384 | > | RealType rab[3],  fpab[3]; | 
| 385 |  |  | 
| 386 |  | int    a,       b, | 
| 387 |  | ax,      ay, | 
| 388 |  | az,      bx, | 
| 389 |  | by,      bz; | 
| 390 |  |  | 
| 391 | < | double rma,     rmb; | 
| 391 | > | RealType rma,     rmb; | 
| 392 |  |  | 
| 393 | < | double rvab; | 
| 393 | > | RealType rvab; | 
| 394 |  |  | 
| 395 | < | double gab; | 
| 395 | > | RealType gab; | 
| 396 |  |  | 
| 397 | < | double rabsq; | 
| 397 | > | RealType rabsq; | 
| 398 |  |  | 
| 399 | < | double rfab; | 
| 399 | > | RealType rfab; | 
| 400 |  |  | 
| 401 |  | int    iteration; | 
| 402 |  |  | 
| 524 |  | calcEnergyGradient(curX, curG, curF, egEvalStatus); | 
| 525 |  | } | 
| 526 |  |  | 
| 527 | < | void Minimizer::calcF(std::vector < double > &x, double&f, int&status) { | 
| 528 | < | std::vector < double > tempG; | 
| 527 | > | void Minimizer::calcF(std::vector < RealType > &x, RealType&f, int&status) { | 
| 528 | > | std::vector < RealType > tempG; | 
| 529 |  |  | 
| 530 |  | tempG.resize(x.size()); | 
| 531 |  |  | 
| 538 |  | calcEnergyGradient(curX, curG, curF, egEvalStatus); | 
| 539 |  | } | 
| 540 |  |  | 
| 541 | < | void Minimizer::calcG(std::vector<double>& x, std::vector<double>& g, double&f, int&status) { | 
| 541 | > | void Minimizer::calcG(std::vector<RealType>& x, std::vector<RealType>& g, RealType&f, int&status) { | 
| 542 |  | calcEnergyGradient(x, g, f, status); | 
| 543 |  | } | 
| 544 |  |  | 
| 564 |  | } | 
| 565 |  | } | 
| 566 |  |  | 
| 567 | < | void Minimizer::setX(std::vector < double > &x) { | 
| 567 | > | void Minimizer::setX(std::vector < RealType > &x) { | 
| 568 |  | if (x.size() != ndim) { | 
| 569 |  | sprintf(painCave.errMsg, "Minimizer Error: dimesion of x and curX does not match\n"); | 
| 570 |  | painCave.isFatal = 1; | 
| 574 |  | curX = x; | 
| 575 |  | } | 
| 576 |  |  | 
| 577 | < | void Minimizer::setG(std::vector < double > &g) { | 
| 577 | > | void Minimizer::setG(std::vector < RealType > &g) { | 
| 578 |  | if (g.size() != ndim) { | 
| 579 |  | sprintf(painCave.errMsg, "Minimizer Error: dimesion of g and curG does not match\n"); | 
| 580 |  | painCave.isFatal = 1; | 
| 596 |  | * @todo optimize this line search algorithm | 
| 597 |  | */ | 
| 598 |  |  | 
| 599 | < | int Minimizer::doLineSearch(std::vector<double> &direction, | 
| 600 | < | double stepSize) { | 
| 599 | > | int Minimizer::doLineSearch(std::vector<RealType> &direction, | 
| 600 | > | RealType stepSize) { | 
| 601 |  |  | 
| 602 | < | std::vector<double> xa; | 
| 603 | < | std::vector<double> xb; | 
| 604 | < | std::vector<double> xc; | 
| 605 | < | std::vector<double> ga; | 
| 606 | < | std::vector<double> gb; | 
| 607 | < | std::vector<double> gc; | 
| 608 | < | double fa; | 
| 609 | < | double fb; | 
| 610 | < | double fc; | 
| 611 | < | double a; | 
| 612 | < | double b; | 
| 613 | < | double c; | 
| 602 | > | std::vector<RealType> xa; | 
| 603 | > | std::vector<RealType> xb; | 
| 604 | > | std::vector<RealType> xc; | 
| 605 | > | std::vector<RealType> ga; | 
| 606 | > | std::vector<RealType> gb; | 
| 607 | > | std::vector<RealType> gc; | 
| 608 | > | RealType fa; | 
| 609 | > | RealType fb; | 
| 610 | > | RealType fc; | 
| 611 | > | RealType a; | 
| 612 | > | RealType b; | 
| 613 | > | RealType c; | 
| 614 |  | int    status; | 
| 615 | < | double initSlope; | 
| 616 | < | double slopeA; | 
| 617 | < | double slopeB; | 
| 618 | < | double slopeC; | 
| 615 | > | RealType initSlope; | 
| 616 | > | RealType slopeA; | 
| 617 | > | RealType slopeB; | 
| 618 | > | RealType slopeC; | 
| 619 |  | bool   foundLower; | 
| 620 |  | int    iter; | 
| 621 |  | int    maxLSIter; | 
| 622 | < | double mu; | 
| 623 | < | double eta; | 
| 624 | < | double ftol; | 
| 625 | < | double lsTol; | 
| 622 | > | RealType mu; | 
| 623 | > | RealType eta; | 
| 624 | > | RealType ftol; | 
| 625 | > | RealType lsTol; | 
| 626 |  |  | 
| 627 |  | xa.resize(ndim); | 
| 628 |  | xb.resize(ndim); | 
| 873 |  | if (curIter == nextWriteIter) { | 
| 874 |  | nextWriteIter += writeFrq; | 
| 875 |  | calcF(); | 
| 876 | < | dumpWriter.writeDump(); | 
| 876 | > | dumpWriter.writeDumpAndEor(); | 
| 877 |  | statWriter.writeStat(curSnapshot->statData); | 
| 878 |  | } | 
| 879 |  |  | 
| 902 |  | } | 
| 903 |  |  | 
| 904 |  |  | 
| 905 | < | double Minimizer::calcPotential() { | 
| 905 | > | RealType Minimizer::calcPotential() { | 
| 906 |  | forceMan->calcForces(true, false); | 
| 907 |  |  | 
| 908 |  | Snapshot* curSnapshot = info->getSnapshotManager()->getCurrentSnapshot(); | 
| 909 | < | double potential_local = curSnapshot->statData[Stats::LONG_RANGE_POTENTIAL] + | 
| 909 | > | RealType potential_local = curSnapshot->statData[Stats::LONG_RANGE_POTENTIAL] + | 
| 910 |  | curSnapshot->statData[Stats::SHORT_RANGE_POTENTIAL] ; | 
| 911 | < | double potential; | 
| 911 | > | RealType potential; | 
| 912 |  |  | 
| 913 |  | #ifdef IS_MPI | 
| 914 | < | MPI_Allreduce(&potential_local, &potential, 1, MPI_DOUBLE, MPI_SUM, | 
| 914 | > | MPI_Allreduce(&potential_local, &potential, 1, MPI_REALTYPE, MPI_SUM, | 
| 915 |  | MPI_COMM_WORLD); | 
| 916 |  | #else | 
| 917 |  | potential = potential_local; |