source: src/Atom/atom_bondedparticle.cpp@ 8b58ac

/*
 * Project: MoleCuilder
 * Description: creates and alters molecular systems
 * Copyright (C)  2010-2012 University of Bonn. All rights reserved.
 * 
 *
 *   This file is part of MoleCuilder.
 *
 *    MoleCuilder is free software: you can redistribute it and/or modify
 *    it under the terms of the GNU General Public License as published by
 *    the Free Software Foundation, either version 2 of the License, or
 *    (at your option) any later version.
 *
 *    MoleCuilder is distributed in the hope that it will be useful,
 *    but WITHOUT ANY WARRANTY; without even the implied warranty of
 *    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *    GNU General Public License for more details.
 *
 *    You should have received a copy of the GNU General Public License
 *    along with MoleCuilder.  If not, see <http://www.gnu.org/licenses/>.
 */

/*
 * atom_bondedparticle.cpp
 *
 *  Created on: Oct 19, 2009
 *      Author: heber
 */

// include config.h
#ifdef HAVE_CONFIG_H
#include <config.h>
#endif

#include "CodePatterns/MemDebug.hpp"

#include "atom.hpp"
#include "atom_bondedparticle.hpp"
#include "Bond/bond.hpp"
#include "CodePatterns/Assert.hpp"
#include "CodePatterns/Log.hpp"
#include "CodePatterns/Verbose.hpp"
#include "Element/element.hpp"
#include "WorldTime.hpp"

/** Constructor of class BondedParticle.
 */
BondedParticle::BondedParticle()
{
  ListOfBonds.push_back(BondList());
};

/** Destructor of class BondedParticle.
 */
BondedParticle::~BondedParticle()
{
  removeAllBonds();
};

/** Outputs the current atom::AdaptiveOrder and atom::MaxOrder to \a *file.
 * \param *file output stream
 */
void BondedParticle::OutputOrder(ofstream *file) const
{
  *file << getNr() << "\t" << (int)AdaptiveOrder << "\t" << (int)MaxOrder << endl;
  //LOG(2, "Storing: " << getNr() << "\t" << (int)AdaptiveOrder << "\t" << (int)MaxOrder << ".");
};

/** Prints all bonds of this atom with total degree.
 */
void BondedParticle::OutputBondOfAtom(std::ostream &ost) const
{
  const BondList& ListOfBonds = getListOfBonds();
  ost << "Atom " << getName() << "/" << getNr() << " with " << ListOfBonds.size() << " bonds: ";
  int TotalDegree = 0;
  for (BondList::const_iterator Runner = ListOfBonds.begin(); Runner != ListOfBonds.end(); ++Runner) {
    ost << **Runner << "\t";
    TotalDegree += (*Runner)->BondDegree;
  }
  ost << " -- TotalDegree: " << TotalDegree;
};

/** Output of atom::Nr along with all bond partners.
 * \param *AdjacencyFile output stream
 */
void BondedParticle::OutputAdjacency(ofstream * const AdjacencyFile) const
{
  const BondList& ListOfBonds = getListOfBonds();
  *AdjacencyFile << getNr() << "\t";
  for (BondList::const_iterator Runner = ListOfBonds.begin(); Runner != ListOfBonds.end(); (++Runner))
    *AdjacencyFile << (*Runner)->GetOtherAtom(this)->getNr() << "\t";
  *AdjacencyFile << endl;
};

/** Output of atom::Nr along each bond partner per line.
 * Only bonds are printed where atom::Nr is smaller than the one of the bond partner.
 * \param *AdjacencyFile output stream
 */
void BondedParticle::OutputBonds(ofstream * const BondFile) const
{
  const BondList& ListOfBonds = getListOfBonds();
  for (BondList::const_iterator Runner = ListOfBonds.begin(); Runner != ListOfBonds.end(); (++Runner))
    if (getNr() < (*Runner)->GetOtherAtom(this)->getNr())
      *BondFile << getNr() << "\t" << (*Runner)->GetOtherAtom(this)->getNr() << "\n";
};

/**
 * Adds a bond between this bonded particle and another. Returns present instance if this
 * bond already exists.
 *
 * @param _step time step to access
 * @param bonding partner
 * @return const pointer to created bond or to already present bonds
 */
