| 1 | /* | 
| 2 | * Copyright (c) 2005 The University of Notre Dame. All Rights Reserved. | 
| 3 | * | 
| 4 | * The University of Notre Dame grants you ("Licensee") a | 
| 5 | * non-exclusive, royalty free, license to use, modify and | 
| 6 | * redistribute this software in source and binary code form, provided | 
| 7 | * that the following conditions are met: | 
| 8 | * | 
| 9 | * 1. Redistributions of source code must retain the above copyright | 
| 10 | *    notice, this list of conditions and the following disclaimer. | 
| 11 | * | 
| 12 | * 2. Redistributions in binary form must reproduce the above copyright | 
| 13 | *    notice, this list of conditions and the following disclaimer in the | 
| 14 | *    documentation and/or other materials provided with the | 
| 15 | *    distribution. | 
| 16 | * | 
| 17 | * This software is provided "AS IS," without a warranty of any | 
| 18 | * kind. All express or implied conditions, representations and | 
| 19 | * warranties, including any implied warranty of merchantability, | 
| 20 | * fitness for a particular purpose or non-infringement, are hereby | 
| 21 | * excluded.  The University of Notre Dame and its licensors shall not | 
| 22 | * be liable for any damages suffered by licensee as a result of | 
| 23 | * using, modifying or distributing the software or its | 
| 24 | * derivatives. In no event will the University of Notre Dame or its | 
| 25 | * licensors be liable for any lost revenue, profit or data, or for | 
| 26 | * direct, indirect, special, consequential, incidental or punitive | 
| 27 | * damages, however caused and regardless of the theory of liability, | 
| 28 | * arising out of the use of or inability to use software, even if the | 
| 29 | * University of Notre Dame has been advised of the possibility of | 
| 30 | * such damages. | 
| 31 | * | 
| 32 | * SUPPORT OPEN SCIENCE!  If you use OpenMD or its source code in your | 
| 33 | * research, please cite the appropriate papers when you publish your | 
| 34 | * work.  Good starting points are: | 
| 35 | * | 
| 36 | * [1]  Meineke, et al., J. Comp. Chem. 26, 252-271 (2005). | 
| 37 | * [2]  Fennell & Gezelter, J. Chem. Phys. 124, 234104 (2006). | 
| 38 | * [3]  Sun, Lin & Gezelter, J. Chem. Phys. 128, 24107 (2008). | 
| 39 | * [4]  Vardeman & Gezelter, in progress (2009). | 
| 40 | */ | 
| 41 |  | 
| 42 | /** | 
| 43 | * @file MoleculeCreator.cpp | 
| 44 | * @author tlin | 
| 45 | * @date 11/04/2004 | 
| 46 | * @time 13:44am | 
| 47 | * @version 1.0 | 
| 48 | */ | 
| 49 |  | 
| 50 | #include <cassert> | 
| 51 | #include <typeinfo> | 
| 52 | #include <set> | 
| 53 |  | 
| 54 | #include "brains/MoleculeCreator.hpp" | 
| 55 | #include "primitives/GhostBend.hpp" | 
| 56 | #include "primitives/GhostTorsion.hpp" | 
| 57 | #include "types/DirectionalAtomType.hpp" | 
| 58 | #include "types/FixedBondType.hpp" | 
| 59 | #include "utils/simError.h" | 
| 60 | #include "utils/StringUtils.hpp" | 
| 61 |  | 
| 62 | namespace OpenMD { | 
| 63 |  | 
| 64 | Molecule* MoleculeCreator::createMolecule(ForceField* ff, | 
| 65 | MoleculeStamp *molStamp, | 
| 66 | int stampId, int globalIndex, | 
| 67 | LocalIndexManager* localIndexMan) { | 
| 68 | Molecule* mol = new Molecule(stampId, globalIndex, molStamp->getName()); | 
| 69 |  | 
| 70 | //create atoms | 
| 71 | Atom* atom; | 
| 72 | AtomStamp* currentAtomStamp; | 
| 73 | int nAtom = molStamp->getNAtoms(); | 
| 74 | for (int i = 0; i < nAtom; ++i) { | 
| 75 | currentAtomStamp = molStamp->getAtomStamp(i); | 
| 76 | atom = createAtom(ff, mol, currentAtomStamp, localIndexMan); | 
| 77 | mol->addAtom(atom); | 
| 78 | } | 
| 79 |  | 
| 80 | //create rigidbodies | 
| 81 | RigidBody* rb; | 
| 82 | RigidBodyStamp * currentRigidBodyStamp; | 
| 83 | int nRigidbodies = molStamp->getNRigidBodies(); | 
| 84 |  | 
| 85 | for (int i = 0; i < nRigidbodies; ++i) { | 
| 86 | currentRigidBodyStamp = molStamp->getRigidBodyStamp(i); | 
| 87 | rb = createRigidBody(molStamp, mol, currentRigidBodyStamp, | 
| 88 | localIndexMan); | 
| 89 | mol->addRigidBody(rb); | 
| 90 | } | 
| 91 |  | 
| 92 | //create bonds | 
| 93 | Bond* bond; | 
| 94 | BondStamp* currentBondStamp; | 
| 95 | int nBonds = molStamp->getNBonds(); | 
| 96 |  | 
| 97 | for (int i = 0; i < nBonds; ++i) { | 
| 98 | currentBondStamp = molStamp->getBondStamp(i); | 
| 99 | bond = createBond(ff, mol, currentBondStamp); | 
| 100 | mol->addBond(bond); | 
| 101 | } | 
| 102 |  | 
| 103 | //create bends | 
| 104 | Bend* bend; | 
| 105 | BendStamp* currentBendStamp; | 
| 106 | int nBends = molStamp->getNBends(); | 
| 107 | for (int i = 0; i < nBends; ++i) { | 
| 108 | currentBendStamp = molStamp->getBendStamp(i); | 
| 109 | bend = createBend(ff, mol, currentBendStamp); | 
| 110 | mol->addBend(bend); | 
| 111 | } | 
| 112 |  | 
| 113 | //create torsions | 
| 114 | Torsion* torsion; | 
| 115 | TorsionStamp* currentTorsionStamp; | 
| 116 | int nTorsions = molStamp->getNTorsions(); | 
| 117 | for (int i = 0; i < nTorsions; ++i) { | 
| 118 | currentTorsionStamp = molStamp->getTorsionStamp(i); | 
| 119 | torsion = createTorsion(ff, mol, currentTorsionStamp); | 
| 120 | mol->addTorsion(torsion); | 
| 121 | } | 
| 122 |  | 
| 123 | //create inversions | 
| 124 | Inversion* inversion; | 
| 125 | InversionStamp* currentInversionStamp; | 
| 126 | int nInversions = molStamp->getNInversions(); | 
| 127 | for (int i = 0; i < nInversions; ++i) { | 
| 128 | currentInversionStamp = molStamp->getInversionStamp(i); | 
| 129 | inversion = createInversion(ff, mol, currentInversionStamp); | 
| 130 | if (inversion != NULL ) { | 
| 131 | mol->addInversion(inversion); | 
| 132 | } | 
| 133 | } | 
| 134 |  | 
| 135 | //create cutoffGroups | 
| 136 | CutoffGroup* cutoffGroup; | 
| 137 | CutoffGroupStamp* currentCutoffGroupStamp; | 
| 138 | int nCutoffGroups = molStamp->getNCutoffGroups(); | 
| 139 | for (int i = 0; i < nCutoffGroups; ++i) { | 
| 140 | currentCutoffGroupStamp = molStamp->getCutoffGroupStamp(i); | 
| 141 | cutoffGroup = createCutoffGroup(mol, currentCutoffGroupStamp); | 
| 142 | mol->addCutoffGroup(cutoffGroup); | 
| 143 | } | 
| 144 |  | 
| 145 | //every free atom is a cutoff group | 
| 146 | std::vector<Atom*> freeAtoms; | 
| 147 | std::vector<Atom*>::iterator ai; | 
| 148 | std::vector<Atom*>::iterator fai; | 
| 149 |  | 
| 150 | //add all atoms into allAtoms set | 
| 151 | for(atom = mol->beginAtom(fai); atom != NULL; atom = mol->nextAtom(fai)) { | 
| 152 | freeAtoms.push_back(atom); | 
| 153 | } | 
| 154 |  | 
| 155 | Molecule::CutoffGroupIterator ci; | 
| 156 | CutoffGroup* cg; | 
| 157 |  | 
| 158 | for (cg = mol->beginCutoffGroup(ci); cg != NULL; | 
| 159 | cg = mol->nextCutoffGroup(ci)) { | 
| 160 |  | 
| 161 | for(atom = cg->beginAtom(ai); atom != NULL; atom = cg->nextAtom(ai)) { | 
| 162 | //erase the atoms belong to cutoff groups from freeAtoms vector | 
| 163 | freeAtoms.erase(std::remove(freeAtoms.begin(), freeAtoms.end(), atom), | 
| 164 | freeAtoms.end()); | 
| 165 | } | 
| 166 | } | 
| 167 |  | 
| 168 | // loop over the free atoms and then create one cutoff group for | 
| 169 | // every single free atom | 
| 170 |  | 
| 171 | for (fai = freeAtoms.begin(); fai != freeAtoms.end(); ++fai) { | 
| 172 | cutoffGroup = createCutoffGroup(mol, *fai); | 
| 173 | mol->addCutoffGroup(cutoffGroup); | 
| 174 | } | 
| 175 | //create constraints | 
| 176 | createConstraintPair(mol); | 
| 177 | createConstraintElem(mol); | 
| 178 |  | 
| 179 | //the construction of this molecule is finished | 
| 180 | mol->complete(); | 
| 181 |  | 
| 182 | return mol; | 
| 183 | } | 
| 184 |  | 
| 185 |  | 
| 186 | Atom* MoleculeCreator::createAtom(ForceField* ff, Molecule* mol, | 
| 187 | AtomStamp* stamp, | 
| 188 | LocalIndexManager* localIndexMan) { | 
| 189 | AtomType * atomType; | 
| 190 | Atom* atom; | 
| 191 |  | 
| 192 | atomType =  ff->getAtomType(stamp->getType()); | 
| 193 |  | 
| 194 | if (atomType == NULL) { | 
| 195 | sprintf(painCave.errMsg, "Can not find Matching Atom Type for[%s]", | 
| 196 | stamp->getType().c_str()); | 
| 197 |  | 
| 198 | painCave.isFatal = 1; | 
| 199 | simError(); | 
| 200 | } | 
| 201 |  | 
| 202 | //below code still have some kind of hard-coding smell | 
| 203 | if (atomType->isDirectional()){ | 
| 204 |  | 
| 205 | DirectionalAtomType* dAtomType = dynamic_cast<DirectionalAtomType*>(atomType); | 
| 206 |  | 
| 207 | if (dAtomType == NULL) { | 
| 208 | sprintf(painCave.errMsg, "Can not cast AtomType to DirectionalAtomType"); | 
| 209 |  | 
| 210 | painCave.isFatal = 1; | 
| 211 | simError(); | 
| 212 | } | 
| 213 |  | 
| 214 | DirectionalAtom* dAtom; | 
| 215 | dAtom = new DirectionalAtom(dAtomType); | 
| 216 | atom = dAtom; | 
| 217 | } | 
| 218 | else{ | 
| 219 | atom = new Atom(atomType); | 
| 220 | } | 
| 221 |  | 
| 222 | atom->setLocalIndex(localIndexMan->getNextAtomIndex()); | 
| 223 |  | 
| 224 | return atom; | 
| 225 | } | 
| 226 |  | 
| 227 | RigidBody* MoleculeCreator::createRigidBody(MoleculeStamp *molStamp, | 
| 228 | Molecule* mol, | 
| 229 | RigidBodyStamp* rbStamp, | 
| 230 | LocalIndexManager* localIndexMan) { | 
| 231 | Atom* atom; | 
| 232 | int nAtoms; | 
| 233 | Vector3d refCoor; | 
| 234 | AtomStamp* atomStamp; | 
| 235 |  | 
| 236 | RigidBody* rb = new RigidBody(); | 
| 237 | nAtoms = rbStamp->getNMembers(); | 
| 238 | for (int i = 0; i < nAtoms; ++i) { | 
| 239 | //rbStamp->getMember(i) return the local index of current atom | 
| 240 | //inside the molecule.  It is not the same as local index of | 
| 241 | //atom which is the index of atom at DataStorage class | 
| 242 | atom = mol->getAtomAt(rbStamp->getMemberAt(i)); | 
| 243 | atomStamp= molStamp->getAtomStamp(rbStamp->getMemberAt(i)); | 
| 244 | rb->addAtom(atom, atomStamp); | 
| 245 | } | 
| 246 |  | 
| 247 | //after all of the atoms are added, we need to calculate the | 
| 248 | //reference coordinates | 
| 249 | rb->calcRefCoords(); | 
| 250 |  | 
| 251 | //set the local index of this rigid body, global index will be set later | 
| 252 | rb->setLocalIndex(localIndexMan->getNextRigidBodyIndex()); | 
| 253 |  | 
| 254 | //the rule for naming rigidbody MoleculeName_RB_Integer | 
| 255 | //The first part is the name of the molecule | 
| 256 | //The second part is alway fixed as "RB" | 
| 257 | //The third part is the index of the rigidbody defined in meta-data file | 
| 258 | //For example, Butane_RB_0 is a valid rigid body name of butane molecule | 
| 259 | /**@todo replace itoa by lexi_cast */ | 
| 260 | std::string s = OpenMD_itoa(mol->getNRigidBodies(), 10); | 
| 261 | rb->setType(mol->getType() + "_RB_" + s.c_str()); | 
| 262 |  | 
| 263 | return rb; | 
| 264 | } | 
| 265 |  | 
| 266 | Bond* MoleculeCreator::createBond(ForceField* ff, Molecule* mol, | 
| 267 | BondStamp* stamp) { | 
| 268 | BondType* bondType; | 
| 269 | Atom* atomA; | 
| 270 | Atom* atomB; | 
| 271 |  | 
| 272 | atomA = mol->getAtomAt(stamp->getA()); | 
| 273 | atomB = mol->getAtomAt(stamp->getB()); | 
| 274 |  | 
| 275 | assert( atomA && atomB); | 
| 276 |  | 
| 277 | bondType = ff->getBondType(atomA->getType(), atomB->getType()); | 
| 278 |  | 
| 279 | if (bondType == NULL) { | 
| 280 | sprintf(painCave.errMsg, "Can not find Matching Bond Type for[%s, %s]", | 
| 281 | atomA->getType().c_str(), | 
| 282 | atomB->getType().c_str()); | 
| 283 |  | 
| 284 | painCave.isFatal = 1; | 
| 285 | simError(); | 
| 286 | } | 
| 287 | return new Bond(atomA, atomB, bondType); | 
| 288 | } | 
| 289 |  | 
| 290 | Bend* MoleculeCreator::createBend(ForceField* ff, Molecule* mol, | 
| 291 | BendStamp* stamp) { | 
| 292 | Bend* bend = NULL; | 
| 293 | std::vector<int> bendAtoms = stamp->getMembers(); | 
| 294 | if (bendAtoms.size() == 3) { | 
| 295 | Atom* atomA = mol->getAtomAt(bendAtoms[0]); | 
| 296 | Atom* atomB = mol->getAtomAt(bendAtoms[1]); | 
| 297 | Atom* atomC = mol->getAtomAt(bendAtoms[2]); | 
| 298 |  | 
| 299 | assert( atomA && atomB && atomC); | 
| 300 |  | 
| 301 | BendType* bendType = ff->getBendType(atomA->getType().c_str(), | 
| 302 | atomB->getType().c_str(), | 
| 303 | atomC->getType().c_str()); | 
| 304 |  | 
| 305 | if (bendType == NULL) { | 
| 306 | sprintf(painCave.errMsg, "Can not find Matching Bend Type for[%s, %s, %s]", | 
| 307 | atomA->getType().c_str(), | 
| 308 | atomB->getType().c_str(), | 
| 309 | atomC->getType().c_str()); | 
| 310 |  | 
| 311 | painCave.isFatal = 1; | 
| 312 | simError(); | 
| 313 | } | 
| 314 |  | 
| 315 | bend = new Bend(atomA, atomB, atomC, bendType); | 
| 316 | } else if ( bendAtoms.size() == 2 && stamp->haveGhostVectorSource()) { | 
| 317 | int ghostIndex = stamp->getGhostVectorSource(); | 
| 318 | int normalIndex = ghostIndex != bendAtoms[0] ? bendAtoms[0] : bendAtoms[1]; | 
| 319 | Atom* normalAtom = mol->getAtomAt(normalIndex) ; | 
| 320 | DirectionalAtom* ghostAtom = dynamic_cast<DirectionalAtom*>(mol->getAtomAt(ghostIndex)); | 
| 321 | if (ghostAtom == NULL) { | 
| 322 | sprintf(painCave.errMsg, "Can not cast Atom to DirectionalAtom"); | 
| 323 | painCave.isFatal = 1; | 
| 324 | simError(); | 
| 325 | } | 
| 326 |  | 
| 327 | BendType* bendType = ff->getBendType(normalAtom->getType(), ghostAtom->getType(), "GHOST"); | 
| 328 |  | 
| 329 | if (bendType == NULL) { | 
| 330 | sprintf(painCave.errMsg, "Can not find Matching Bend Type for[%s, %s, %s]", | 
| 331 | normalAtom->getType().c_str(), | 
| 332 | ghostAtom->getType().c_str(), | 
| 333 | "GHOST"); | 
| 334 |  | 
| 335 | painCave.isFatal = 1; | 
| 336 | simError(); | 
| 337 | } | 
| 338 |  | 
| 339 | bend = new GhostBend(normalAtom, ghostAtom, bendType); | 
| 340 |  | 
| 341 | } | 
| 342 |  | 
| 343 | return bend; | 
| 344 | } | 
| 345 |  | 
| 346 | Torsion* MoleculeCreator::createTorsion(ForceField* ff, Molecule* mol, | 
| 347 | TorsionStamp* stamp) { | 
| 348 |  | 
| 349 | Torsion* torsion = NULL; | 
| 350 | std::vector<int> torsionAtoms = stamp->getMembers(); | 
| 351 | if (torsionAtoms.size() < 3) { | 
| 352 | return torsion; | 
| 353 | } | 
| 354 |  | 
| 355 | Atom* atomA = mol->getAtomAt(torsionAtoms[0]); | 
| 356 | Atom* atomB = mol->getAtomAt(torsionAtoms[1]); | 
| 357 | Atom* atomC = mol->getAtomAt(torsionAtoms[2]); | 
| 358 |  | 
| 359 | if (torsionAtoms.size() == 4) { | 
| 360 | Atom* atomD = mol->getAtomAt(torsionAtoms[3]); | 
| 361 |  | 
| 362 | assert(atomA && atomB && atomC && atomD); | 
| 363 |  | 
| 364 | TorsionType* torsionType = ff->getTorsionType(atomA->getType(), | 
| 365 | atomB->getType(), | 
| 366 | atomC->getType(), | 
| 367 | atomD->getType()); | 
| 368 | if (torsionType == NULL) { | 
| 369 | sprintf(painCave.errMsg, "Can not find Matching Torsion Type for[%s, %s, %s, %s]", | 
| 370 | atomA->getType().c_str(), | 
| 371 | atomB->getType().c_str(), | 
| 372 | atomC->getType().c_str(), | 
| 373 | atomD->getType().c_str()); | 
| 374 |  | 
| 375 | painCave.isFatal = 1; | 
| 376 | simError(); | 
| 377 | } | 
| 378 |  | 
| 379 | torsion = new Torsion(atomA, atomB, atomC, atomD, torsionType); | 
| 380 | } | 
| 381 | else { | 
| 382 |  | 
| 383 | DirectionalAtom* dAtom = dynamic_cast<DirectionalAtom*>(mol->getAtomAt(stamp->getGhostVectorSource())); | 
| 384 | if (dAtom == NULL) { | 
| 385 | sprintf(painCave.errMsg, "Can not cast Atom to DirectionalAtom"); | 
| 386 | painCave.isFatal = 1; | 
| 387 | simError(); | 
| 388 | } | 
| 389 |  | 
| 390 | TorsionType* torsionType = ff->getTorsionType(atomA->getType(), atomB->getType(), | 
| 391 | atomC->getType(), "GHOST"); | 
| 392 |  | 
| 393 | if (torsionType == NULL) { | 
| 394 | sprintf(painCave.errMsg, "Can not find Matching Torsion Type for[%s, %s, %s, %s]", | 
| 395 | atomA->getType().c_str(), | 
| 396 | atomB->getType().c_str(), | 
| 397 | atomC->getType().c_str(), | 
| 398 | "GHOST"); | 
| 399 |  | 
| 400 | painCave.isFatal = 1; | 
| 401 | simError(); | 
| 402 | } | 
| 403 |  | 
| 404 | torsion = new GhostTorsion(atomA, atomB, dAtom, torsionType); | 
| 405 | } | 
| 406 |  | 
| 407 | return torsion; | 
| 408 | } | 
| 409 |  | 
| 410 | Inversion* MoleculeCreator::createInversion(ForceField* ff, Molecule* mol, | 
| 411 | InversionStamp* stamp) { | 
| 412 |  | 
| 413 | Inversion* inversion = NULL; | 
| 414 | int center = stamp->getCenter(); | 
| 415 | std::vector<int> satellites = stamp->getSatellites(); | 
| 416 | if (satellites.size() != 3) { | 
| 417 | return inversion; | 
| 418 | } | 
| 419 |  | 
| 420 | Atom* atomA = mol->getAtomAt(center); | 
| 421 | Atom* atomB = mol->getAtomAt(satellites[0]); | 
| 422 | Atom* atomC = mol->getAtomAt(satellites[1]); | 
| 423 | Atom* atomD = mol->getAtomAt(satellites[2]); | 
| 424 |  | 
| 425 | assert(atomA && atomB && atomC && atomD); | 
| 426 |  | 
| 427 | InversionType* inversionType = ff->getInversionType(atomA->getType(), | 
| 428 | atomB->getType(), | 
| 429 | atomC->getType(), | 
| 430 | atomD->getType()); | 
| 431 |  | 
| 432 | if (inversionType == NULL) { | 
| 433 | sprintf(painCave.errMsg, "No Matching Inversion Type for[%s, %s, %s, %s]\n" | 
| 434 | "\t(May not be a problem: not all inversions are parametrized)\n", | 
| 435 | atomA->getType().c_str(), | 
| 436 | atomB->getType().c_str(), | 
| 437 | atomC->getType().c_str(), | 
| 438 | atomD->getType().c_str()); | 
| 439 |  | 
| 440 | painCave.isFatal = 0; | 
| 441 | painCave.severity = OPENMD_INFO; | 
| 442 | simError(); | 
| 443 | return NULL; | 
| 444 | } else { | 
| 445 |  | 
| 446 | inversion = new Inversion(atomA, atomB, atomC, atomD, inversionType); | 
| 447 | return inversion; | 
| 448 | } | 
| 449 | } | 
| 450 |  | 
| 451 |  | 
| 452 | CutoffGroup* MoleculeCreator::createCutoffGroup(Molecule* mol, CutoffGroupStamp* stamp) { | 
| 453 | int nAtoms; | 
| 454 | CutoffGroup* cg; | 
| 455 | Atom* atom; | 
| 456 | cg = new CutoffGroup(); | 
| 457 |  | 
| 458 | nAtoms = stamp->getNMembers(); | 
| 459 | for (int i =0; i < nAtoms; ++i) { | 
| 460 | atom = mol->getAtomAt(stamp->getMemberAt(i)); | 
| 461 | assert(atom); | 
| 462 | cg->addAtom(atom); | 
| 463 | } | 
| 464 |  | 
| 465 | return cg; | 
| 466 | } | 
| 467 |  | 
| 468 | CutoffGroup* MoleculeCreator::createCutoffGroup(Molecule * mol, Atom* atom) { | 
| 469 | CutoffGroup* cg; | 
| 470 | cg  = new CutoffGroup(); | 
| 471 | cg->addAtom(atom); | 
| 472 | return cg; | 
| 473 | } | 
| 474 |  | 
| 475 | void MoleculeCreator::createConstraintPair(Molecule* mol) { | 
| 476 |  | 
| 477 | //add bond constraints | 
| 478 | Molecule::BondIterator bi; | 
| 479 | Bond* bond; | 
| 480 | for (bond = mol->beginBond(bi); bond != NULL; bond = mol->nextBond(bi)) { | 
| 481 |  | 
| 482 | BondType* bt = bond->getBondType(); | 
| 483 |  | 
| 484 | //class Parent1 {}; | 
| 485 | //class Child1 : public Parent {}; | 
| 486 | //class Child2 : public Parent {}; | 
| 487 | //Child1* ch1 = new Child1(); | 
| 488 | //Child2* ch2 = dynamic_cast<Child2*>(ch1); | 
| 489 | //the dynamic_cast is succeed in above line. A compiler bug? | 
| 490 |  | 
| 491 | if (typeid(FixedBondType) == typeid(*bt)) { | 
| 492 | FixedBondType* fbt = dynamic_cast<FixedBondType*>(bt); | 
| 493 |  | 
| 494 | ConstraintElem* consElemA = new ConstraintElem(bond->getAtomA()); | 
| 495 | ConstraintElem* consElemB = new ConstraintElem(bond->getAtomB()); | 
| 496 | ConstraintPair* consPair = new ConstraintPair(consElemA, consElemB, fbt->getEquilibriumBondLength()); | 
| 497 | mol->addConstraintPair(consPair); | 
| 498 | } | 
| 499 | } | 
| 500 |  | 
| 501 | //rigidbody -- rigidbody constraint is not support yet | 
| 502 | } | 
| 503 |  | 
| 504 | void MoleculeCreator::createConstraintElem(Molecule* mol) { | 
| 505 |  | 
| 506 | ConstraintPair* consPair; | 
| 507 | Molecule::ConstraintPairIterator cpi; | 
| 508 | std::set<StuntDouble*> sdSet; | 
| 509 | for (consPair = mol->beginConstraintPair(cpi); consPair != NULL; consPair = mol->nextConstraintPair(cpi)) { | 
| 510 |  | 
| 511 | StuntDouble* sdA = consPair->getConsElem1()->getStuntDouble(); | 
| 512 | if (sdSet.find(sdA) == sdSet.end()){ | 
| 513 | sdSet.insert(sdA); | 
| 514 | mol->addConstraintElem(new ConstraintElem(sdA)); | 
| 515 | } | 
| 516 |  | 
| 517 | StuntDouble* sdB = consPair->getConsElem2()->getStuntDouble(); | 
| 518 | if (sdSet.find(sdB) == sdSet.end()){ | 
| 519 | sdSet.insert(sdB); | 
| 520 | mol->addConstraintElem(new ConstraintElem(sdB)); | 
| 521 | } | 
| 522 |  | 
| 523 | } | 
| 524 |  | 
| 525 | } | 
| 526 |  | 
| 527 | } |