Ignore:
File:
1 edited

Legend:

Unmodified
Added
Removed
  • src/molecule.cpp

    ra67d19 r68f03d  
    66
    77#include <cstring>
    8 
     8#include <boost/bind.hpp>
     9
     10#include "World.hpp"
    911#include "atom.hpp"
    1012#include "bond.hpp"
     
    2426#include "vector.hpp"
    2527#include "World.hpp"
     28#include "Plane.hpp"
     29#include "Exceptions/LinearDependenceException.hpp"
     30
    2631
    2732/************************************* Functions for class molecule *********************************/
     
    3035 * Initialises molecule list with correctly referenced start and end, and sets molecule::last_atom to zero.
    3136 */
    32 molecule::molecule(const periodentafel * const teil) : elemente(teil), start(new atom), end(new atom),
     37molecule::molecule(const periodentafel * const teil) : elemente(teil), start(World::getInstance().createAtom()), end(World::getInstance().createAtom()),
    3338  first(new bond(start, end, 1, -1)), last(new bond(start, end, 1, -1)), MDSteps(0), AtomCount(0),
    3439  BondCount(0), ElementCount(0), NoNonHydrogen(0), NoNonBonds(0), NoCyclicBonds(0), BondDistance(0.),
    35   ActiveFlag(false), IndexNr(-1), last_atom(0), InternalPointer(start)
     40  ActiveFlag(false), IndexNr(-1),
     41  formula(this,boost::bind(&molecule::calcFormula,this)),
     42  last_atom(0),
     43  InternalPointer(start)
    3644{
    3745  // init atom chain list
     
    4654  for(int i=MAX_ELEMENTS;i--;)
    4755    ElementsInMolecule[i] = 0;
    48   strcpy(name,World::get()->DefaultName);
    49 };
     56  strcpy(name,World::getInstance().getDefaultName());
     57};
     58
     59molecule *NewMolecule(){
     60  return new molecule(World::getInstance().getPeriode());
     61}
    5062
    5163/** Destructor of class molecule.
     
    5769  delete(first);
    5870  delete(last);
    59   delete(end);
    60   delete(start);
    61 };
     71  end->getWorld()->destroyAtom(end);
     72  start->getWorld()->destroyAtom(start);
     73};
     74
     75
     76void DeleteMolecule(molecule *mol){
     77  delete mol;
     78}
     79
     80// getter and setter
     81const std::string molecule::getName(){
     82  return std::string(name);
     83}
     84
     85void molecule::setName(const std::string _name){
     86  OBSERVE;
     87  strncpy(name,_name.c_str(),MAXSTRINGSIZE);
     88}
     89
     90moleculeId_t molecule::getId(){
     91  return id;
     92}
     93
     94void molecule::setId(moleculeId_t _id){
     95  id =_id;
     96}
     97
     98const std::string molecule::getFormula(){
     99  return *formula;
     100}
     101
     102std::string molecule::calcFormula(){
     103  std::map<atomicNumber_t,unsigned int> counts;
     104  stringstream sstr;
     105  periodentafel *periode = World::getInstance().getPeriode();
     106  for(atom *Walker = start; Walker != end; Walker = Walker->next) {
     107    counts[Walker->type->getNumber()]++;
     108  }
     109  std::map<atomicNumber_t,unsigned int>::reverse_iterator iter;
     110  for(iter = counts.rbegin(); iter != counts.rend(); ++iter) {
     111    atomicNumber_t Z = (*iter).first;
     112    sstr << periode->FindElement(Z)->symbol << (*iter).second;
     113  }
     114  return sstr.str();
     115}
    62116
    63117
     
    69123bool molecule::AddAtom(atom *pointer)
    70124{
     125  bool retval = false;
     126  OBSERVE;
    71127  if (pointer != NULL) {
    72128    pointer->sort = &pointer->nr;
     
    79135      if (pointer->type->Z != 1)
    80136        NoNonHydrogen++;
    81       if (pointer->Name == NULL) {
    82         Free(&pointer->Name);
    83         pointer->Name = Malloc<char>(6, "molecule::AddAtom: *pointer->Name");
    84         sprintf(pointer->Name, "%2s%02d", pointer->type->symbol, pointer->nr+1);
    85       }
    86     }
    87     return add(pointer, end);
    88   } else
    89     return false;
     137      if(pointer->getName() == "Unknown"){
     138        stringstream sstr;
     139        sstr << pointer->type->symbol << pointer->nr+1;
     140        pointer->setName(sstr.str());
     141      }
     142    }
     143    retval = add(pointer, end);
     144  }
     145  return retval;
    90146};
    91147
     
    97153atom * molecule::AddCopyAtom(atom *pointer)
    98154{
     155  atom *retval = NULL;
     156  OBSERVE;
    99157  if (pointer != NULL) {
    100     atom *walker = new atom(pointer);
    101     walker->Name = Malloc<char>(strlen(pointer->Name) + 1, "atom::atom: *Name");
    102     strcpy (walker->Name, pointer->Name);
     158    atom *walker = pointer->clone();
     159    stringstream sstr;
     160    sstr << pointer->getName();
     161    walker->setName(sstr.str());
    103162    walker->nr = last_atom++;  // increase number within molecule
    104163    add(walker, end);
     
    106165      NoNonHydrogen++;
    107166    AtomCount++;
    108     return walker;
    109   } else
    110     return NULL;
     167    retval=walker;
     168  }
     169  return retval;
    111170};
    112171
     
    147206bool molecule::AddHydrogenReplacementAtom(bond *TopBond, atom *BottomOrigin, atom *TopOrigin, atom *TopReplacement, bool IsAngstroem)
    148207{
     208  bool AllWentWell = true;    // flag gathering the boolean return value of molecule::AddAtom and other functions, as return value on exit
     209  OBSERVE;
    149210  double bondlength;  // bond length of the bond to be replaced/cut
    150211  double bondangle;  // bond angle of the bond to be replaced/cut
    151212  double BondRescale;   // rescale value for the hydrogen bond length
    152   bool AllWentWell = true;    // flag gathering the boolean return value of molecule::AddAtom and other functions, as return value on exit
    153213  bond *FirstBond = NULL, *SecondBond = NULL; // Other bonds in double bond case to determine "other" plane
    154214  atom *FirstOtherAtom = NULL, *SecondOtherAtom = NULL, *ThirdOtherAtom = NULL; // pointer to hydrogen atoms to be added
     
    158218  double *matrix = NULL;
    159219  bond *Binder = NULL;
    160   double * const cell_size = World::get()->cell_size;
     220  double * const cell_size = World::getInstance().getDomain();
    161221
    162222//  Log() << Verbose(3) << "Begin of AddHydrogenReplacementAtom." << endl;
    163223  // create vector in direction of bond
    164   InBondvector.CopyVector(&TopReplacement->x);
    165   InBondvector.SubtractVector(&TopOrigin->x);
     224  InBondvector = TopReplacement->x - TopOrigin->x;
    166225  bondlength = InBondvector.Norm();
    167226
     
    175234    Orthovector1.Zero();
    176235    for (int i=NDIM;i--;) {
    177       l = TopReplacement->x.x[i] - TopOrigin->x.x[i];
     236      l = TopReplacement->x[i] - TopOrigin->x[i];
    178237      if (fabs(l) > BondDistance) { // is component greater than bond distance
    179         Orthovector1.x[i] = (l < 0) ? -1. : +1.;
     238        Orthovector1[i] = (l < 0) ? -1. : +1.;
    180239      } // (signs are correct, was tested!)
    181240    }
    182241    matrix = ReturnFullMatrixforSymmetric(cell_size);
    183242    Orthovector1.MatrixMultiplication(matrix);
    184     InBondvector.SubtractVector(&Orthovector1); // subtract just the additional translation
     243    InBondvector -= Orthovector1; // subtract just the additional translation
    185244    Free(&matrix);
    186245    bondlength = InBondvector.Norm();
     
    194253  BondRescale = TopOrigin->type->HBondDistance[TopBond->BondDegree-1];
    195254  if (BondRescale == -1) {
    196     DoeLog(1) && (eLog()<< Verbose(1) << "There is no typical hydrogen bond distance in replacing bond (" << TopOrigin->Name << "<->" << TopReplacement->Name << ") of degree " << TopBond->BondDegree << "!" << endl);
     255    DoeLog(1) && (eLog()<< Verbose(1) << "There is no typical hydrogen bond distance in replacing bond (" << TopOrigin->getName() << "<->" << TopReplacement->getName() << ") of degree " << TopBond->BondDegree << "!" << endl);
    197256    return false;
    198257    BondRescale = bondlength;
     
    205264  switch(TopBond->BondDegree) {
    206265    case 1:
    207       FirstOtherAtom = new atom();    // new atom
     266      FirstOtherAtom = World::getInstance().createAtom();    // new atom
    208267      FirstOtherAtom->type = elemente->FindElement(1);  // element is Hydrogen
    209       FirstOtherAtom->v.CopyVector(&TopReplacement->v); // copy velocity
     268      FirstOtherAtom->v = TopReplacement->v; // copy velocity
    210269      FirstOtherAtom->FixedIon = TopReplacement->FixedIon;
    211270      if (TopReplacement->type->Z == 1) { // neither rescale nor replace if it's already hydrogen
     
    215274        FirstOtherAtom->father = NULL;  // if we replace hydrogen, we mark it as our father, otherwise we are just an added hydrogen with no father
    216275      }
    217       InBondvector.Scale(&BondRescale);   // rescale the distance vector to Hydrogen bond length
    218       FirstOtherAtom->x.CopyVector(&TopOrigin->x); // set coordination to origin ...
    219       FirstOtherAtom->x.AddVector(&InBondvector);  // ... and add distance vector to replacement atom
     276      InBondvector *= BondRescale;   // rescale the distance vector to Hydrogen bond length
     277      FirstOtherAtom->x = TopOrigin->x; // set coordination to origin ...
     278      FirstOtherAtom->x = InBondvector;  // ... and add distance vector to replacement atom
    220279      AllWentWell = AllWentWell && AddAtom(FirstOtherAtom);
    221280//      Log() << Verbose(4) << "Added " << *FirstOtherAtom << " at: ";
     
    237296            SecondOtherAtom = (*Runner)->GetOtherAtom(TopOrigin);
    238297          } else {
    239             DoeLog(2) && (eLog()<< Verbose(2) << "Detected more than four bonds for atom " << TopOrigin->Name);
     298            DoeLog(2) && (eLog()<< Verbose(2) << "Detected more than four bonds for atom " << TopOrigin->getName());
    240299          }
    241300        }
     
    249308
    250309        // determine the plane of these two with the *origin
    251         AllWentWell = AllWentWell && Orthovector1.MakeNormalVector(&TopOrigin->x, &FirstOtherAtom->x, &SecondOtherAtom->x);
     310        try {
     311          Orthovector1 =Plane(TopOrigin->x, FirstOtherAtom->x, SecondOtherAtom->x).getNormal();
     312        }
     313        catch(LinearDependenceException &excp){
     314          Log() << Verbose(0) << excp;
     315          // TODO: figure out what to do with the Orthovector in this case
     316          AllWentWell = false;
     317        }
    252318      } else {
    253         Orthovector1.GetOneNormalVector(&InBondvector);
     319        Orthovector1.GetOneNormalVector(InBondvector);
    254320      }
    255321      //Log() << Verbose(3)<< "Orthovector1: ";
     
    257323      //Log() << Verbose(0) << endl;
    258324      // orthogonal vector and bond vector between origin and replacement form the new plane
    259       Orthovector1.MakeNormalVector(&InBondvector);
     325      Orthovector1.MakeNormalTo(InBondvector);
    260326      Orthovector1.Normalize();
    261327      //Log() << Verbose(3) << "ReScaleCheck: " << Orthovector1.Norm() << " and " << InBondvector.Norm() << "." << endl;
    262328
    263329      // create the two Hydrogens ...
    264       FirstOtherAtom = new atom();
    265       SecondOtherAtom = new atom();
     330      FirstOtherAtom = World::getInstance().createAtom();
     331      SecondOtherAtom = World::getInstance().createAtom();
    266332      FirstOtherAtom->type = elemente->FindElement(1);
    267333      SecondOtherAtom->type = elemente->FindElement(1);
    268       FirstOtherAtom->v.CopyVector(&TopReplacement->v); // copy velocity
     334      FirstOtherAtom->v = TopReplacement->v; // copy velocity
    269335      FirstOtherAtom->FixedIon = TopReplacement->FixedIon;
    270       SecondOtherAtom->v.CopyVector(&TopReplacement->v); // copy velocity
     336      SecondOtherAtom->v = TopReplacement->v; // copy velocity
    271337      SecondOtherAtom->FixedIon = TopReplacement->FixedIon;
    272338      FirstOtherAtom->father = NULL;  // we are just an added hydrogen with no father
     
    274340      bondangle = TopOrigin->type->HBondAngle[1];
    275341      if (bondangle == -1) {
    276         DoeLog(1) && (eLog()<< Verbose(1) << "There is no typical hydrogen bond angle in replacing bond (" << TopOrigin->Name << "<->" << TopReplacement->Name << ") of degree " << TopBond->BondDegree << "!" << endl);
     342        DoeLog(1) && (eLog()<< Verbose(1) << "There is no typical hydrogen bond angle in replacing bond (" << TopOrigin->getName() << "<->" << TopReplacement->getName() << ") of degree " << TopBond->BondDegree << "!" << endl);
    277343        return false;
    278344        bondangle = 0;
     
    289355      SecondOtherAtom->x.Zero();
    290356      for(int i=NDIM;i--;) { // rotate by half the bond angle in both directions (InBondvector is bondangle = 0 direction)
    291         FirstOtherAtom->x.x[i] = InBondvector.x[i] * cos(bondangle) + Orthovector1.x[i] * (sin(bondangle));
    292         SecondOtherAtom->x.x[i] = InBondvector.x[i] * cos(bondangle) + Orthovector1.x[i] * (-sin(bondangle));
    293       }
    294       FirstOtherAtom->x.Scale(&BondRescale);  // rescale by correct BondDistance
    295       SecondOtherAtom->x.Scale(&BondRescale);
     357        FirstOtherAtom->x[i] = InBondvector[i] * cos(bondangle) + Orthovector1[i] * (sin(bondangle));
     358        SecondOtherAtom->x[i] = InBondvector[i] * cos(bondangle) + Orthovector1[i] * (-sin(bondangle));
     359      }
     360      FirstOtherAtom->x *= BondRescale;  // rescale by correct BondDistance
     361      SecondOtherAtom->x *= BondRescale;
    296362      //Log() << Verbose(3) << "ReScaleCheck: " << FirstOtherAtom->x.Norm() << " and " << SecondOtherAtom->x.Norm() << "." << endl;
    297363      for(int i=NDIM;i--;) { // and make relative to origin atom
    298         FirstOtherAtom->x.x[i] += TopOrigin->x.x[i];
    299         SecondOtherAtom->x.x[i] += TopOrigin->x.x[i];
     364        FirstOtherAtom->x[i] += TopOrigin->x[i];
     365        SecondOtherAtom->x[i] += TopOrigin->x[i];
    300366      }
    301367      // ... and add to molecule
     
    317383    case 3:
    318384      // take the "usual" tetraoidal angle and add the three Hydrogen in direction of the bond (height of the tetraoid)
    319       FirstOtherAtom = new atom();
    320       SecondOtherAtom = new atom();
    321       ThirdOtherAtom = new atom();
     385      FirstOtherAtom = World::getInstance().createAtom();
     386      SecondOtherAtom = World::getInstance().createAtom();
     387      ThirdOtherAtom = World::getInstance().createAtom();
    322388      FirstOtherAtom->type = elemente->FindElement(1);
    323389      SecondOtherAtom->type = elemente->FindElement(1);
    324390      ThirdOtherAtom->type = elemente->FindElement(1);
    325       FirstOtherAtom->v.CopyVector(&TopReplacement->v); // copy velocity
     391      FirstOtherAtom->v = TopReplacement->v; // copy velocity
    326392      FirstOtherAtom->FixedIon = TopReplacement->FixedIon;
    327       SecondOtherAtom->v.CopyVector(&TopReplacement->v); // copy velocity
     393      SecondOtherAtom->v = TopReplacement->v; // copy velocity
    328394      SecondOtherAtom->FixedIon = TopReplacement->FixedIon;
    329       ThirdOtherAtom->v.CopyVector(&TopReplacement->v); // copy velocity
     395      ThirdOtherAtom->v = TopReplacement->v; // copy velocity
    330396      ThirdOtherAtom->FixedIon = TopReplacement->FixedIon;
    331397      FirstOtherAtom->father = NULL;  //  we are just an added hydrogen with no father
     
    334400
    335401      // we need to vectors orthonormal the InBondvector
    336       AllWentWell = AllWentWell && Orthovector1.GetOneNormalVector(&InBondvector);
     402      AllWentWell = AllWentWell && Orthovector1.GetOneNormalVector(InBondvector);
    337403//      Log() << Verbose(3) << "Orthovector1: ";
    338404//      Orthovector1.Output(out);
    339405//      Log() << Verbose(0) << endl;
    340       AllWentWell = AllWentWell && Orthovector2.MakeNormalVector(&InBondvector, &Orthovector1);
     406      try{
     407        Orthovector2 = Plane(InBondvector, Orthovector1,0).getNormal();
     408      }
     409      catch(LinearDependenceException &excp) {
     410        Log() << Verbose(0) << excp;
     411        AllWentWell = false;
     412      }
    341413//      Log() << Verbose(3) << "Orthovector2: ";
    342414//      Orthovector2.Output(out);
     
    355427      factors[1] = f;
    356428      factors[2] = 0.;
    357       FirstOtherAtom->x.LinearCombinationOfVectors(&InBondvector, &Orthovector1, &Orthovector2, factors);
     429      FirstOtherAtom->x.LinearCombinationOfVectors(InBondvector, Orthovector1, Orthovector2, factors);
    358430      factors[1] = -0.5*f;
    359431      factors[2] = g;
    360       SecondOtherAtom->x.LinearCombinationOfVectors(&InBondvector, &Orthovector1, &Orthovector2, factors);
     432      SecondOtherAtom->x.LinearCombinationOfVectors(InBondvector, Orthovector1, Orthovector2, factors);
    361433      factors[2] = -g;
    362       ThirdOtherAtom->x.LinearCombinationOfVectors(&InBondvector, &Orthovector1, &Orthovector2, factors);
     434      ThirdOtherAtom->x.LinearCombinationOfVectors(InBondvector, Orthovector1, Orthovector2, factors);
    363435
    364436      // rescale each to correct BondDistance
     
    368440
    369441      // and relative to *origin atom
    370       FirstOtherAtom->x.AddVector(&TopOrigin->x);
    371       SecondOtherAtom->x.AddVector(&TopOrigin->x);
    372       ThirdOtherAtom->x.AddVector(&TopOrigin->x);
     442      FirstOtherAtom->x += TopOrigin->x;
     443      SecondOtherAtom->x += TopOrigin->x;
     444      ThirdOtherAtom->x += TopOrigin->x;
    373445
    374446      // ... and add to molecule
     
    413485bool molecule::AddXYZFile(string filename)
    414486{
     487
    415488  istringstream *input = NULL;
    416489  int NumberOfAtoms = 0; // atom number in xyz read
     
    426499    return false;
    427500
     501  OBSERVE;
    428502  getline(xyzfile,line,'\n'); // Read numer of atoms in file
    429503  input = new istringstream(line);
     
    436510    MDSteps++;
    437511  for(i=0;i<NumberOfAtoms;i++){
    438     Walker = new atom;
     512    Walker = World::getInstance().createAtom();
    439513    getline(xyzfile,line,'\n');
    440514    istringstream *item = new istringstream(line);
     
    456530    }
    457531    for(j=NDIM;j--;) {
    458       Walker->x.x[j] = x[j];
    459       Walker->Trajectory.R.at(MDSteps-1).x[j] = x[j];
    460       Walker->Trajectory.U.at(MDSteps-1).x[j] = 0;
    461       Walker->Trajectory.F.at(MDSteps-1).x[j] = 0;
     532      Walker->x[j] = x[j];
     533      Walker->Trajectory.R.at(MDSteps-1)[j] = x[j];
     534      Walker->Trajectory.U.at(MDSteps-1)[j] = 0;
     535      Walker->Trajectory.F.at(MDSteps-1)[j] = 0;
    462536    }
    463537    AddAtom(Walker);  // add to molecule
     
    474548molecule *molecule::CopyMolecule()
    475549{
    476   molecule *copy = new molecule(elemente);
     550  molecule *copy = World::getInstance().createMolecule();
    477551  atom *LeftAtom = NULL, *RightAtom = NULL;
    478552
     
    517591 */
    518592molecule* molecule::CopyMoleculeFromSubRegion(const Vector offset, const double *parallelepiped) const {
    519   molecule *copy = new molecule(elemente);
     593  molecule *copy = World::getInstance().createMolecule();
    520594
    521595  ActOnCopyWithEachAtomIfTrue ( &molecule::AddCopyAtom, copy, &atom::IsInParallelepiped, offset, parallelepiped );
     
    543617    add(Binder, last);
    544618  } else {
    545     DoeLog(1) && (eLog()<< Verbose(1) << "Could not add bond between " << atom1->Name << " and " << atom2->Name << " as one or both are not present in the molecule." << endl);
     619    DoeLog(1) && (eLog()<< Verbose(1) << "Could not add bond between " << atom1->getName() << " and " << atom2->getName() << " as one or both are not present in the molecule." << endl);
    546620  }
    547621  return Binder;
     
    603677void molecule::SetBoxDimension(Vector *dim)
    604678{
    605   double * const cell_size = World::get()->cell_size;
    606   cell_size[0] = dim->x[0];
     679  double * const cell_size = World::getInstance().getDomain();
     680  cell_size[0] = dim->at(0);
    607681  cell_size[1] = 0.;
    608   cell_size[2] = dim->x[1];
     682  cell_size[2] = dim->at(1);
    609683  cell_size[3] = 0.;
    610684  cell_size[4] = 0.;
    611   cell_size[5] = dim->x[2];
     685  cell_size[5] = dim->at(2);
    612686};
    613687
     
    622696    AtomCount--;
    623697  } else
    624     DoeLog(1) && (eLog()<< Verbose(1) << "Atom " << pointer->Name << " is of element " << pointer->type->Z << " but the entry in the table of the molecule is 0!" << endl);
     698    DoeLog(1) && (eLog()<< Verbose(1) << "Atom " << pointer->getName() << " is of element " << pointer->type->Z << " but the entry in the table of the molecule is 0!" << endl);
    625699  if (ElementsInMolecule[pointer->type->Z] == 0)  // was last atom of this element?
    626700    ElementCount--;
     
    640714    ElementsInMolecule[pointer->type->Z]--; // decrease number of atom of this element
    641715  else
    642     DoeLog(1) && (eLog()<< Verbose(1) << "Atom " << pointer->Name << " is of element " << pointer->type->Z << " but the entry in the table of the molecule is 0!" << endl);
     716    DoeLog(1) && (eLog()<< Verbose(1) << "Atom " << pointer->getName() << " is of element " << pointer->type->Z << " but the entry in the table of the molecule is 0!" << endl);
    643717  if (ElementsInMolecule[pointer->type->Z] == 0)  // was last atom of this element?
    644718    ElementCount--;
     
    694768bool molecule::CheckBounds(const Vector *x) const
    695769{
    696   double * const cell_size = World::get()->cell_size;
     770  double * const cell_size = World::getInstance().getDomain();
    697771  bool result = true;
    698772  int j =-1;
    699773  for (int i=0;i<NDIM;i++) {
    700774    j += i+1;
    701     result = result && ((x->x[i] >= 0) && (x->x[i] < cell_size[j]));
     775    result = result && ((x->at(i) >= 0) && (x->at(i) < cell_size[j]));
    702776  }
    703777  //return result;
     
    842916        if (Walker->type->Z != 1) // count non-hydrogen atoms whilst at it
    843917          NoNonHydrogen++;
    844         Free(&Walker->Name);
    845         Walker->Name = Malloc<char>(6, "molecule::CountAtoms: *walker->Name");
    846         sprintf(Walker->Name, "%2s%02d", Walker->type->symbol, Walker->nr+1);
    847         DoLog(3) && (Log() << Verbose(3) << "Naming atom nr. " << Walker->nr << " " << Walker->Name << "." << endl);
     918        stringstream sstr;
     919        sstr << Walker->type->symbol << Walker->nr+1;
     920        Walker->setName(sstr.str());
     921        DoLog(3) && (Log() << Verbose(3) << "Naming atom nr. " << Walker->nr << " " << Walker->getName() << "." << endl);
    848922        i++;
    849923      }
     
    9491023    DeterminePeriodicCenter(CenterOfGravity);
    9501024    OtherMolecule->DeterminePeriodicCenter(OtherCenterOfGravity);
    951     DoLog(5) && (Log() << Verbose(5) << "Center of Gravity: ");
    952     CenterOfGravity.Output();
    953     DoLog(0) && (Log() << Verbose(0) << endl << Verbose(5) << "Other Center of Gravity: ");
    954     OtherCenterOfGravity.Output();
    955     DoLog(0) && (Log() << Verbose(0) << endl);
    956     if (CenterOfGravity.DistanceSquared(&OtherCenterOfGravity) > threshold*threshold) {
     1025    DoLog(5) && (Log() << Verbose(5) << "Center of Gravity: " << CenterOfGravity << endl);
     1026    DoLog(5) && (Log() << Verbose(5) << "Other Center of Gravity: " << OtherCenterOfGravity << endl);
     1027    if (CenterOfGravity.DistanceSquared(OtherCenterOfGravity) > threshold*threshold) {
    9571028      DoLog(4) && (Log() << Verbose(4) << "Centers of gravity don't match." << endl);
    9581029      result = false;
     
    10851156  }
    10861157};
     1158
     1159void molecule::flipActiveFlag(){
     1160  ActiveFlag = !ActiveFlag;
     1161}
Note: See TracChangeset for help on using the changeset viewer.