bond * const BondedParticle::addBond(const unsigned int _step, BondedParticle* Partner)
{
  const BondList &bondlist = getListOfBondsAtStep(_step);
  for (BondList::const_iterator runner = bondlist.begin();
      runner != bondlist.end();
      runner++) {
    if ((*runner)->Contains(Partner))
      return *runner;
  }

  bond* newBond = new bond((atom*) this, (atom*) Partner, 1);
  RegisterBond(_step, newBond);
  Partner->RegisterBond(_step, newBond);

  return newBond;
}

/** Removes a bond for this atom.
 *
 * @param Binder bond to remove
 */
void BondedParticle::removeBond(bond * binder)
{
  UnregisterBond(binder);
}

/** Removes all bonds in all timesteps and their instances, too.
 *
 */
void BondedParticle::removeAllBonds()
{
  for (size_t index = 0; index < ListOfBonds.size(); ++index)
    removeAllBonds(index);
}

/** Removes all bonds for a given \a _step and their instances, too.
 *
 * @param _step time step to access
 */
void BondedParticle::removeAllBonds(const unsigned int _step)
{
  for (BondList::iterator iter = ListOfBonds[_step].begin();
      !ListOfBonds[_step].empty();
      iter = ListOfBonds[_step].begin()) {
    delete (*iter);
    // unregister/NOTIFY is done by bond::~bond()
  }
}

/** Puts a given bond into atom::ListOfBonds.
 * @param _step time step to access
 * \param *Binder bond to insert
 */
bool BondedParticle::RegisterBond(const unsigned int _step, bond * const Binder)
{
  OBSERVE;
  bool status = false;
  if (Binder != NULL) {
    if (Binder->Contains(this)) {
      //LOG(3,"INFO: Registering bond "<< *Binder << " with atom " << *this << " at step " << _step);
      if (ListOfBonds.size() <= _step)
        ListOfBonds.resize(_step+1);
      ListOfBonds[_step].push_back(Binder);
      if (WorldTime::getTime() == _step)
        NOTIFY(AtomObservable::BondsAdded);
      status = true;
    } else {
      ELOG(1, *Binder << " does not contain " << *this << ".");
    }
  } else {
    ELOG(1, "Binder is " << Binder << ".");
  }
  return status;
};

/** Removes a given bond from atom::ListOfBonds.
 * @param _step time step to access
 * \param *Binder bond to remove
 */
bool BondedParticle::UnregisterBond(bond * const Binder)
{
  OBSERVE;
  bool status = false;
  ASSERT(Binder != NULL, "BondedParticle::UnregisterBond() - Binder is NULL.");
  const int step = ContainsBondAtStep(Binder);
  if (step != -1) {
    //LOG(0,"INFO: Unregistering bond "<< *Binder << " from list " << &ListOfBonds << " of atom " << *this << " at step " << step);
    ListOfBonds[step].remove(Binder);
    if (WorldTime::getTime() == step)
      NOTIFY(AtomObservable::BondsRemoved);
    status = true;
  } else {
    ELOG(1, *Binder << " does not contain " << *this << ".");
  }
  return status;
};

/** Removes all bonds from atom::ListOfBonds.
 * \note Does not do any memory de-allocation.
 */
void BondedParticle::UnregisterAllBond(const unsigned int _step)
{
  OBSERVE;
  NOTIFY(AtomObservable::BondsRemoved);
  ListOfBonds[_step].clear();
}

/** Removes all bonds of given \a _step with freeing memory.
 *
 * @param _step time step whose bonds to free
 */
void BondedParticle::ClearBondsAtStep(const unsigned int _step)
{
  OBSERVE;
  NOTIFY(AtomObservable::BondsRemoved);
  //LOG(3,"INFO: Clearing all bonds of " << *this << ": " << ListOfBonds[_step]);
  for (BondList::iterator iter = (ListOfBonds[_step]).begin();
      !(ListOfBonds[_step]).empty();
      iter = (ListOfBonds[_step]).begin()) {
    //LOG(3,"INFO: Clearing bond (" << *iter << ") " << *(*iter) << " of list " << &ListOfBonds);
    delete((*iter)); // will also unregister with us and remove from list
  }
}

/** Searches for the time step where the given bond \a *Binder is a bond of this particle.
 *
 * @param Binder bond to check
 * @return >=0 - first time step where bond appears, -1 - bond not present in lists
 */
