Changeset 300220
- Timestamp:
- Mar 1, 2011, 1:17:08 PM (14 years ago)
- Branches:
- Action_Thermostats, Add_AtomRandomPerturbation, Add_FitFragmentPartialChargesAction, Add_RotateAroundBondAction, Add_SelectAtomByNameAction, Added_ParseSaveFragmentResults, AddingActions_SaveParseParticleParameters, Adding_Graph_to_ChangeBondActions, Adding_MD_integration_tests, Adding_ParticleName_to_Atom, Adding_StructOpt_integration_tests, AtomFragments, Automaking_mpqc_open, AutomationFragmentation_failures, Candidate_v1.5.4, Candidate_v1.6.0, Candidate_v1.6.1, ChangeBugEmailaddress, ChangingTestPorts, ChemicalSpaceEvaluator, CombiningParticlePotentialParsing, Combining_Subpackages, Debian_Package_split, Debian_package_split_molecuildergui_only, Disabling_MemDebug, Docu_Python_wait, EmpiricalPotential_contain_HomologyGraph, EmpiricalPotential_contain_HomologyGraph_documentation, Enable_parallel_make_install, Enhance_userguide, Enhanced_StructuralOptimization, Enhanced_StructuralOptimization_continued, Example_ManyWaysToTranslateAtom, Exclude_Hydrogens_annealWithBondGraph, FitPartialCharges_GlobalError, Fix_BoundInBox_CenterInBox_MoleculeActions, Fix_ChargeSampling_PBC, Fix_ChronosMutex, Fix_FitPartialCharges, Fix_FitPotential_needs_atomicnumbers, Fix_ForceAnnealing, Fix_IndependentFragmentGrids, Fix_ParseParticles, Fix_ParseParticles_split_forward_backward_Actions, Fix_PopActions, Fix_QtFragmentList_sorted_selection, Fix_Restrictedkeyset_FragmentMolecule, Fix_StatusMsg, Fix_StepWorldTime_single_argument, Fix_Verbose_Codepatterns, Fix_fitting_potentials, Fixes, ForceAnnealing_goodresults, ForceAnnealing_oldresults, ForceAnnealing_tocheck, ForceAnnealing_with_BondGraph, ForceAnnealing_with_BondGraph_continued, ForceAnnealing_with_BondGraph_continued_betteresults, ForceAnnealing_with_BondGraph_contraction-expansion, FragmentAction_writes_AtomFragments, FragmentMolecule_checks_bonddegrees, GeometryObjects, Gui_Fixes, Gui_displays_atomic_force_velocity, ImplicitCharges, IndependentFragmentGrids, IndependentFragmentGrids_IndividualZeroInstances, IndependentFragmentGrids_IntegrationTest, IndependentFragmentGrids_Sole_NN_Calculation, JobMarket_RobustOnKillsSegFaults, JobMarket_StableWorkerPool, JobMarket_unresolvable_hostname_fix, MoreRobust_FragmentAutomation, ODR_violation_mpqc_open, PartialCharges_OrthogonalSummation, PdbParser_setsAtomName, PythonUI_with_named_parameters, QtGui_reactivate_TimeChanged_changes, Recreated_GuiChecks, Rewrite_FitPartialCharges, RotateToPrincipalAxisSystem_UndoRedo, SaturateAtoms_findBestMatching, SaturateAtoms_singleDegree, StoppableMakroAction, Subpackage_CodePatterns, Subpackage_JobMarket, Subpackage_LinearAlgebra, Subpackage_levmar, Subpackage_mpqc_open, Subpackage_vmg, Switchable_LogView, ThirdParty_MPQC_rebuilt_buildsystem, TrajectoryDependenant_MaxOrder, TremoloParser_IncreasedPrecision, TremoloParser_MultipleTimesteps, TremoloParser_setsAtomName, Ubuntu_1604_changes, stable
- Children:
- 0ec7fe
- Parents:
- 0cbad2
- git-author:
- Frederik Heber <heber@…> (02/28/11 15:00:34)
- git-committer:
- Frederik Heber <heber@…> (03/01/11 13:17:08)
- Files:
-
- 7 edited
Legend:
- Unmodified
- Added
- Removed
-
configure.ac
r0cbad2 r300220 115 115 #AC_MSG_NOTICE(["GLU_CFLAGS: $GLU_CFLAGS, GLU_CXXFLAGS: $GLU_CXXFLAGS, GLU_LDFLAGS: $GLU_LDFLAGS, GLU_LIBS: $GLU_LIBS"]) 116 116 117 # CodePatterns library (needs ObservedIterator<> operator==,!= to be const)118 AM_PATH_CODEPATTERNS([1.0.1 1], $have_debug,,[AC_MSG_ERROR([Missing CodePatterns library, please specify PKG_CONFIG_PATH or use --with-codepatterns...])])117 # CodePatterns library (needs operator<<(.., range<>) ) 118 AM_PATH_CODEPATTERNS([1.0.13], $have_debug,,[AC_MSG_ERROR([Missing CodePatterns library, please specify PKG_CONFIG_PATH or use --with-codepatterns...])]) 119 119 120 120 # Checks for header files. -
src/Actions/FragmentationAction/CreateAdjacencyAction.cpp
r0cbad2 r300220 27 27 #include "bondgraph.hpp" 28 28 #include "CodePatterns/Log.hpp" 29 #include "CodePatterns/Range.hpp" 29 30 #include "CodePatterns/Verbose.hpp" 30 31 #include "config.hpp" … … 65 66 atom *OtherWalker = NULL; 66 67 int n[NDIM]; 67 double MinDistance, MaxDistance;68 68 LinkedCell *LC = NULL; 69 69 Box &domain = World::getInstance().getDomain(); … … 136 136 ASSERT(OtherWalker,"TesselPoint that was not an atom retrieved from LinkedNode"); 137 137 //Log() << Verbose(1) << "Checking distance " << OtherWalker->x.PeriodicDistanceSquared(&(Walker->x), cell_size) << " against typical bond length of " << bonddistance*bonddistance << "." << endl; 138 BG->getMinMaxDistance(Walker, OtherWalker, MinDistance, MaxDistance, IsAngstroem); 138 range<double> MinMaxDistanceSquared(0.,0.); 139 BG->getMinMaxDistanceSquared(Walker, OtherWalker, MinMaxDistanceSquared, IsAngstroem); 139 140 const double distance = domain.periodicDistanceSquared(OtherWalker->getPosition(),Walker->getPosition()); 140 const bool status = (distance <= MaxDistance * MaxDistance) && (distance >= MinDistance * MinDistance);141 // L og() << Verbose(1) << "MinDistance is " << MinDistance << " and MaxDistance is " << MaxDistance << "." << endl;141 const bool status = MinMaxDistanceSquared.isInRange(distance); 142 // LOG3, "INFO: MinMaxDistanceSquared interval is " << MinMaxDistanceSquared << "."); 142 143 if (AtomIds[OtherWalker->father] > AtomIds[Walker->father]) { 143 144 if (status) { // create bond if distance is smaller -
src/bondgraph.cpp
r0cbad2 r300220 93 93 } 94 94 95 double BondGraph::GetBondLength(int firstZ, int secondZ) 95 double BondGraph::GetBondLength( 96 int firstZ, 97 int secondZ) const 96 98 { 97 99 std::cout << "Request for length between " << firstZ << " and " << secondZ << ": "; … … 110 112 } 111 113 112 113 void BondGraph::CovalentMinMaxDistance(const BondedParticle * const Walker, const BondedParticle * const OtherWalker, double &MinDistance, double &MaxDistance, bool IsAngstroem) 114 { 115 MinDistance = OtherWalker->getType()->getCovalentRadius() + Walker->getType()->getCovalentRadius(); 116 MinDistance *= (IsAngstroem) ? 1. : 1. / AtomicLengthToAngstroem; 117 MaxDistance = MinDistance + BondThreshold; 118 MinDistance -= BondThreshold; 119 } 120 121 void BondGraph::BondLengthMatrixMinMaxDistance(const BondedParticle * const Walker, const BondedParticle * const OtherWalker, double &MinDistance, double &MaxDistance, bool IsAngstroem) 122 { 123 ASSERT(BondLengthMatrix != NULL, 124 "BondGraph::BondLengthMatrixMinMaxDistance() called without NULL BondLengthMatrix."); 125 MinDistance = GetBondLength(Walker->getType()->getAtomicNumber()-1, OtherWalker->getType()->getAtomicNumber()-1); 126 MinDistance *= (IsAngstroem) ? 1. : 1. / AtomicLengthToAngstroem; 127 MaxDistance = MinDistance + BondThreshold; 128 MinDistance -= BondThreshold; 129 } 130 131 void BondGraph::getMinMaxDistance(const BondedParticle * const Walker, const BondedParticle * const OtherWalker, double &MinDistance, double &MaxDistance, bool IsAngstroem) 114 void BondGraph::CovalentMinMaxDistance( 115 const element * const Walker, 116 const element * const OtherWalker, 117 range<double> &MinMaxDistance, 118 bool IsAngstroem) const 119 { 120 MinMaxDistance.first = OtherWalker->getCovalentRadius() + Walker->getCovalentRadius(); 121 MinMaxDistance.first *= (IsAngstroem) ? 1. : 1. / AtomicLengthToAngstroem; 122 MinMaxDistance.last = MinMaxDistance.first + BondThreshold; 123 MinMaxDistance.first -= BondThreshold; 124 } 125 126 void BondGraph::BondLengthMatrixMinMaxDistance( 127 const element * const Walker, 128 const element * const OtherWalker, 129 range<double> &MinMaxDistance, 130 bool IsAngstroem) const 131 { 132 ASSERT(BondLengthMatrix, "BondGraph::BondLengthMatrixMinMaxDistance() called without NULL BondLengthMatrix."); 133 ASSERT(Walker, "BondGraph::BondLengthMatrixMinMaxDistance() - illegal element given."); 134 ASSERT(OtherWalker, "BondGraph::BondLengthMatrixMinMaxDistance() - illegal other element given."); 135 MinMaxDistance.first = GetBondLength(Walker->getAtomicNumber()-1, OtherWalker->getAtomicNumber()-1); 136 MinMaxDistance.first *= (IsAngstroem) ? 1. : 1. / AtomicLengthToAngstroem; 137 MinMaxDistance.last = MinMaxDistance.first + BondThreshold; 138 MinMaxDistance.first -= BondThreshold; 139 } 140 141 void BondGraph::getMinMaxDistance( 142 const element * const Walker, 143 const element * const OtherWalker, 144 range<double> &MinMaxDistance, 145 bool IsAngstroem) const 132 146 { 133 147 if (BondLengthMatrix == NULL) {// safety measure if no matrix has been parsed yet 134 LOG(2, "INFO: Using Covalent radii criterion for [min,max ]distances.");135 CovalentMinMaxDistance(Walker, OtherWalker, Min Distance,MaxDistance, IsAngstroem);148 LOG(2, "INFO: Using Covalent radii criterion for [min,max) distances."); 149 CovalentMinMaxDistance(Walker, OtherWalker, MinMaxDistance, IsAngstroem); 136 150 } else { 137 LOG(2, "INFO: Using Covalent radii criterion for [min,max] distances."); 138 BondLengthMatrixMinMaxDistance(Walker, OtherWalker, MinDistance, MaxDistance, IsAngstroem); 139 } 151 LOG(2, "INFO: Using Covalent radii criterion for [min,max) distances."); 152 BondLengthMatrixMinMaxDistance(Walker, OtherWalker, MinMaxDistance, IsAngstroem); 153 } 154 } 155 156 void BondGraph::getMinMaxDistance( 157 const BondedParticle * const Walker, 158 const BondedParticle * const OtherWalker, 159 range<double> &MinMaxDistance, 160 bool IsAngstroem) const 161 { 162 getMinMaxDistance(Walker->getType(), OtherWalker->getType(),MinMaxDistance, IsAngstroem); 163 } 164 165 void BondGraph::getMinMaxDistanceSquared( 166 const BondedParticle * const Walker, 167 const BondedParticle * const OtherWalker, 168 range<double> &MinMaxDistance, 169 bool IsAngstroem) const 170 { 171 // use non-squared version 172 getMinMaxDistance(Walker, OtherWalker,MinMaxDistance, IsAngstroem); 173 // and square 174 MinMaxDistance.first *= MinMaxDistance.first; 175 MinMaxDistance.last *= MinMaxDistance.last; 140 176 } 141 177 … … 145 181 atom *OtherWalker = NULL; 146 182 int n[NDIM]; 147 double MinDistance, MaxDistance;148 183 Box &domain = World::getInstance().getDomain(); 149 184 … … 174 209 ASSERT(OtherWalker != NULL, 175 210 "BondGraph::CreateAdjacency() - TesselPoint that was not an atom retrieved from LinkedNode"); 176 getMinMaxDistance(Walker, OtherWalker, MinDistance, MaxDistance, IsAngstroem); 211 range<double> MinMaxDistanceSquared(0.,0.); 212 getMinMaxDistanceSquared(Walker, OtherWalker, MinMaxDistanceSquared, IsAngstroem); 177 213 const double distance = domain.periodicDistanceSquared(OtherWalker->getPosition(),Walker->getPosition()); 178 LOG(2, "INFO: Checking distance " << distance << " against typical bond length of [" << MinDistance * MinDistance << "," << MaxDistance * MaxDistance << "]."); 179 const bool status = (distance <= MaxDistance * MaxDistance) && (distance >= MinDistance * MinDistance); 180 LOG(3, "INFO: MinDistance is " << MinDistance << " and MaxDistance is " << MaxDistance << "."); 214 LOG(2, "INFO: Checking squared distance " << distance << " against typical bond length of " << MinMaxDistanceSquared << "."); 215 const bool status = MinMaxDistanceSquared.isInRange(distance); 181 216 if (OtherWalker->father > Walker->father ) { // just to not add bonds from both sides 182 217 if (status) { // create bond if distance is smaller … … 188 223 } else { 189 224 LOG(1, "REJECT: Squared distance " 190 << distance << " is out of covalent bounds [" 191 << MinDistance*MinDistance << "," 192 << MaxDistance * MaxDistance << "]."); 225 << distance << " is out of squared covalent bounds " 226 << MinMaxDistanceSquared << "."); 193 227 } 194 228 } else { -
src/bondgraph.hpp
r0cbad2 r300220 44 44 45 45 class BondGraph { 46 //!> analysis bonds unit test should be friend to access private parts. 47 friend class AnalysisBondsTest; 48 //!> own bond graph unit test should be friend to access private parts. 49 friend class BondGraphTest; 46 50 public: 47 51 /** Constructor of class BondGraph. … … 62 66 bool LoadBondLengthTable(std::istream &input); 63 67 64 /** Returns the entry for a given index pair.65 * \param firstelement index/atom number of first element (row index)66 * \param secondelement index/atom number of second element (column index)67 * \note matrix is of course symmetric.68 */69 double GetBondLength(int firstelement, int secondelement);70 71 68 /** Determines the maximum of all element::CovalentRadius for elements present in \a &Set. 72 69 * … … 103 100 104 101 /** Returns bond criterion for given pair based on a bond length matrix. 105 * This calls e ither the covalent or the bond matrix criterion.102 * This calls element-version of getMinMaxDistance(). 106 103 * \param *Walker first BondedParticle 107 104 * \param *OtherWalker second BondedParticle 108 * \param &MinDistance lower bond bound on return 109 * \param &MaxDistance upper bond bound on return 105 * \param &MinMaxDistance Range for interval on return 110 106 * \param IsAngstroem whether units are in angstroem or bohr radii 111 107 */ 112 void getMinMaxDistance(const BondedParticle * const Walker, const BondedParticle * const OtherWalker, double &MinDistance, double &MaxDistance, bool IsAngstroem); 108 void getMinMaxDistance( 109 const BondedParticle * const Walker, 110 const BondedParticle * const OtherWalker, 111 range<double> &MinMaxDistance, 112 bool IsAngstroem) const; 113 114 /** Returns SQUARED bond criterion for given pair based on a bond length matrix. 115 * This calls element-version of getMinMaxDistance() and squares the values 116 * of either interval end. 117 * \param *Walker first BondedParticle 118 * \param *OtherWalker second BondedParticle 119 * \param &MinMaxDistance Range for interval on return 120 * \param IsAngstroem whether units are in angstroem or bohr radii 121 */ 122 void getMinMaxDistanceSquared( 123 const BondedParticle * const Walker, 124 const BondedParticle * const OtherWalker, 125 range<double> &MinMaxDistance, 126 bool IsAngstroem) const; 113 127 114 128 /** Creates the adjacency list for a given \a Range of iterable atoms. … … 291 305 static const double BondThreshold; 292 306 307 /** Returns the BondLengthMatrix entry for a given index pair. 308 * \param firstelement index/atom number of first element (row index) 309 * \param secondelement index/atom number of second element (column index) 310 * \note matrix is of course symmetric. 311 */ 312 double GetBondLength( 313 int firstelement, 314 int secondelement) const; 315 293 316 /** Returns bond criterion for given pair based on a bond length matrix. 294 317 * The matrix should be contained in \a this BondGraph and contain an element- … … 296 319 * \param *Walker first BondedParticle 297 320 * \param *OtherWalker second BondedParticle 298 * \param &MinDistance lower bond bound on return 299 * \param &MaxDistance upper bond bound on return 321 * \param &MinMaxDistance Range for interval on return 300 322 * \param IsAngstroem whether units are in angstroem or bohr radii 301 323 */ 302 void BondLengthMatrixMinMaxDistance(const BondedParticle * const Walker, const BondedParticle * const OtherWalker, double &MinDistance, double &MaxDistance, bool IsAngstroem); 303 304 /** Returns bond criterion for given pair based on covalent radius. 305 * \param *Walker first BondedParticle 306 * \param *OtherWalker second BondedParticle 307 * \param &MinDistance lower bond bound on return 308 * \param &MaxDistance upper bond bound on return 324 void getMinMaxDistance( 325 const element * const Walker, 326 const element * const OtherWalker, 327 range<double> &MinMaxDistance, 328 bool IsAngstroem) const; 329 330 /** Returns bond criterion for given pair of elements based on a bond length matrix. 331 * The matrix should be contained in \a this BondGraph and contain an element- 332 * to-element length. 333 * \param *Walker first element 334 * \param *OtherWalker second element 335 * @param MinMaxDistance reference to range type set on return 309 336 * \param IsAngstroem whether units are in angstroem or bohr radii 310 337 */ 311 void CovalentMinMaxDistance(const BondedParticle * const Walker, const BondedParticle * const OtherWalker, double &MinDistance, double &MaxDistance, bool IsAngstroem); 338 void BondLengthMatrixMinMaxDistance( 339 const element * const Walker, 340 const element * const OtherWalker, 341 range<double> &DistanceInterval, 342 bool IsAngstroem) const; 343 344 /** Returns bond criterion for given pair of elements based on covalent radius. 345 * \param *Walker first element 346 * \param *OtherWalker second element 347 * @param MinMaxDistance reference to range type set on return 348 * \param IsAngstroem whether units are in angstroem or bohr radii 349 */ 350 void CovalentMinMaxDistance( 351 const element * const Walker, 352 const element * const OtherWalker, 353 range<double> &DistanceInterval, 354 bool IsAngstroem) const; 312 355 313 356 -
src/molecule.cpp
r0cbad2 r300220 344 344 // the problem is not the H being out of the box, but InBondvector have the wrong direction 345 345 // due to TopReplacement or Origin being on the wrong side! 346 double MinBondDistance; 347 double MaxBondDistance; 348 World::getInstance().getBondGraph()->getMinMaxDistance( 349 TopOrigin, 350 TopReplacement, 351 MinBondDistance,MaxBondDistance,IsAngstroem); 352 std::cout << "bondlength > MaxBondDistance: " << bondlength << " > " << MaxBondDistance << std::endl; 353 std::cout << "bondlength < MinBondDistance: " << bondlength << " < " << MinBondDistance << std::endl; 354 if ((bondlength < MinBondDistance) || (bondlength > MaxBondDistance)) { 346 range<double> MinMaxBondDistance(0.,0.); 347 const BondGraph * const BG = World::getInstance().getBondGraph(); 348 BG->getMinMaxDistance(TopOrigin,TopReplacement,MinMaxBondDistance,IsAngstroem); 349 if (!MinMaxBondDistance.isInRange(bondlength)) { 355 350 // Log() << Verbose(4) << "InBondvector is: "; 356 351 // InBondvector.Output(out); … … 359 354 for (int i=NDIM;i--;) { 360 355 l = TopReplacement->at(i) - TopOrigin->at(i); 361 if (fabs(l) > M axBondDistance) { // is component greater than bond distance (check against min not useful here)356 if (fabs(l) > MinMaxBondDistance.last) { // is component greater than bond distance (check against min not useful here) 362 357 Orthovector1[i] = (l < 0) ? -1. : +1.; 363 358 } // (signs are correct, was tested!) -
src/molecule_fragmentation.cpp
r0cbad2 r300220 1811 1811 tmp = fabs(Binder->leftatom->at(i) - Binder->rightatom->at(i)); 1812 1812 //Log() << Verbose(3) << "Checking " << i << "th distance of " << *Binder->leftatom << " to " << *Binder->rightatom << ": " << tmp << "." << endl; 1813 double MinBondDistance; 1814 double MaxBondDistance; 1815 BG->getMinMaxDistance(Binder->leftatom, Binder->rightatom, MinBondDistance, MaxBondDistance, true); 1816 if ((tmp < MinBondDistance) || ((tmp > MaxBondDistance))) { 1813 range<double> MinMaxDistance(0.,0.); 1814 BG->getMinMaxDistance(Binder->leftatom, Binder->rightatom, MinMaxDistance, true); 1815 if (!MinMaxDistance.isInRange(tmp)) { 1817 1816 DoLog(2) && (Log() << Verbose(2) << "Correcting at bond " << *Binder << "." << endl); 1818 1817 flag = true; … … 1827 1826 for (int i=NDIM;i--;) { 1828 1827 tmp = Binder->leftatom->at(i) - Binder->rightatom->at(i); 1829 double MinBondDistance; 1830 double MaxBondDistance; 1831 BG->getMinMaxDistance(Binder->leftatom, Binder->rightatom, MinBondDistance, MaxBondDistance, true); 1832 if (fabs(tmp) > MaxBondDistance) // check against Min is not useful for components 1828 range<double> MinMaxDistance(0.,0.); 1829 BG->getMinMaxDistance(Binder->leftatom, Binder->rightatom, MinMaxDistance, true); 1830 if (fabs(tmp) > MinMaxDistance.last) // check against Min is not useful for components 1833 1831 Translationvector[i] = (tmp < 0) ? +1. : -1.; 1834 1832 } -
src/molecule_geometry.cpp
r0cbad2 r300220 393 393 for (int j=0;j<NDIM;j++) { 394 394 tmp = (*iter)->at(j) - (*Runner)->GetOtherAtom(*iter)->at(j); 395 double MinBondDistance; 396 double MaxBondDistance; 397 BG->getMinMaxDistance((*iter), (*Runner)->GetOtherAtom(*iter), MinBondDistance, MaxBondDistance, true); 398 if (fabs(tmp) > MaxBondDistance) { // check against Min is not useful for components 395 range<double> MinMaxBondDistance(0.,0.); 396 BG->getMinMaxDistance((*iter), (*Runner)->GetOtherAtom(*iter), MinMaxBondDistance, true); 397 if (fabs(tmp) > MinMaxBondDistance.last) { // check against Min is not useful for components 399 398 flag = false; 400 399 DoLog(0) && (Log() << Verbose(0) << "Hit: atom " << (*iter)->getName() << " in bond " << *(*Runner) << " has to be shifted due to " << tmp << "." << endl);
Note:
See TracChangeset
for help on using the changeset viewer.