| 1 | /* | 
|---|
| 2 | * Project: MoleCuilder | 
|---|
| 3 | * Description: creates and alters molecular systems | 
|---|
| 4 | * Copyright (C)  2010-2012 University of Bonn. All rights reserved. | 
|---|
| 5 | * Copyright (C)  2013 Frederik Heber. All rights reserved. | 
|---|
| 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/>. | 
|---|
| 22 | */ | 
|---|
| 23 |  | 
|---|
| 24 | /** \file atom.cpp | 
|---|
| 25 | * | 
|---|
| 26 | * Function implementations for the class atom. | 
|---|
| 27 | * | 
|---|
| 28 | */ | 
|---|
| 29 |  | 
|---|
| 30 | // include config.h | 
|---|
| 31 | #ifdef HAVE_CONFIG_H | 
|---|
| 32 | #include <config.h> | 
|---|
| 33 | #endif | 
|---|
| 34 |  | 
|---|
| 35 | #include "CodePatterns/MemDebug.hpp" | 
|---|
| 36 |  | 
|---|
| 37 | #include "atom.hpp" | 
|---|
| 38 | #include "AtomObserver.hpp" | 
|---|
| 39 | #include "Bond/bond.hpp" | 
|---|
| 40 | #include "CodePatterns/Log.hpp" | 
|---|
| 41 | #include "config.hpp" | 
|---|
| 42 | #include "Element/element.hpp" | 
|---|
| 43 | #include "LinearAlgebra/Vector.hpp" | 
|---|
| 44 | #include "World.hpp" | 
|---|
| 45 | #include "WorldTime.hpp" | 
|---|
| 46 | #include "molecule.hpp" | 
|---|
| 47 | #include "Shapes/Shape.hpp" | 
|---|
| 48 |  | 
|---|
| 49 | #include <iomanip> | 
|---|
| 50 | #include <iostream> | 
|---|
| 51 |  | 
|---|
| 52 | /************************************* Functions for class atom *************************************/ | 
|---|
| 53 |  | 
|---|
| 54 |  | 
|---|
| 55 | atom::atom() : | 
|---|
| 56 | father(this), | 
|---|
| 57 | sort(&Nr), | 
|---|
| 58 | mol(0), | 
|---|
| 59 | selected(false) | 
|---|
| 60 | { | 
|---|
| 61 | // note AtomObserver about inserted atom | 
|---|
| 62 | AtomObserver::getInstance().Inserted(this); | 
|---|
| 63 | } | 
|---|
| 64 |  | 
|---|
| 65 | atom::atom(atom *pointer) : | 
|---|
| 66 | ParticleInfo(*pointer), | 
|---|
| 67 | AtomInfo(*pointer), | 
|---|
| 68 | father(pointer), | 
|---|
| 69 | sort(&Nr), | 
|---|
| 70 | mol(0), | 
|---|
| 71 | selected(false) | 
|---|
| 72 | { | 
|---|
| 73 | // sign on to father atom to be notified when it is removed | 
|---|
| 74 | father->signOn(this); | 
|---|
| 75 |  | 
|---|
| 76 | // note AtomObserver about inserted atom | 
|---|
| 77 | AtomObserver::getInstance().Inserted(this); | 
|---|
| 78 | }; | 
|---|
| 79 |  | 
|---|
| 80 | atom *atom::clone(){ | 
|---|
| 81 | atom *res = new atom(this); | 
|---|
| 82 | World::getInstance().registerAtom(res); | 
|---|
| 83 | return res; | 
|---|
| 84 | } | 
|---|
| 85 |  | 
|---|
| 86 |  | 
|---|
| 87 | /** Destructor of class atom. | 
|---|
| 88 | */ | 
|---|
| 89 | atom::~atom() | 
|---|
| 90 | { | 
|---|
| 91 | // sign off from possible father | 
|---|
| 92 | if ((father != this) && (father != NULL)) | 
|---|
| 93 | father->signOff(this); | 
|---|
| 94 |  | 
|---|
| 95 | removeFromMolecule(); | 
|---|
| 96 | // note AtomObserver about removed atom | 
|---|
| 97 | AtomObserver::getInstance().Removed(this); | 
|---|
| 98 | } | 
|---|
| 99 |  | 
|---|
| 100 |  | 
|---|
| 101 | void atom::UpdateStep(const unsigned int _step) | 
|---|
| 102 | { | 
|---|
| 103 | LOG(4,"atom::UpdateStep() called."); | 
|---|
| 104 | // append to position, velocity and force vector | 
|---|
| 105 | AtomInfo::AppendTrajectoryStep(WorldTime::getTime()+1); | 
|---|
| 106 | // append to ListOfBonds vector | 
|---|
| 107 | BondedParticleInfo::AppendTrajectoryStep(WorldTime::getTime()+1); | 
|---|
| 108 | } | 
|---|
| 109 |  | 
|---|
| 110 | void atom::removeStep(const unsigned int _step) | 
|---|
| 111 | { | 
|---|
| 112 | LOG(4,"atom::removeStep() called."); | 
|---|
| 113 | // append to position, velocity and force vector | 
|---|
| 114 | AtomInfo::removeTrajectoryStep(_step); | 
|---|
| 115 | // append to ListOfBonds vector | 
|---|
| 116 | BondedParticleInfo::removeTrajectoryStep(_step); | 
|---|
| 117 | } | 
|---|
| 118 |  | 
|---|
| 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 | 
|---|
| 126 | { | 
|---|
| 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 | } | 
|---|
| 136 | } | 
|---|
| 137 |  | 
|---|
| 138 | void atom::setFather(atom * const _father) | 
|---|
| 139 | { | 
|---|
| 140 | // sign off from old father | 
|---|
| 141 | if ((father != this) && (father != NULL)) | 
|---|
| 142 | father->signOff(this); | 
|---|
| 143 |  | 
|---|
| 144 | father = _father; | 
|---|
| 145 | father->signOn(this); | 
|---|
| 146 | } | 
|---|
| 147 |  | 
|---|
| 148 | /** Sets father to itself or its father in case of copying a molecule. | 
|---|
| 149 | */ | 
|---|
| 150 | void atom::CorrectFather() | 
|---|
| 151 | { | 
|---|
| 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 | 
|---|
| 155 | father = father->father;  // set father to original's father | 
|---|
| 156 |  | 
|---|
| 157 | }; | 
|---|
| 158 |  | 
|---|
| 159 | void atom::EqualsFather ( const atom *ptr, const atom **res ) const | 
|---|
| 160 | { | 
|---|
| 161 | if ( ptr == father ) | 
|---|
| 162 | *res = this; | 
|---|
| 163 | }; | 
|---|
| 164 |  | 
|---|
| 165 | bool atom::isFather(const atom *ptr){ | 
|---|
| 166 | return ptr==father; | 
|---|
| 167 | } | 
|---|
| 168 |  | 
|---|
| 169 | bool atom::OutputIndexed(ofstream * const out, const int ElementNo, const int AtomNo, const char *comment) const | 
|---|
| 170 | { | 
|---|
| 171 | if (out != NULL) { | 
|---|
| 172 | *out << "Ion_Type" << ElementNo << "_" << AtomNo << "\t"  << fixed << setprecision(9) << showpoint; | 
|---|
| 173 | *out << at(0) << "\t" << at(1) << "\t" << at(2); | 
|---|
| 174 | *out << "\t" << (int)(getFixedIon()); | 
|---|
| 175 | if (getAtomicVelocity().Norm() > MYEPSILON) | 
|---|
| 176 | *out << "\t" << scientific << setprecision(6) << getAtomicVelocity()[0] << "\t" << getAtomicVelocity()[1] << "\t" << getAtomicVelocity()[2] << "\t"; | 
|---|
| 177 | if (comment != NULL) | 
|---|
| 178 | *out << " # " << comment << endl; | 
|---|
| 179 | else | 
|---|
| 180 | *out << " # molecule nr " << getNr() << endl; | 
|---|
| 181 | return true; | 
|---|
| 182 | } else | 
|---|
| 183 | return false; | 
|---|
| 184 | }; | 
|---|
| 185 |  | 
|---|
| 186 | bool atom::OutputArrayIndexed(ostream * const out,const enumeration<const element*> &elementLookup, int *AtomNo, const char *comment) const | 
|---|
| 187 | { | 
|---|
| 188 | AtomNo[getType()->getAtomicNumber()]++;  // increment number | 
|---|
| 189 | if (out != NULL) { | 
|---|
| 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"); | 
|---|
| 192 | *out << "Ion_Type" << elementLookup.there.find(elemental)->second << "_" << AtomNo[elemental->getAtomicNumber()] << "\t"  << fixed << setprecision(9) << showpoint; | 
|---|
| 193 | *out << at(0) << "\t" << at(1) << "\t" << at(2); | 
|---|
| 194 | *out << "\t" << getFixedIon(); | 
|---|
| 195 | if (getAtomicVelocity().Norm() > MYEPSILON) | 
|---|
| 196 | *out << "\t" << scientific << setprecision(6) << getAtomicVelocity()[0] << "\t" << getAtomicVelocity()[1] << "\t" << getAtomicVelocity()[2] << "\t"; | 
|---|
| 197 | if (comment != NULL) | 
|---|
| 198 | *out << " # " << comment << endl; | 
|---|
| 199 | else | 
|---|
| 200 | *out << " # molecule nr " << getNr() << endl; | 
|---|
| 201 | return true; | 
|---|
| 202 | } else | 
|---|
| 203 | return false; | 
|---|
| 204 | }; | 
|---|
| 205 |  | 
|---|
| 206 | bool atom::Compare(const atom &ptr) const | 
|---|
| 207 | { | 
|---|
| 208 | if (getNr() < ptr.getNr()) | 
|---|
| 209 | return true; | 
|---|
| 210 | else | 
|---|
| 211 | return false; | 
|---|
| 212 | }; | 
|---|
| 213 |  | 
|---|
| 214 | double atom::DistanceSquaredToVector(const Vector &origin) const | 
|---|
| 215 | { | 
|---|
| 216 | return DistanceSquared(origin); | 
|---|
| 217 | }; | 
|---|
| 218 |  | 
|---|
| 219 | double atom::DistanceToVector(const Vector &origin) const | 
|---|
| 220 | { | 
|---|
| 221 | return distance(origin); | 
|---|
| 222 | }; | 
|---|
| 223 |  | 
|---|
| 224 | void atom::InitComponentNr() | 
|---|
| 225 | { | 
|---|
| 226 | if (ComponentNr != NULL) | 
|---|
| 227 | delete[](ComponentNr); | 
|---|
| 228 | const BondList& ListOfBonds = getListOfBonds(); | 
|---|
| 229 | ComponentNr = new int[ListOfBonds.size()+1]; | 
|---|
| 230 | for (int i=ListOfBonds.size()+1;i--;) | 
|---|
| 231 | ComponentNr[i] = -1; | 
|---|
| 232 | }; | 
|---|
| 233 |  | 
|---|
| 234 | void atom::resetGraphNr(){ | 
|---|
| 235 | GraphNr=-1; | 
|---|
| 236 | } | 
|---|
| 237 |  | 
|---|
| 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 | } | 
|---|
| 251 |  | 
|---|
| 252 | bool operator < (atom &a, atom &b) | 
|---|
| 253 | { | 
|---|
| 254 | return a.Compare(b); | 
|---|
| 255 | }; | 
|---|
| 256 |  | 
|---|
| 257 | World *atom::getWorld(){ | 
|---|
| 258 | return world; | 
|---|
| 259 | } | 
|---|
| 260 |  | 
|---|
| 261 | void atom::setWorld(World* _world){ | 
|---|
| 262 | world = _world; | 
|---|
| 263 | } | 
|---|
| 264 |  | 
|---|
| 265 | bool atom::changeId(atomId_t newId){ | 
|---|
| 266 | // first we move ourselves in the world | 
|---|
| 267 | // the world lets us know if that succeeded | 
|---|
| 268 | atomId_t oldid = id; | 
|---|
| 269 | if(world->changeAtomId(id,newId,this)){ | 
|---|
| 270 | OBSERVE; | 
|---|
| 271 | id = newId; | 
|---|
| 272 | if (mol != NULL) | 
|---|
| 273 | mol->changeAtomId(oldid, newId); | 
|---|
| 274 | NOTIFY(IndexChanged); | 
|---|
| 275 | return true; | 
|---|
| 276 | } | 
|---|
| 277 | else{ | 
|---|
| 278 | return false; | 
|---|
| 279 | } | 
|---|
| 280 | } | 
|---|
| 281 |  | 
|---|
| 282 | void atom::setId(atomId_t _id) { | 
|---|
| 283 | id=_id; | 
|---|
| 284 | } | 
|---|
| 285 |  | 
|---|
| 286 | atomId_t atom::getId() const { | 
|---|
| 287 | return id; | 
|---|
| 288 | } | 
|---|
| 289 |  | 
|---|
| 290 | void atom::setMolecule(molecule *_mol){ | 
|---|
| 291 | OBSERVE; | 
|---|
| 292 | NOTIFY(MoleculeChanged); | 
|---|
| 293 | // take this atom from the old molecule | 
|---|
| 294 | removeFromMolecule(); | 
|---|
| 295 | mol = _mol; | 
|---|
| 296 | if ((mol) && (!mol->containsAtom(this))) { | 
|---|
| 297 | mol->associateAtomWithMolecule(this); | 
|---|
| 298 | } | 
|---|
| 299 | } | 
|---|
| 300 |  | 
|---|
| 301 | const molecule* atom::getMolecule() const { | 
|---|
| 302 | return mol; | 
|---|
| 303 | } | 
|---|
| 304 |  | 
|---|
| 305 | void atom::removeFromMolecule(){ | 
|---|
| 306 | if(mol){ | 
|---|
| 307 | if(mol->containsAtom(this)){ | 
|---|
| 308 | mol->disassociateAtomWithMolecule(this); | 
|---|
| 309 | } | 
|---|
| 310 | mol=0; | 
|---|
| 311 | } | 
|---|
| 312 | } | 
|---|
| 313 |  | 
|---|
| 314 | bool atom::changeNr(const int newNr) | 
|---|
| 315 | { | 
|---|
| 316 | if ((mol) && (mol->changeAtomNr(getNr(),newNr,this))) { | 
|---|
| 317 | return true; | 
|---|
| 318 | } else{ | 
|---|
| 319 | return false; | 
|---|
| 320 | } | 
|---|
| 321 | } | 
|---|
| 322 |  | 
|---|
| 323 | int atom::getNr() const{ | 
|---|
| 324 | return ParticleInfo::getNr(); | 
|---|
| 325 | } | 
|---|
| 326 |  | 
|---|
| 327 | atom* NewAtom(atomId_t _id){ | 
|---|
| 328 | atom * res = new atom(); | 
|---|
| 329 | res->setId(_id); | 
|---|
| 330 | return res; | 
|---|
| 331 | } | 
|---|
| 332 |  | 
|---|
| 333 | void DeleteAtom(atom* atom){ | 
|---|
| 334 | delete atom; | 
|---|
| 335 | } | 
|---|
| 336 |  | 
|---|
| 337 | bool compareAtomElements(atom* atom1,atom* atom2){ | 
|---|
| 338 | return atom1->getType()->getAtomicNumber() < atom2->getType()->getAtomicNumber(); | 
|---|
| 339 | } | 
|---|
| 340 | /* | 
|---|
| 341 | void atom::update(Observable *publisher) | 
|---|
| 342 | {} | 
|---|
| 343 |  | 
|---|
| 344 | void atom::recieveNotification(Observable *publisher, Notification_ptr notification) | 
|---|
| 345 | { | 
|---|
| 346 | ASSERT(0, "atom::recieveNotification() - we are not signed on to any notifications."); | 
|---|
| 347 | } | 
|---|
| 348 | */ | 
|---|
| 349 | void atom::subjectKilled(Observable *publisher) | 
|---|
| 350 | { | 
|---|
| 351 | // as publisher has been half-deallocated (Observable is one of the base classes, hence | 
|---|
| 352 | // becomes destroyed latest), we cannot senibly cast it anymore. | 
|---|
| 353 | // Hence, we simply have to check here whether it is NOT one of the other instances | 
|---|
| 354 | // we are signed on to. | 
|---|
| 355 | father = this; | 
|---|
| 356 | // no need to sign off | 
|---|
| 357 | } | 
|---|
| 358 |  | 
|---|
| 359 | void atom::select() | 
|---|
| 360 | { | 
|---|
| 361 | OBSERVE; | 
|---|
| 362 | selected = true; | 
|---|
| 363 | NOTIFY(SelectionChanged); | 
|---|
| 364 | } | 
|---|
| 365 |  | 
|---|
| 366 | void atom::unselect() | 
|---|
| 367 | { | 
|---|
| 368 | OBSERVE; | 
|---|
| 369 | selected = false; | 
|---|
| 370 | NOTIFY(SelectionChanged); | 
|---|
| 371 | } | 
|---|