| [bcf653] | 1 | /* | 
|---|
|  | 2 | * Project: MoleCuilder | 
|---|
|  | 3 | * Description: creates and alters molecular systems | 
|---|
| [0aa122] | 4 | * Copyright (C)  2010-2012 University of Bonn. All rights reserved. | 
|---|
| [5aaa43] | 5 | * Copyright (C)  2013 Frederik Heber. All rights reserved. | 
|---|
| [94d5ac6] | 6 | * | 
|---|
|  | 7 | * | 
|---|
|  | 8 | *   This file is part of MoleCuilder. | 
|---|
|  | 9 | * | 
|---|
|  | 10 | *    MoleCuilder is free software: you can redistribute it and/or modify | 
|---|
|  | 11 | *    it under the terms of the GNU General Public License as published by | 
|---|
|  | 12 | *    the Free Software Foundation, either version 2 of the License, or | 
|---|
|  | 13 | *    (at your option) any later version. | 
|---|
|  | 14 | * | 
|---|
|  | 15 | *    MoleCuilder is distributed in the hope that it will be useful, | 
|---|
|  | 16 | *    but WITHOUT ANY WARRANTY; without even the implied warranty of | 
|---|
|  | 17 | *    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the | 
|---|
|  | 18 | *    GNU General Public License for more details. | 
|---|
|  | 19 | * | 
|---|
|  | 20 | *    You should have received a copy of the GNU General Public License | 
|---|
|  | 21 | *    along with MoleCuilder.  If not, see <http://www.gnu.org/licenses/>. | 
|---|
| [bcf653] | 22 | */ | 
|---|
|  | 23 |  | 
|---|
| [14de469] | 24 | /** \file atom.cpp | 
|---|
| [1907a7] | 25 | * | 
|---|
| [14de469] | 26 | * Function implementations for the class atom. | 
|---|
| [1907a7] | 27 | * | 
|---|
| [14de469] | 28 | */ | 
|---|
|  | 29 |  | 
|---|
| [bf3817] | 30 | // include config.h | 
|---|
|  | 31 | #ifdef HAVE_CONFIG_H | 
|---|
|  | 32 | #include <config.h> | 
|---|
|  | 33 | #endif | 
|---|
|  | 34 |  | 
|---|
| [9eb71b3] | 35 | //#include "CodePatterns/MemDebug.hpp" | 
|---|
| [112b09] | 36 |  | 
|---|
| [357fba] | 37 | #include "atom.hpp" | 
|---|
| [708277] | 38 | #include "AtomObserver.hpp" | 
|---|
| [129204] | 39 | #include "Bond/bond.hpp" | 
|---|
| [e2373df] | 40 | #include "CodePatterns/Log.hpp" | 
|---|
| [4a7776a] | 41 | #include "config.hpp" | 
|---|
| [3bdb6d] | 42 | #include "Element/element.hpp" | 
|---|
| [57f243] | 43 | #include "LinearAlgebra/Vector.hpp" | 
|---|
| [d346b6] | 44 | #include "World.hpp" | 
|---|
| [11f0fa] | 45 | #include "WorldTime.hpp" | 
|---|
| [6cfa36] | 46 | #include "molecule.hpp" | 
|---|
| [c550dd] | 47 | #include "Shapes/Shape.hpp" | 
|---|
| [a0064e] | 48 |  | 
|---|
| [36166d] | 49 | #include <iomanip> | 
|---|
| [0ba410] | 50 | #include <iostream> | 
|---|
| [36166d] | 51 |  | 
|---|
| [14de469] | 52 | /************************************* Functions for class atom *************************************/ | 
|---|
|  | 53 |  | 
|---|
| [70ff32] | 54 |  | 
|---|
| [d05088] | 55 | atom::atom() : | 
|---|
| [97b825] | 56 | father(this), | 
|---|
| [5309ba] | 57 | sort(&Nr), | 
|---|
| [b71881] | 58 | mol(0), | 
|---|
|  | 59 | selected(false) | 
|---|
| [708277] | 60 | { | 
|---|
| [d05088] | 61 | // note AtomObserver about inserted atom | 
|---|
| [d4ba3f] | 62 | AtomObserver::getInstance().Inserted(this); | 
|---|
| [708277] | 63 | } | 
|---|
| [14de469] | 64 |  | 
|---|
| [ea0c8b] | 65 | atom::atom(atom *pointer) : | 
|---|
|  | 66 | ParticleInfo(*pointer), | 
|---|
| [c742fe] | 67 | AtomInfo(*pointer), | 
|---|
| [97b825] | 68 | father(pointer), | 
|---|
| [5309ba] | 69 | sort(&Nr), | 
|---|
| [b71881] | 70 | mol(0), | 
|---|
|  | 71 | selected(false) | 
|---|
| [2319ed] | 72 | { | 
|---|
| [d05088] | 73 | // sign on to father atom to be notified when it is removed | 
|---|
|  | 74 | father->signOn(this); | 
|---|
|  | 75 |  | 
|---|
|  | 76 | // note AtomObserver about inserted atom | 
|---|
| [d4ba3f] | 77 | AtomObserver::getInstance().Inserted(this); | 
|---|
| [b453f9] | 78 | }; | 
|---|
| [2319ed] | 79 |  | 
|---|
| [46d958] | 80 | atom *atom::clone(){ | 
|---|
| [68f03d] | 81 | atom *res = new atom(this); | 
|---|
| [23b547] | 82 | World::getInstance().registerAtom(res); | 
|---|
| [46d958] | 83 | return res; | 
|---|
|  | 84 | } | 
|---|
|  | 85 |  | 
|---|
| [2319ed] | 86 |  | 
|---|
| [14de469] | 87 | /** Destructor of class atom. | 
|---|
|  | 88 | */ | 
|---|
| [1907a7] | 89 | atom::~atom() | 
|---|
| [14de469] | 90 | { | 
|---|
| [d05088] | 91 | // sign off from possible father | 
|---|
|  | 92 | if ((father != this) && (father != NULL)) | 
|---|
|  | 93 | father->signOff(this); | 
|---|
|  | 94 |  | 
|---|
| [6cfa36] | 95 | removeFromMolecule(); | 
|---|
| [d05088] | 96 | // note AtomObserver about removed atom | 
|---|
| [d4ba3f] | 97 | AtomObserver::getInstance().Removed(this); | 
|---|
| [708277] | 98 | } | 
|---|
| [e2373df] | 99 |  | 
|---|
|  | 100 |  | 
|---|
| [8cc22f] | 101 | void atom::UpdateStep(const unsigned int _step) | 
|---|
| [e2373df] | 102 | { | 
|---|
| [8cc22f] | 103 | LOG(4,"atom::UpdateStep() called."); | 
|---|
| [e2373df] | 104 | // append to position, velocity and force vector | 
|---|
| [05a2c5] | 105 | AtomInfo::AppendTrajectoryStep(WorldTime::getTime()); | 
|---|
| [1e6249] | 106 | // append to ListOfBonds vector | 
|---|
| [05a2c5] | 107 | BondedParticleInfo::AppendTrajectoryStep(WorldTime::getTime()); | 
|---|
| [e2373df] | 108 | } | 
|---|
|  | 109 |  | 
|---|
| [8cc22f] | 110 | void atom::removeStep(const unsigned int _step) | 
|---|
| [7e51e1] | 111 | { | 
|---|
| [8cc22f] | 112 | LOG(4,"atom::removeStep() called."); | 
|---|
| [7e51e1] | 113 | // append to position, velocity and force vector | 
|---|
| [8cc22f] | 114 | AtomInfo::removeTrajectoryStep(_step); | 
|---|
| [7e51e1] | 115 | // append to ListOfBonds vector | 
|---|
| [8cc22f] | 116 | BondedParticleInfo::removeTrajectoryStep(_step); | 
|---|
| [7e51e1] | 117 | } | 
|---|
|  | 118 |  | 
|---|
| [59fff1] | 119 | atom *atom::GetTrueFather() | 
|---|
|  | 120 | { | 
|---|
|  | 121 | const atom *father = const_cast<const atom *>(this)->GetTrueFather(); | 
|---|
|  | 122 | return const_cast<atom *>(father); | 
|---|
|  | 123 | } | 
|---|
|  | 124 |  | 
|---|
|  | 125 | const atom *atom::GetTrueFather() const | 
|---|
| [14de469] | 126 | { | 
|---|
| [215df0] | 127 | if(father == this){ // top most father is the one that points on itself | 
|---|
|  | 128 | return this; | 
|---|
|  | 129 | } | 
|---|
|  | 130 | else if(!father) { | 
|---|
|  | 131 | return 0; | 
|---|
|  | 132 | } | 
|---|
|  | 133 | else { | 
|---|
|  | 134 | return father->GetTrueFather(); | 
|---|
|  | 135 | } | 
|---|
| [910a5d] | 136 | } | 
|---|
|  | 137 |  | 
|---|
|  | 138 | void atom::setFather(atom * const _father) | 
|---|
|  | 139 | { | 
|---|
| [d05088] | 140 | // sign off from old father | 
|---|
|  | 141 | if ((father != this) && (father != NULL)) | 
|---|
|  | 142 | father->signOff(this); | 
|---|
|  | 143 |  | 
|---|
| [910a5d] | 144 | father = _father; | 
|---|
| [d05088] | 145 | father->signOn(this); | 
|---|
| [910a5d] | 146 | } | 
|---|
| [14de469] | 147 |  | 
|---|
| [e65246] | 148 | /** Sets father to itself or its father in case of copying a molecule. | 
|---|
|  | 149 | */ | 
|---|
|  | 150 | void atom::CorrectFather() | 
|---|
|  | 151 | { | 
|---|
| [2e352f] | 152 | if (father->father != father)   // same atom in copy's father points to itself | 
|---|
|  | 153 | //    father = this;  // set father to itself (copy of a whole molecule) | 
|---|
|  | 154 | //  else | 
|---|
| [e65246] | 155 | father = father->father;  // set father to original's father | 
|---|
|  | 156 |  | 
|---|
|  | 157 | }; | 
|---|
|  | 158 |  | 
|---|
| [b453f9] | 159 | void atom::EqualsFather ( const atom *ptr, const atom **res ) const | 
|---|
| [e65246] | 160 | { | 
|---|
|  | 161 | if ( ptr == father ) | 
|---|
|  | 162 | *res = this; | 
|---|
|  | 163 | }; | 
|---|
|  | 164 |  | 
|---|
| [00abfc] | 165 | bool atom::isFather(const atom *ptr){ | 
|---|
|  | 166 | return ptr==father; | 
|---|
|  | 167 | } | 
|---|
|  | 168 |  | 
|---|
| [e138de] | 169 | bool atom::OutputIndexed(ofstream * const out, const int ElementNo, const int AtomNo, const char *comment) const | 
|---|
| [14de469] | 170 | { | 
|---|
|  | 171 | if (out != NULL) { | 
|---|
|  | 172 | *out << "Ion_Type" << ElementNo << "_" << AtomNo << "\t"  << fixed << setprecision(9) << showpoint; | 
|---|
| [d74077] | 173 | *out << at(0) << "\t" << at(1) << "\t" << at(2); | 
|---|
| [6625c3] | 174 | *out << "\t" << (int)(getFixedIon()); | 
|---|
| [bce72c] | 175 | if (getAtomicVelocity().Norm() > MYEPSILON) | 
|---|
|  | 176 | *out << "\t" << scientific << setprecision(6) << getAtomicVelocity()[0] << "\t" << getAtomicVelocity()[1] << "\t" << getAtomicVelocity()[2] << "\t"; | 
|---|
| [437922] | 177 | if (comment != NULL) | 
|---|
|  | 178 | *out << " # " << comment << endl; | 
|---|
| [e9f8f9] | 179 | else | 
|---|
| [735b1c] | 180 | *out << " # molecule nr " << getNr() << endl; | 
|---|
| [e9f8f9] | 181 | return true; | 
|---|
|  | 182 | } else | 
|---|
|  | 183 | return false; | 
|---|
|  | 184 | }; | 
|---|
| [b453f9] | 185 |  | 
|---|
| [0ba410] | 186 | bool atom::OutputArrayIndexed(ostream * const out,const enumeration<const element*> &elementLookup, int *AtomNo, const char *comment) const | 
|---|
| [e9f8f9] | 187 | { | 
|---|
| [83f176] | 188 | AtomNo[getType()->getAtomicNumber()]++;  // increment number | 
|---|
| [e9f8f9] | 189 | if (out != NULL) { | 
|---|
| [8f4df1] | 190 | const element *elemental = getType(); | 
|---|
|  | 191 | ASSERT(elementLookup.there.find(elemental)!=elementLookup.there.end(),"Type of this atom was not in the formula upon enumeration"); | 
|---|
| [83f176] | 192 | *out << "Ion_Type" << elementLookup.there.find(elemental)->second << "_" << AtomNo[elemental->getAtomicNumber()] << "\t"  << fixed << setprecision(9) << showpoint; | 
|---|
| [d74077] | 193 | *out << at(0) << "\t" << at(1) << "\t" << at(2); | 
|---|
| [6625c3] | 194 | *out << "\t" << getFixedIon(); | 
|---|
| [bce72c] | 195 | if (getAtomicVelocity().Norm() > MYEPSILON) | 
|---|
|  | 196 | *out << "\t" << scientific << setprecision(6) << getAtomicVelocity()[0] << "\t" << getAtomicVelocity()[1] << "\t" << getAtomicVelocity()[2] << "\t"; | 
|---|
| [e9f8f9] | 197 | if (comment != NULL) | 
|---|
|  | 198 | *out << " # " << comment << endl; | 
|---|
| [437922] | 199 | else | 
|---|
| [735b1c] | 200 | *out << " # molecule nr " << getNr() << endl; | 
|---|
| [14de469] | 201 | return true; | 
|---|
|  | 202 | } else | 
|---|
|  | 203 | return false; | 
|---|
|  | 204 | }; | 
|---|
|  | 205 |  | 
|---|
| [b453f9] | 206 | bool atom::Compare(const atom &ptr) const | 
|---|
| [4455f4] | 207 | { | 
|---|
| [735b1c] | 208 | if (getNr() < ptr.getNr()) | 
|---|
| [4455f4] | 209 | return true; | 
|---|
|  | 210 | else | 
|---|
|  | 211 | return false; | 
|---|
|  | 212 | }; | 
|---|
|  | 213 |  | 
|---|
| [b453f9] | 214 | double atom::DistanceSquaredToVector(const Vector &origin) const | 
|---|
| [4455f4] | 215 | { | 
|---|
| [d74077] | 216 | return DistanceSquared(origin); | 
|---|
| [4455f4] | 217 | }; | 
|---|
|  | 218 |  | 
|---|
| [b453f9] | 219 | double atom::DistanceToVector(const Vector &origin) const | 
|---|
| [4455f4] | 220 | { | 
|---|
| [d74077] | 221 | return distance(origin); | 
|---|
| [4455f4] | 222 | }; | 
|---|
|  | 223 |  | 
|---|
|  | 224 | void atom::InitComponentNr() | 
|---|
|  | 225 | { | 
|---|
|  | 226 | if (ComponentNr != NULL) | 
|---|
| [920c70] | 227 | delete[](ComponentNr); | 
|---|
| [9d83b6] | 228 | const BondList& ListOfBonds = getListOfBonds(); | 
|---|
| [920c70] | 229 | ComponentNr = new int[ListOfBonds.size()+1]; | 
|---|
| [4455f4] | 230 | for (int i=ListOfBonds.size()+1;i--;) | 
|---|
|  | 231 | ComponentNr[i] = -1; | 
|---|
| [14b65e] | 232 | }; | 
|---|
|  | 233 |  | 
|---|
|  | 234 | void atom::resetGraphNr(){ | 
|---|
|  | 235 | GraphNr=-1; | 
|---|
|  | 236 | } | 
|---|
| [4455f4] | 237 |  | 
|---|
| [d74077] | 238 | std::ostream & atom::operator << (std::ostream &ost) const | 
|---|
|  | 239 | { | 
|---|
|  | 240 | ParticleInfo::operator<<(ost); | 
|---|
|  | 241 | ost << "," << getPosition(); | 
|---|
|  | 242 | return ost; | 
|---|
|  | 243 | } | 
|---|
|  | 244 |  | 
|---|
|  | 245 | std::ostream & operator << (std::ostream &ost, const atom &a) | 
|---|
|  | 246 | { | 
|---|
|  | 247 | a.ParticleInfo::operator<<(ost); | 
|---|
|  | 248 | ost << "," << a.getPosition(); | 
|---|
|  | 249 | return ost; | 
|---|
|  | 250 | } | 
|---|
| [4455f4] | 251 |  | 
|---|
|  | 252 | bool operator < (atom &a, atom &b) | 
|---|
|  | 253 | { | 
|---|
|  | 254 | return a.Compare(b); | 
|---|
|  | 255 | }; | 
|---|
|  | 256 |  | 
|---|
| [46d958] | 257 | World *atom::getWorld(){ | 
|---|
|  | 258 | return world; | 
|---|
|  | 259 | } | 
|---|
|  | 260 |  | 
|---|
|  | 261 | void atom::setWorld(World* _world){ | 
|---|
|  | 262 | world = _world; | 
|---|
|  | 263 | } | 
|---|
|  | 264 |  | 
|---|
| [88d586] | 265 | bool atom::changeId(atomId_t newId){ | 
|---|
|  | 266 | // first we move ourselves in the world | 
|---|
| [ceaab1] | 267 | // the world lets us know if that succeeded | 
|---|
|  | 268 | atomId_t oldid = id; | 
|---|
| [4f7f0bf] | 269 | if(world->changeAtomId(id,newId,this)){ | 
|---|
|  | 270 | OBSERVE; | 
|---|
|  | 271 | id = newId; | 
|---|
| [ceaab1] | 272 | if (mol != NULL) | 
|---|
|  | 273 | mol->changeAtomId(oldid, newId); | 
|---|
| [4f7f0bf] | 274 | NOTIFY(IndexChanged); | 
|---|
| [88d586] | 275 | return true; | 
|---|
|  | 276 | } | 
|---|
|  | 277 | else{ | 
|---|
|  | 278 | return false; | 
|---|
|  | 279 | } | 
|---|
|  | 280 | } | 
|---|
|  | 281 |  | 
|---|
|  | 282 | void atom::setId(atomId_t _id) { | 
|---|
| [46d958] | 283 | id=_id; | 
|---|
|  | 284 | } | 
|---|
|  | 285 |  | 
|---|
| [ad2b411] | 286 | atomId_t atom::getId() const { | 
|---|
| [46d958] | 287 | return id; | 
|---|
|  | 288 | } | 
|---|
|  | 289 |  | 
|---|
| [918a6c] | 290 | void atom::setMolecule(molecule *_mol){ | 
|---|
| [1187c5] | 291 | { | 
|---|
|  | 292 | // we must inform about molecule changed before the molecule itself tells about it | 
|---|
|  | 293 | OBSERVE; | 
|---|
|  | 294 | NOTIFY(MoleculeChanged); | 
|---|
|  | 295 | // take this atom from the old molecule | 
|---|
|  | 296 | removeFromMolecule(); | 
|---|
|  | 297 | mol = _mol; | 
|---|
|  | 298 | } | 
|---|
| [c32d21] | 299 | if ((mol) && (!mol->containsAtom(this))) { | 
|---|
| [c0f2fc] | 300 | mol->associateAtomWithMolecule(this); | 
|---|
| [c32d21] | 301 | } | 
|---|
| [6cfa36] | 302 | } | 
|---|
| [0d9546] | 303 |  | 
|---|
| [270bdf] | 304 | const molecule* atom::getMolecule() const { | 
|---|
| [c084cc] | 305 | return mol; | 
|---|
|  | 306 | } | 
|---|
|  | 307 |  | 
|---|
| [6cfa36] | 308 | void atom::removeFromMolecule(){ | 
|---|
|  | 309 | if(mol){ | 
|---|
| [1187c5] | 310 | OBSERVE; | 
|---|
|  | 311 | NOTIFY(MoleculeChanged); | 
|---|
| [c32d21] | 312 | if(mol->containsAtom(this)){ | 
|---|
| [c0f2fc] | 313 | mol->disassociateAtomWithMolecule(this); | 
|---|
| [6cfa36] | 314 | } | 
|---|
|  | 315 | mol=0; | 
|---|
|  | 316 | } | 
|---|
| [1f8337] | 317 | } | 
|---|
|  | 318 |  | 
|---|
| [560bbe] | 319 | bool atom::changeNr(const int newNr) | 
|---|
|  | 320 | { | 
|---|
|  | 321 | if ((mol) && (mol->changeAtomNr(getNr(),newNr,this))) { | 
|---|
|  | 322 | return true; | 
|---|
|  | 323 | } else{ | 
|---|
|  | 324 | return false; | 
|---|
|  | 325 | } | 
|---|
|  | 326 | } | 
|---|
|  | 327 |  | 
|---|
| [e8a21f] | 328 | int atom::getNr() const{ | 
|---|
| [735b1c] | 329 | return ParticleInfo::getNr(); | 
|---|
| [e8a21f] | 330 | } | 
|---|
| [6cfa36] | 331 |  | 
|---|
| [88d586] | 332 | atom* NewAtom(atomId_t _id){ | 
|---|
| [11f0fa] | 333 | atom * res = new atom(); | 
|---|
| [88d586] | 334 | res->setId(_id); | 
|---|
|  | 335 | return res; | 
|---|
| [46d958] | 336 | } | 
|---|
|  | 337 |  | 
|---|
| [88d586] | 338 | void DeleteAtom(atom* atom){ | 
|---|
| [46d958] | 339 | delete atom; | 
|---|
| [e5f64de] | 340 | } | 
|---|
|  | 341 |  | 
|---|
|  | 342 | bool compareAtomElements(atom* atom1,atom* atom2){ | 
|---|
| [ed26ae] | 343 | return atom1->getType()->getAtomicNumber() < atom2->getType()->getAtomicNumber(); | 
|---|
| [46d958] | 344 | } | 
|---|
| [d05088] | 345 | /* | 
|---|
|  | 346 | void atom::update(Observable *publisher) | 
|---|
|  | 347 | {} | 
|---|
|  | 348 |  | 
|---|
|  | 349 | void atom::recieveNotification(Observable *publisher, Notification_ptr notification) | 
|---|
|  | 350 | { | 
|---|
|  | 351 | ASSERT(0, "atom::recieveNotification() - we are not signed on to any notifications."); | 
|---|
|  | 352 | } | 
|---|
|  | 353 | */ | 
|---|
|  | 354 | void atom::subjectKilled(Observable *publisher) | 
|---|
|  | 355 | { | 
|---|
|  | 356 | // as publisher has been half-deallocated (Observable is one of the base classes, hence | 
|---|
|  | 357 | // becomes destroyed latest), we cannot senibly cast it anymore. | 
|---|
|  | 358 | // Hence, we simply have to check here whether it is NOT one of the other instances | 
|---|
|  | 359 | // we are signed on to. | 
|---|
|  | 360 | father = this; | 
|---|
|  | 361 | // no need to sign off | 
|---|
|  | 362 | } | 
|---|
| [b71881] | 363 |  | 
|---|
|  | 364 | void atom::select() | 
|---|
|  | 365 | { | 
|---|
|  | 366 | OBSERVE; | 
|---|
|  | 367 | selected = true; | 
|---|
|  | 368 | NOTIFY(SelectionChanged); | 
|---|
|  | 369 | } | 
|---|
|  | 370 |  | 
|---|
|  | 371 | void atom::unselect() | 
|---|
|  | 372 | { | 
|---|
|  | 373 | OBSERVE; | 
|---|
|  | 374 | selected = false; | 
|---|
|  | 375 | NOTIFY(SelectionChanged); | 
|---|
|  | 376 | } | 
|---|