int BondedParticle::ContainsBondAtStep(bond *Binder) const
{
  int step = -1;
  int tempstep = 0;
  for(std::vector<BondList>::const_iterator iter = ListOfBonds.begin();
      iter != ListOfBonds.end();
      ++iter,++tempstep) {
    for (BondList::const_iterator bonditer = iter->begin();
        bonditer != iter->end();
        ++bonditer) {
      if ((*bonditer) == Binder) {
        step = tempstep;
        break;
      }
    }
    if (step != -1)
      break;
  }

  return step;
}

/** Corrects the bond degree by one at most if necessary.
 * \return number of corrections done
 */
int BondedParticle::CorrectBondDegree()
{
  OBSERVE;
  NOTIFY(AtomObservable::BondDegreeChanged);
  int NoBonds = 0;
  int OtherNoBonds = 0;
  int FalseBondDegree = 0;
  atom *OtherWalker = NULL;
  bond *CandidateBond = NULL;

  NoBonds = CountBonds();
  //LOG(3, "Walker " << *this << ": " << (int)this->type->NoValenceOrbitals << " > " << NoBonds << "?");
  if ((int)(getType()->getNoValenceOrbitals()) > NoBonds) { // we have a mismatch, check all bonding partners for mismatch
    const BondList& ListOfBonds = getListOfBonds();
    for (BondList::const_iterator Runner = ListOfBonds.begin(); Runner != ListOfBonds.end(); (++Runner)) {
      OtherWalker = (*Runner)->GetOtherAtom(this);
      OtherNoBonds = OtherWalker->CountBonds();
      //LOG(3, "OtherWalker " << *OtherWalker << ": " << (int)OtherWalker->type->NoValenceOrbitals << " > " << OtherNoBonds << "?");
      if ((int)(OtherWalker->getType()->getNoValenceOrbitals()) > OtherNoBonds) { // check if possible candidate
        const BondList& OtherListOfBonds = OtherWalker->getListOfBonds();
        if ((CandidateBond == NULL) || (ListOfBonds.size() > OtherListOfBonds.size())) { // pick the one with fewer number of bonds first
          CandidateBond = (*Runner);
          //LOG(3, "New candidate is " << *CandidateBond << ".");
        }
      }
    }
    if ((CandidateBond != NULL)) {
      CandidateBond->BondDegree++;
      //LOG(2, "Increased bond degree for bond " << *CandidateBond << ".");
    } else {
      ELOG(2, "Could not find correct degree for atom " << *this << ".");
      FalseBondDegree++;
    }
  }
  return FalseBondDegree;
};

/** Sets the weight of all connected bonds to one.
 */
void BondedParticle::resetBondDegree()
{
  OBSERVE;
  NOTIFY(BondedParticle::BondDegreeChanged);
  for (std::vector<BondList>::iterator Runner = ListOfBonds.begin();
      Runner != ListOfBonds.end();
      ++Runner)
    for (BondList::iterator BondRunner = (*Runner).begin();
        BondRunner != (*Runner).end();
        ++BondRunner)
      (*BondRunner)->BondDegree = 1;
};

/** Counts the number of bonds weighted by bond::BondDegree.
   * @param _step time step to access
 * \param bonds times bond::BondDegree
 */
int BondedParticle::CountBonds() const
{
  int NoBonds = 0;
  const BondList& ListOfBonds = getListOfBonds();
  for (BondList::const_iterator Runner = ListOfBonds.begin();
      Runner != ListOfBonds.end();
      (++Runner))
    NoBonds += (*Runner)->BondDegree;
  return NoBonds;
};

/** Checks whether there is a bond between \a this atom and the given \a *BondPartner.
 * @param _step time step to access
 * \param *BondPartner atom to check for
 * \return true - bond exists, false - bond does not exist
 */
bool BondedParticle::IsBondedTo(const unsigned int _step, BondedParticle * const BondPartner) const
{
  bool status = false;

  const BondList& ListOfBonds = getListOfBondsAtStep(_step);
  for (BondList::const_iterator runner = ListOfBonds.begin();
      runner != ListOfBonds.end();
      runner++) {
    status = status || ((*runner)->Contains(BondPartner));
  }
  return status;
};

std::ostream & BondedParticle::operator << (std::ostream &ost) const
{
  ParticleInfo::operator<<(ost);
  ost << "," << getPosition();
  return ost;
}

std::ostream & operator << (std::ostream &ost, const BondedParticle &a)
{
  a.ParticleInfo::operator<<(ost);
  ost << "," << a.getPosition();
  return ost;
}