| 1 | /*
|
|---|
| 2 | * Project: MoleCuilder
|
|---|
| 3 | * Description: creates and alters molecular systems
|
|---|
| 4 | * Copyright (C) 2010 University of Bonn. All rights reserved.
|
|---|
| 5 | * Please see the LICENSE file or "Copyright notice" in builder.cpp for details.
|
|---|
| 6 | */
|
|---|
| 7 |
|
|---|
| 8 | /*
|
|---|
| 9 | * atom_bondedparticle.cpp
|
|---|
| 10 | *
|
|---|
| 11 | * Created on: Oct 19, 2009
|
|---|
| 12 | * Author: heber
|
|---|
| 13 | */
|
|---|
| 14 |
|
|---|
| 15 | // include config.h
|
|---|
| 16 | #ifdef HAVE_CONFIG_H
|
|---|
| 17 | #include <config.h>
|
|---|
| 18 | #endif
|
|---|
| 19 |
|
|---|
| 20 | #include "CodePatterns/MemDebug.hpp"
|
|---|
| 21 |
|
|---|
| 22 | #include "atom.hpp"
|
|---|
| 23 | #include "atom_bondedparticle.hpp"
|
|---|
| 24 | #include "Bond/bond.hpp"
|
|---|
| 25 | #include "CodePatterns/Assert.hpp"
|
|---|
| 26 | #include "CodePatterns/Log.hpp"
|
|---|
| 27 | #include "CodePatterns/Verbose.hpp"
|
|---|
| 28 | #include "element.hpp"
|
|---|
| 29 |
|
|---|
| 30 | /** Constructor of class BondedParticle.
|
|---|
| 31 | */
|
|---|
| 32 | BondedParticle::BondedParticle()
|
|---|
| 33 | {
|
|---|
| 34 | ListOfBonds.push_back(BondList());
|
|---|
| 35 | };
|
|---|
| 36 |
|
|---|
| 37 | /** Destructor of class BondedParticle.
|
|---|
| 38 | */
|
|---|
| 39 | BondedParticle::~BondedParticle()
|
|---|
| 40 | {
|
|---|
| 41 | const size_t max = ListOfBonds.size();
|
|---|
| 42 | for (size_t i = 0; i < max; ++i) {
|
|---|
| 43 | ClearBondsAtStep(i);
|
|---|
| 44 | }
|
|---|
| 45 | };
|
|---|
| 46 |
|
|---|
| 47 | /** Outputs the current atom::AdaptiveOrder and atom::MaxOrder to \a *file.
|
|---|
| 48 | * \param *file output stream
|
|---|
| 49 | */
|
|---|
| 50 | void BondedParticle::OutputOrder(ofstream *file) const
|
|---|
| 51 | {
|
|---|
| 52 | *file << getNr() << "\t" << (int)AdaptiveOrder << "\t" << (int)MaxOrder << endl;
|
|---|
| 53 | //Log() << Verbose(2) << "Storing: " << getNr() << "\t" << (int)AdaptiveOrder << "\t" << (int)MaxOrder << "." << endl;
|
|---|
| 54 | };
|
|---|
| 55 |
|
|---|
| 56 | /** Prints all bonds of this atom with total degree.
|
|---|
| 57 | */
|
|---|
| 58 | void BondedParticle::OutputBondOfAtom(std::ostream &ost) const
|
|---|
| 59 | {
|
|---|
| 60 | const BondList& ListOfBonds = getListOfBonds();
|
|---|
| 61 | ost << "Atom " << getName() << "/" << getNr() << " with " << ListOfBonds.size() << " bonds: ";
|
|---|
| 62 | int TotalDegree = 0;
|
|---|
| 63 | for (BondList::const_iterator Runner = ListOfBonds.begin(); Runner != ListOfBonds.end(); ++Runner) {
|
|---|
| 64 | ost << **Runner << "\t";
|
|---|
| 65 | TotalDegree += (*Runner)->BondDegree;
|
|---|
| 66 | }
|
|---|
| 67 | ost << " -- TotalDegree: " << TotalDegree;
|
|---|
| 68 | };
|
|---|
| 69 |
|
|---|
| 70 | /** Output of atom::Nr along with all bond partners.
|
|---|
| 71 | * \param *AdjacencyFile output stream
|
|---|
| 72 | */
|
|---|
| 73 | void BondedParticle::OutputAdjacency(ofstream * const AdjacencyFile) const
|
|---|
| 74 | {
|
|---|
| 75 | const BondList& ListOfBonds = getListOfBonds();
|
|---|
| 76 | *AdjacencyFile << getNr() << "\t";
|
|---|
| 77 | for (BondList::const_iterator Runner = ListOfBonds.begin(); Runner != ListOfBonds.end(); (++Runner))
|
|---|
| 78 | *AdjacencyFile << (*Runner)->GetOtherAtom(this)->getNr() << "\t";
|
|---|
| 79 | *AdjacencyFile << endl;
|
|---|
| 80 | };
|
|---|
| 81 |
|
|---|
| 82 | /** Output of atom::Nr along each bond partner per line.
|
|---|
| 83 | * Only bonds are printed where atom::Nr is smaller than the one of the bond partner.
|
|---|
| 84 | * \param *AdjacencyFile output stream
|
|---|
| 85 | */
|
|---|
| 86 | void BondedParticle::OutputBonds(ofstream * const BondFile) const
|
|---|
| 87 | {
|
|---|
| 88 | const BondList& ListOfBonds = getListOfBonds();
|
|---|
| 89 | for (BondList::const_iterator Runner = ListOfBonds.begin(); Runner != ListOfBonds.end(); (++Runner))
|
|---|
| 90 | if (getNr() < (*Runner)->GetOtherAtom(this)->getNr())
|
|---|
| 91 | *BondFile << getNr() << "\t" << (*Runner)->GetOtherAtom(this)->getNr() << "\n";
|
|---|
| 92 | };
|
|---|
| 93 |
|
|---|
| 94 | /**
|
|---|
| 95 | * Adds a bond between this bonded particle and another. Does nothing if this
|
|---|
| 96 | * bond already exists.
|
|---|
| 97 | *
|
|---|
| 98 | * @param _step time step to access
|
|---|
| 99 | * \param bonding partner
|
|---|
| 100 | */
|
|---|
| 101 | void BondedParticle::addBond(const unsigned int _step, BondedParticle* Partner) {
|
|---|
| 102 | if (IsBondedTo(_step, Partner)) {
|
|---|
| 103 | return;
|
|---|
| 104 | }
|
|---|
| 105 |
|
|---|
| 106 | bond* newBond = new bond((atom*) this, (atom*) Partner, 1, 0);
|
|---|
| 107 | RegisterBond(_step, newBond);
|
|---|
| 108 | Partner->RegisterBond(_step, newBond);
|
|---|
| 109 | }
|
|---|
| 110 |
|
|---|
| 111 | /** Puts a given bond into atom::ListOfBonds.
|
|---|
| 112 | * @param _step time step to access
|
|---|
| 113 | * \param *Binder bond to insert
|
|---|
| 114 | */
|
|---|
| 115 | bool BondedParticle::RegisterBond(const unsigned int _step, bond *Binder)
|
|---|
| 116 | {
|
|---|
| 117 | bool status = false;
|
|---|
| 118 | if (Binder != NULL) {
|
|---|
| 119 | if (Binder->Contains(this)) {
|
|---|
| 120 | //LOG(3,"INFO: Registering bond "<< *Binder << " with atom " << *this << " at step " << _step);
|
|---|
| 121 | BondList& ListOfBonds = getListOfBondsAtStep(_step);
|
|---|
| 122 | ListOfBonds.push_back(Binder);
|
|---|
| 123 | status = true;
|
|---|
| 124 | } else {
|
|---|
| 125 | DoeLog(1) && (eLog()<< Verbose(1) << *Binder << " does not contain " << *this << "." << endl);
|
|---|
| 126 | }
|
|---|
| 127 | } else {
|
|---|
| 128 | DoeLog(1) && (eLog()<< Verbose(1) << "Binder is " << Binder << "." << endl);
|
|---|
| 129 | }
|
|---|
| 130 | return status;
|
|---|
| 131 | };
|
|---|
| 132 |
|
|---|
| 133 | /** Removes a given bond from atom::ListOfBonds.
|
|---|
| 134 | * @param _step time step to access
|
|---|
| 135 | * \param *Binder bond to remove
|
|---|
| 136 | */
|
|---|
| 137 | bool BondedParticle::UnregisterBond(bond *Binder)
|
|---|
| 138 | {
|
|---|
| 139 | bool status = false;
|
|---|
| 140 | ASSERT(Binder != NULL, "BondedParticle::UnregisterBond() - Binder is NULL.");
|
|---|
| 141 | const int step = ContainsBondAtStep(Binder);
|
|---|
| 142 | if (step != -1) {
|
|---|
| 143 | //LOG(3,"INFO: Unregistering bond "<< *Binder << " from list " << &ListOfBonds << " of atom " << *this << " at step " << step);
|
|---|
| 144 | ListOfBonds[step].remove(Binder);
|
|---|
| 145 | status = true;
|
|---|
| 146 | } else {
|
|---|
| 147 | DoeLog(1) && (eLog()<< Verbose(1) << *Binder << " does not contain " << *this << "." << endl);
|
|---|
| 148 | }
|
|---|
| 149 | return status;
|
|---|
| 150 | };
|
|---|
| 151 |
|
|---|
| 152 | /** Removes all bonds from atom::ListOfBonds.
|
|---|
| 153 | * \note Does not do any memory de-allocation.
|
|---|
| 154 | */
|
|---|
| 155 | void BondedParticle::UnregisterAllBond(const unsigned int _step)
|
|---|
| 156 | {
|
|---|
| 157 | ListOfBonds[_step].clear();
|
|---|
| 158 | }
|
|---|
| 159 |
|
|---|
| 160 | /** Removes all bonds of given \a _step with freeing memory.
|
|---|
| 161 | *
|
|---|
| 162 | * @param _step time step whose bonds to free
|
|---|
| 163 | */
|
|---|
| 164 | void BondedParticle::ClearBondsAtStep(const unsigned int _step)
|
|---|
| 165 | {
|
|---|
| 166 | //LOG(3,"INFO: Clearing all bonds of " << *this << ": " << ListOfBonds[_step]);
|
|---|
| 167 | for (BondList::iterator iter = (ListOfBonds[_step]).begin();
|
|---|
| 168 | !(ListOfBonds[_step]).empty();
|
|---|
| 169 | iter = (ListOfBonds[_step]).begin()) {
|
|---|
| 170 | //LOG(3,"INFO: Clearing bond (" << *iter << ") " << *(*iter) << " of list " << &ListOfBonds);
|
|---|
| 171 | delete((*iter)); // will also unregister with us and remove from list
|
|---|
| 172 | }
|
|---|
| 173 | }
|
|---|
| 174 |
|
|---|
| 175 | /** Searches for the time step where the given bond \a *Binder is a bond of this particle.
|
|---|
| 176 | *
|
|---|
| 177 | * @param Binder bond to check
|
|---|
| 178 | * @return >=0 - first time step where bond appears, -1 - bond not present in lists
|
|---|
| 179 | */
|
|---|
| 180 | int BondedParticle::ContainsBondAtStep(bond *Binder)
|
|---|
| 181 | {
|
|---|
| 182 | int step = -1;
|
|---|
| 183 | int tempstep = 0;
|
|---|
| 184 | for(std::vector<BondList>::const_iterator iter = ListOfBonds.begin();
|
|---|
| 185 | iter != ListOfBonds.end();
|
|---|
| 186 | ++iter,++tempstep) {
|
|---|
| 187 | for (BondList::const_iterator bonditer = iter->begin();
|
|---|
| 188 | bonditer != iter->end();
|
|---|
| 189 | ++bonditer) {
|
|---|
| 190 | if ((*bonditer) == Binder) {
|
|---|
| 191 | step = tempstep;
|
|---|
| 192 | break;
|
|---|
| 193 | }
|
|---|
| 194 | }
|
|---|
| 195 | if (step != -1)
|
|---|
| 196 | break;
|
|---|
| 197 | }
|
|---|
| 198 |
|
|---|
| 199 | return step;
|
|---|
| 200 | }
|
|---|
| 201 |
|
|---|
| 202 | /** Corrects the bond degree by one at most if necessary.
|
|---|
| 203 | * \return number of corrections done
|
|---|
| 204 | */
|
|---|
| 205 | int BondedParticle::CorrectBondDegree()
|
|---|
| 206 | {
|
|---|
| 207 | int NoBonds = 0;
|
|---|
| 208 | int OtherNoBonds = 0;
|
|---|
| 209 | int FalseBondDegree = 0;
|
|---|
| 210 | atom *OtherWalker = NULL;
|
|---|
| 211 | bond *CandidateBond = NULL;
|
|---|
| 212 |
|
|---|
| 213 | NoBonds = CountBonds();
|
|---|
| 214 | //Log() << Verbose(3) << "Walker " << *this << ": " << (int)this->type->NoValenceOrbitals << " > " << NoBonds << "?" << endl;
|
|---|
| 215 | if ((int)(getType()->getNoValenceOrbitals()) > NoBonds) { // we have a mismatch, check all bonding partners for mismatch
|
|---|
| 216 | const BondList& ListOfBonds = getListOfBonds();
|
|---|
| 217 | for (BondList::const_iterator Runner = ListOfBonds.begin(); Runner != ListOfBonds.end(); (++Runner)) {
|
|---|
| 218 | OtherWalker = (*Runner)->GetOtherAtom(this);
|
|---|
| 219 | OtherNoBonds = OtherWalker->CountBonds();
|
|---|
| 220 | //Log() << Verbose(3) << "OtherWalker " << *OtherWalker << ": " << (int)OtherWalker->type->NoValenceOrbitals << " > " << OtherNoBonds << "?" << endl;
|
|---|
| 221 | if ((int)(OtherWalker->getType()->getNoValenceOrbitals()) > OtherNoBonds) { // check if possible candidate
|
|---|
| 222 | const BondList& OtherListOfBonds = OtherWalker->getListOfBonds();
|
|---|
| 223 | if ((CandidateBond == NULL) || (ListOfBonds.size() > OtherListOfBonds.size())) { // pick the one with fewer number of bonds first
|
|---|
| 224 | CandidateBond = (*Runner);
|
|---|
| 225 | //Log() << Verbose(3) << "New candidate is " << *CandidateBond << "." << endl;
|
|---|
| 226 | }
|
|---|
| 227 | }
|
|---|
| 228 | }
|
|---|
| 229 | if ((CandidateBond != NULL)) {
|
|---|
| 230 | CandidateBond->BondDegree++;
|
|---|
| 231 | //Log() << Verbose(2) << "Increased bond degree for bond " << *CandidateBond << "." << endl;
|
|---|
| 232 | } else {
|
|---|
| 233 | DoeLog(2) && (eLog()<< Verbose(2) << "Could not find correct degree for atom " << *this << "." << endl);
|
|---|
| 234 | FalseBondDegree++;
|
|---|
| 235 | }
|
|---|
| 236 | }
|
|---|
| 237 | return FalseBondDegree;
|
|---|
| 238 | };
|
|---|
| 239 |
|
|---|
| 240 | /** Counts the number of bonds weighted by bond::BondDegree.
|
|---|
| 241 | * @param _step time step to access
|
|---|
| 242 | * \param bonds times bond::BondDegree
|
|---|
| 243 | */
|
|---|
| 244 | int BondedParticle::CountBonds() const
|
|---|
| 245 | {
|
|---|
| 246 | int NoBonds = 0;
|
|---|
| 247 | const BondList& ListOfBonds = getListOfBonds();
|
|---|
| 248 | for (BondList::const_iterator Runner = ListOfBonds.begin();
|
|---|
| 249 | Runner != ListOfBonds.end();
|
|---|
| 250 | (++Runner))
|
|---|
| 251 | NoBonds += (*Runner)->BondDegree;
|
|---|
| 252 | return NoBonds;
|
|---|
| 253 | };
|
|---|
| 254 |
|
|---|
| 255 | /** Checks whether there is a bond between \a this atom and the given \a *BondPartner.
|
|---|
| 256 | * @param _step time step to access
|
|---|
| 257 | * \param *BondPartner atom to check for
|
|---|
| 258 | * \return true - bond exists, false - bond does not exist
|
|---|
| 259 | */
|
|---|
| 260 | bool BondedParticle::IsBondedTo(const unsigned int _step, BondedParticle * const BondPartner) const
|
|---|
| 261 | {
|
|---|
| 262 | bool status = false;
|
|---|
| 263 |
|
|---|
| 264 | const BondList& ListOfBonds = getListOfBondsAtStep(_step);
|
|---|
| 265 | for (BondList::const_iterator runner = ListOfBonds.begin();
|
|---|
| 266 | runner != ListOfBonds.end();
|
|---|
| 267 | runner++) {
|
|---|
| 268 | status = status || ((*runner)->Contains(BondPartner));
|
|---|
| 269 | }
|
|---|
| 270 | return status;
|
|---|
| 271 | };
|
|---|
| 272 |
|
|---|
| 273 | std::ostream & BondedParticle::operator << (std::ostream &ost) const
|
|---|
| 274 | {
|
|---|
| 275 | ParticleInfo::operator<<(ost);
|
|---|
| 276 | ost << "," << getPosition();
|
|---|
| 277 | return ost;
|
|---|
| 278 | }
|
|---|
| 279 |
|
|---|
| 280 | std::ostream & operator << (std::ostream &ost, const BondedParticle &a)
|
|---|
| 281 | {
|
|---|
| 282 | a.ParticleInfo::operator<<(ost);
|
|---|
| 283 | ost << "," << a.getPosition();
|
|---|
| 284 | return ost;
|
|---|
| 285 | }
|
|---|
| 286 |
|
|---|