Changeset 300220


Ignore:
Timestamp:
Mar 1, 2011, 1:17:08 PM (14 years ago)
Author:
Frederik Heber <heber@…>
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)
Message:

BondGraph::..MinMaxDistance() changed to use range<double>.

Files:
7 edited

Legend:

Unmodified
Added
Removed
  • configure.ac

    r0cbad2 r300220  
    115115#AC_MSG_NOTICE(["GLU_CFLAGS: $GLU_CFLAGS, GLU_CXXFLAGS: $GLU_CXXFLAGS, GLU_LDFLAGS: $GLU_LDFLAGS, GLU_LIBS: $GLU_LIBS"])
    116116
    117 # CodePatterns library (needs ObservedIterator<> operator==,!= to be const)
    118 AM_PATH_CODEPATTERNS([1.0.11], $have_debug,,[AC_MSG_ERROR([Missing CodePatterns library, please specify PKG_CONFIG_PATH or use --with-codepatterns...])])
     117# CodePatterns library (needs operator<<(.., range<>) )
     118AM_PATH_CODEPATTERNS([1.0.13], $have_debug,,[AC_MSG_ERROR([Missing CodePatterns library, please specify PKG_CONFIG_PATH or use --with-codepatterns...])])
    119119
    120120# Checks for header files.
  • src/Actions/FragmentationAction/CreateAdjacencyAction.cpp

    r0cbad2 r300220  
    2727#include "bondgraph.hpp"
    2828#include "CodePatterns/Log.hpp"
     29#include "CodePatterns/Range.hpp"
    2930#include "CodePatterns/Verbose.hpp"
    3031#include "config.hpp"
     
    6566  atom *OtherWalker = NULL;
    6667  int n[NDIM];
    67   double MinDistance, MaxDistance;
    6868  LinkedCell *LC = NULL;
    6969  Box &domain = World::getInstance().getDomain();
     
    136136                          ASSERT(OtherWalker,"TesselPoint that was not an atom retrieved from LinkedNode");
    137137                          //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);
    139140                          const double distance = domain.periodicDistanceSquared(OtherWalker->getPosition(),Walker->getPosition());
    140                           const bool status = (distance <= MaxDistance * MaxDistance) && (distance >= MinDistance * MinDistance);
    141 //                          Log() << Verbose(1) << "MinDistance is " << MinDistance << " and MaxDistance is " << MaxDistance << "." << endl;
     141                          const bool status = MinMaxDistanceSquared.isInRange(distance);
     142//                          LOG3, "INFO: MinMaxDistanceSquared interval is " << MinMaxDistanceSquared << ".");
    142143                          if (AtomIds[OtherWalker->father] > AtomIds[Walker->father]) {
    143144                            if (status) { // create bond if distance is smaller
  • src/bondgraph.cpp

    r0cbad2 r300220  
    9393}
    9494
    95 double BondGraph::GetBondLength(int firstZ, int secondZ)
     95double BondGraph::GetBondLength(
     96    int firstZ,
     97    int secondZ) const
    9698{
    9799  std::cout << "Request for length between " << firstZ << " and " << secondZ << ": ";
     
    110112}
    111113
    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)
     114void 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
     126void 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
     141void BondGraph::getMinMaxDistance(
     142    const element * const Walker,
     143    const element * const OtherWalker,
     144    range<double> &MinMaxDistance,
     145    bool IsAngstroem) const
    132146{
    133147  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, MinDistance, MaxDistance, IsAngstroem);
     148    LOG(2, "INFO: Using Covalent radii criterion for [min,max) distances.");
     149    CovalentMinMaxDistance(Walker, OtherWalker, MinMaxDistance, IsAngstroem);
    136150  } 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
     156void 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
     165void 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;
    140176}
    141177
     
    145181  atom *OtherWalker = NULL;
    146182  int n[NDIM];
    147   double MinDistance, MaxDistance;
    148183  Box &domain = World::getInstance().getDomain();
    149184
     
    174209                        ASSERT(OtherWalker != NULL,
    175210                            "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);
    177213                        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);
    181216                        if (OtherWalker->father > Walker->father ) { // just to not add bonds from both sides
    182217                          if (status) { // create bond if distance is smaller
     
    188223                          } else {
    189224                            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 << ".");
    193227                          }
    194228                        } else {
  • src/bondgraph.hpp

    r0cbad2 r300220  
    4444
    4545class 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;
    4650public:
    4751  /** Constructor of class BondGraph.
     
    6266  bool LoadBondLengthTable(std::istream &input);
    6367
    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 
    7168  /** Determines the maximum of all element::CovalentRadius for elements present in \a &Set.
    7269   *
     
    103100
    104101  /** Returns bond criterion for given pair based on a bond length matrix.
    105    * This calls either the covalent or the bond matrix criterion.
     102   * This calls element-version of getMinMaxDistance().
    106103   * \param *Walker first BondedParticle
    107104   * \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
    110106   * \param IsAngstroem whether units are in angstroem or bohr radii
    111107   */
    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;
    113127
    114128  /** Creates the adjacency list for a given \a Range of iterable atoms.
     
    291305  static const double BondThreshold;
    292306
     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
    293316  /** Returns bond criterion for given pair based on a bond length matrix.
    294317   * The matrix should be contained in \a this BondGraph and contain an element-
     
    296319   * \param *Walker first BondedParticle
    297320   * \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
    300322   * \param IsAngstroem whether units are in angstroem or bohr radii
    301323   */
    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
    309336   * \param IsAngstroem whether units are in angstroem or bohr radii
    310337   */
    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;
    312355
    313356
  • src/molecule.cpp

    r0cbad2 r300220  
    344344   // the problem is not the H being out of the box, but InBondvector have the wrong direction
    345345   // 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)) {
    355350//    Log() << Verbose(4) << "InBondvector is: ";
    356351//    InBondvector.Output(out);
     
    359354    for (int i=NDIM;i--;) {
    360355      l = TopReplacement->at(i) - TopOrigin->at(i);
    361       if (fabs(l) > MaxBondDistance) { // 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)
    362357        Orthovector1[i] = (l < 0) ? -1. : +1.;
    363358      } // (signs are correct, was tested!)
  • src/molecule_fragmentation.cpp

    r0cbad2 r300220  
    18111811          tmp = fabs(Binder->leftatom->at(i) - Binder->rightatom->at(i));
    18121812          //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)) {
    18171816            DoLog(2) && (Log() << Verbose(2) << "Correcting at bond " << *Binder << "." << endl);
    18181817            flag = true;
     
    18271826      for (int i=NDIM;i--;) {
    18281827        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
    18331831          Translationvector[i] = (tmp < 0) ? +1. : -1.;
    18341832      }
  • src/molecule_geometry.cpp

    r0cbad2 r300220  
    393393            for (int j=0;j<NDIM;j++) {
    394394              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
    399398                flag = false;
    400399                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.