Changeset 99db9b for src/Analysis/analysis_correlation.cpp
- Timestamp:
- Jul 24, 2015, 4:44:35 PM (10 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, Candidate_v1.7.0, 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:
- fac58f
- Parents:
- a58c16
- git-author:
- Frederik Heber <heber@…> (06/02/15 08:32:10)
- git-committer:
- Frederik Heber <heber@…> (07/24/15 16:44:35)
- File:
-
- 1 edited
-
src/Analysis/analysis_correlation.cpp (modified) (15 diffs)
Legend:
- Unmodified
- Added
- Removed
-
src/Analysis/analysis_correlation.cpp
ra58c16 r99db9b 125 125 * \return range with [min, max] 126 126 */ 127 range<size_t> getMaximumTrajectoryBounds(const std::vector< atom *> &atoms)127 range<size_t> getMaximumTrajectoryBounds(const std::vector<const atom *> &atoms) 128 128 { 129 129 // get highest trajectory size … … 133 133 size_t max_timesteps = std::numeric_limits<size_t>::min(); 134 134 size_t min_timesteps = std::numeric_limits<size_t>::max(); 135 BOOST_FOREACH( atom *_atom, atoms) {135 BOOST_FOREACH(const atom *_atom, atoms) { 136 136 if (_atom->getTrajectorySize() > max_timesteps) 137 137 max_timesteps = _atom->getTrajectorySize(); … … 149 149 * \return map with orientation vector for each atomic id given in \a atoms. 150 150 */ 151 std::map<atomId_t, Vector> CalculateZeroAngularDipole(const std::vector< molecule *> &molecules)151 std::map<atomId_t, Vector> CalculateZeroAngularDipole(const std::vector<const molecule *> &molecules) 152 152 { 153 153 // get zero orientation for each molecule. 154 154 LOG(0,"STATUS: Calculating dipoles for current time step ..."); 155 155 std::map<atomId_t, Vector> ZeroVector; 156 BOOST_FOREACH(molecule *_mol, molecules) { 157 const Vector Dipole = 158 getDipole( 159 const_cast<const molecule *>(_mol)->begin(), 160 const_cast<const molecule *>(_mol)->end()); 161 for(molecule::const_iterator iter = const_cast<const molecule *>(_mol)->begin(); 162 iter != const_cast<const molecule *>(_mol)->end(); 163 ++iter) 156 BOOST_FOREACH(const molecule *_mol, molecules) { 157 const Vector Dipole = getDipole(_mol->begin(),_mol->end()); 158 for(molecule::const_iterator iter = _mol->begin(); iter != _mol->end(); ++iter) 164 159 ZeroVector[(*iter)->getId()] = Dipole; 165 160 LOG(2,"INFO: Zero alignment for molecule " << _mol->getId() << " is " << Dipole); … … 205 200 World::getInstance().clearMoleculeSelection(); 206 201 World::getInstance().selectAllMolecules(MoleculeByFormula(DipoleFormula)); 207 std::vector<molecule *> molecules = World::getInstance().getSelectedMolecules(); 202 std::vector<const molecule *> molecules = 203 const_cast<const World &>(World::getInstance()).getSelectedMolecules(); 208 204 LOG(1,"INFO: There are " << molecules.size() << " molecules for time step " << timestep << "."); 209 205 … … 212 208 size_t i=0; 213 209 size_t Counter_rejections = 0; 214 BOOST_FOREACH(molecule *_mol, molecules) { 215 const Vector Dipole = 216 getDipole( 217 const_cast<const molecule *>(_mol)->begin(), 218 const_cast<const molecule *>(_mol)->end()); 210 BOOST_FOREACH(const molecule *_mol, molecules) { 211 const Vector Dipole = getDipole(_mol->begin(),_mol->end()); 219 212 LOG(3,"INFO: Dipole vector at time step " << timestep << " for for molecule " 220 213 << _mol->getId() << " is " << Dipole); 221 214 // check that all atoms are valid (zeroVector known) 222 molecule::const_iterator iter = const_cast<const molecule *>(_mol)->begin();223 for(; iter != const_cast<const molecule *>(_mol)->end(); ++iter) {215 molecule::const_iterator iter = _mol->begin(); 216 for(; iter != _mol->end(); ++iter) { 224 217 if (!ZeroVector.count((*iter)->getId())) 225 218 break; 226 219 } 227 if (iter != const_cast<const molecule *>(_mol)->end()) {220 if (iter != _mol->end()) { 228 221 ELOG(2, "Skipping molecule " << _mol->getName() << " as not all atoms have a valid zeroVector."); 229 222 ++Counter_rejections; 230 223 continue; 231 224 } else 232 iter = const_cast<const molecule *>(_mol)->begin();225 iter = _mol->begin(); 233 226 std::map<atomId_t, Vector>::const_iterator zeroValue = ZeroVector.find((*iter)->getId()); //due to iter is const 234 227 double angle = 0.; … … 274 267 * \return Map of doubles with values the pair of the two atoms. 275 268 */ 276 DipoleCorrelationMap *DipoleCorrelation(std::vector<molecule *> &molecules) 269 DipoleCorrelationMap *DipoleCorrelation( 270 const std::vector<const molecule *> &molecules) 277 271 { 278 272 Info FunctionInfo(__func__); … … 286 280 } 287 281 288 for (std::vector< molecule *>::const_iterator MolWalker = molecules.begin();282 for (std::vector<const molecule *>::const_iterator MolWalker = molecules.begin(); 289 283 MolWalker != molecules.end(); ++MolWalker) { 290 284 LOG(2, "INFO: Current molecule is " << (*MolWalker)->getId() << "."); 291 const Vector Dipole = 292 getDipole( 293 const_cast<const molecule *>(*MolWalker)->begin(), 294 const_cast<const molecule *>(*MolWalker)->end()); 295 std::vector<molecule *>::const_iterator MolOtherWalker = MolWalker; 296 for (++MolOtherWalker; 297 MolOtherWalker != molecules.end(); 285 const Vector Dipole = getDipole((*MolWalker)->begin(), (*MolWalker)->end()); 286 std::vector<const molecule *>::const_iterator MolOtherWalker = MolWalker; 287 for (++MolOtherWalker; MolOtherWalker != molecules.end(); 298 288 ++MolOtherWalker) { 299 289 LOG(2, "INFO: Current other molecule is " << (*MolOtherWalker)->getId() << "."); 300 const Vector OtherDipole = getDipole( 301 const_cast<const molecule *>(*MolOtherWalker)->begin(), 302 const_cast<const molecule *>(*MolOtherWalker)->end()); 290 const Vector OtherDipole = getDipole((*MolOtherWalker)->begin(), (*MolOtherWalker)->end()); 303 291 const double angle = Dipole.Angle(OtherDipole) * (180./M_PI); 304 292 LOG(1, "Angle is " << angle << "."); … … 409 397 * \return Map of dobules with values as pairs of atom and the vector 410 398 */ 411 CorrelationToPointMap *CorrelationToPoint(std::vector<molecule *> &molecules, const std::vector<const element *> &elements, const Vector *point ) 399 CorrelationToPointMap *CorrelationToPoint( 400 const std::vector<const molecule *> &molecules, 401 const std::vector<const element *> &elements, 402 const Vector *point ) 412 403 { 413 404 Info FunctionInfo(__func__); … … 421 412 } 422 413 423 for (std::vector<molecule *>::const_iterator MolWalker = molecules.begin(); MolWalker != molecules.end(); MolWalker++) { 414 for (std::vector<const molecule *>::const_iterator MolWalker = molecules.begin(); 415 MolWalker != molecules.end(); 416 MolWalker++) { 424 417 LOG(2, "Current molecule is " << *MolWalker << "."); 425 for (molecule::const_iterator iter = const_cast<const molecule *>(*MolWalker)->begin(); 426 iter != const_cast<const molecule *>(*MolWalker)->end(); 427 ++iter) { 418 for (molecule::const_iterator iter = (*MolWalker)->begin(); iter != (*MolWalker)->end(); ++iter) { 428 419 LOG(3, "Current atom is " << **iter << "."); 429 420 for (vector<const element *>::const_iterator type = elements.begin(); type != elements.end(); ++type) … … 453 444 * \return Map of dobules with values as pairs of atom and the vector 454 445 */ 455 CorrelationToPointMap *PeriodicCorrelationToPoint(std::vector<molecule *> &molecules, const std::vector<const element *> &elements, const Vector *point, const int ranges[NDIM] ) 446 CorrelationToPointMap *PeriodicCorrelationToPoint( 447 const std::vector<const molecule *> &molecules, 448 const std::vector<const element *> &elements, 449 const Vector *point, 450 const int ranges[NDIM] ) 456 451 { 457 452 Info FunctionInfo(__func__); … … 467 462 } 468 463 469 for (std::vector<molecule *>::const_iterator MolWalker = molecules.begin(); MolWalker != molecules.end(); MolWalker++) { 464 for (std::vector<const molecule *>::const_iterator MolWalker = molecules.begin(); 465 MolWalker != molecules.end(); 466 MolWalker++) { 470 467 RealSpaceMatrix FullMatrix = World::getInstance().getDomain().getM(); 471 468 RealSpaceMatrix FullInverseMatrix = World::getInstance().getDomain().getMinv(); 472 469 LOG(2, "Current molecule is " << *MolWalker << "."); 473 for (molecule::const_iterator iter = const_cast<const molecule *>(*MolWalker)->begin(); 474 iter != const_cast<const molecule *>(*MolWalker)->end(); 475 ++iter) { 470 for (molecule::const_iterator iter = (*MolWalker)->begin(); iter != (*MolWalker)->end(); ++iter) { 476 471 LOG(3, "Current atom is " << **iter << "."); 477 472 for (vector<const element *>::const_iterator type = elements.begin(); type != elements.end(); ++type) … … 509 504 * \return Map of doubles with values as pairs of atom and the BoundaryTriangleSet that's closest 510 505 */ 511 CorrelationToSurfaceMap *CorrelationToSurface(std::vector<molecule *> &molecules, const std::vector<const element *> &elements, const Tesselation * const Surface, const LinkedCell_deprecated *LC ) 506 CorrelationToSurfaceMap *CorrelationToSurface( 507 const std::vector<const molecule *> &molecules, 508 const std::vector<const element *> &elements, 509 const Tesselation * const Surface, 510 const LinkedCell_deprecated *LC ) 512 511 { 513 512 Info FunctionInfo(__func__); … … 522 521 } 523 522 524 for (std::vector<molecule *>::const_iterator MolWalker = molecules.begin(); MolWalker != molecules.end(); MolWalker++) { 523 for (std::vector<const molecule *>::const_iterator MolWalker = molecules.begin(); 524 MolWalker != molecules.end(); 525 MolWalker++) { 525 526 LOG(2, "Current molecule is " << (*MolWalker)->name << "."); 526 527 if ((*MolWalker)->empty()) 527 528 LOG(2, "\t is empty."); 528 for (molecule::const_iterator iter = const_cast<const molecule *>(*MolWalker)->begin(); 529 iter != const_cast<const molecule *>(*MolWalker)->end(); 530 ++iter) { 529 for (molecule::const_iterator iter = (*MolWalker)->begin(); iter != (*MolWalker)->end(); ++iter) { 531 530 LOG(3, "\tCurrent atom is " << *(*iter) << "."); 532 531 for (vector<const element *>::const_iterator type = elements.begin(); type != elements.end(); ++type) … … 563 562 * \return Map of doubles with values as pairs of atom and the BoundaryTriangleSet that's closest 564 563 */ 565 CorrelationToSurfaceMap *PeriodicCorrelationToSurface(std::vector<molecule *> &molecules, const std::vector<const element *> &elements, const Tesselation * const Surface, const LinkedCell_deprecated *LC, const int ranges[NDIM] ) 564 CorrelationToSurfaceMap *PeriodicCorrelationToSurface( 565 const std::vector<const molecule *> &molecules, 566 const std::vector<const element *> &elements, 567 const Tesselation * const Surface, 568 const LinkedCell_deprecated *LC, 569 const int ranges[NDIM] ) 566 570 { 567 571 Info FunctionInfo(__func__); … … 581 585 double ShortestDistance = 0.; 582 586 BoundaryTriangleSet *ShortestTriangle = NULL; 583 for (std::vector<molecule *>::const_iterator MolWalker = molecules.begin(); MolWalker != molecules.end(); MolWalker++) { 587 for (std::vector<const molecule *>::const_iterator MolWalker = molecules.begin(); 588 MolWalker != molecules.end(); 589 MolWalker++) { 584 590 RealSpaceMatrix FullMatrix = World::getInstance().getDomain().getM(); 585 591 RealSpaceMatrix FullInverseMatrix = World::getInstance().getDomain().getMinv(); 586 592 LOG(2, "Current molecule is " << *MolWalker << "."); 587 for (molecule::const_iterator iter = const_cast<const molecule *>(*MolWalker)->begin(); 588 iter != const_cast<const molecule *>(*MolWalker)->end(); 589 ++iter) { 593 for (molecule::const_iterator iter = (*MolWalker)->begin(); iter != (*MolWalker)->end(); ++iter) { 590 594 LOG(3, "Current atom is " << **iter << "."); 591 595 for (vector<const element *>::const_iterator type = elements.begin(); type != elements.end(); ++type)
Note:
See TracChangeset
for help on using the changeset viewer.
