| 48 | 
  | 
 | 
| 49 | 
  | 
#include <algorithm> | 
| 50 | 
  | 
#include <set> | 
| 51 | 
+ | 
#include <map> | 
| 52 | 
  | 
 | 
| 53 | 
  | 
#include "brains/SimInfo.hpp" | 
| 54 | 
  | 
#include "math/Vector3.hpp" | 
| 71 | 
  | 
#endif  | 
| 72 | 
  | 
 | 
| 73 | 
  | 
namespace oopse { | 
| 74 | 
+ | 
  std::set<int> getRigidSet(int index, std::map<int, std::set<int> >& container) { | 
| 75 | 
+ | 
    std::map<int, std::set<int> >::iterator i = container.find(index); | 
| 76 | 
+ | 
    std::set<int> result; | 
| 77 | 
+ | 
    if (i != container.end()) { | 
| 78 | 
+ | 
        result = i->second; | 
| 79 | 
+ | 
    } | 
| 80 | 
  | 
 | 
| 81 | 
+ | 
    return result; | 
| 82 | 
+ | 
  } | 
| 83 | 
+ | 
   | 
| 84 | 
  | 
  SimInfo::SimInfo(MakeStamps* stamps, std::vector<std::pair<MoleculeStamp*, int> >& molStampPairs,  | 
| 85 | 
  | 
                   ForceField* ff, Globals* simParams) :  | 
| 86 | 
  | 
    stamps_(stamps), forceField_(ff), simParams_(simParams),  | 
| 355 | 
  | 
    int b; | 
| 356 | 
  | 
    int c; | 
| 357 | 
  | 
    int d; | 
| 358 | 
+ | 
 | 
| 359 | 
+ | 
    std::map<int, std::set<int> > atomGroups; | 
| 360 | 
+ | 
 | 
| 361 | 
+ | 
    Molecule::RigidBodyIterator rbIter; | 
| 362 | 
+ | 
    RigidBody* rb; | 
| 363 | 
+ | 
    Molecule::IntegrableObjectIterator ii; | 
| 364 | 
+ | 
    StuntDouble* integrableObject; | 
| 365 | 
+ | 
     | 
| 366 | 
+ | 
    for (integrableObject = mol->beginIntegrableObject(ii); integrableObject != NULL; | 
| 367 | 
+ | 
           integrableObject = mol->nextIntegrableObject(ii)) { | 
| 368 | 
+ | 
 | 
| 369 | 
+ | 
      if (integrableObject->isRigidBody()) { | 
| 370 | 
+ | 
          rb = static_cast<RigidBody*>(integrableObject); | 
| 371 | 
+ | 
          std::vector<Atom*> atoms = rb->getAtoms(); | 
| 372 | 
+ | 
          std::set<int> rigidAtoms; | 
| 373 | 
+ | 
          for (int i = 0; i < atoms.size(); ++i) { | 
| 374 | 
+ | 
            rigidAtoms.insert(atoms[i]->getGlobalIndex()); | 
| 375 | 
+ | 
          } | 
| 376 | 
+ | 
          for (int i = 0; i < atoms.size(); ++i) { | 
| 377 | 
+ | 
            atomGroups.insert(std::map<int, std::set<int> >::value_type(atoms[i]->getGlobalIndex(), rigidAtoms)); | 
| 378 | 
+ | 
          }       | 
| 379 | 
+ | 
      } else { | 
| 380 | 
+ | 
        std::set<int> oneAtomSet; | 
| 381 | 
+ | 
        oneAtomSet.insert(integrableObject->getGlobalIndex()); | 
| 382 | 
+ | 
        atomGroups.insert(std::map<int, std::set<int> >::value_type(integrableObject->getGlobalIndex(), oneAtomSet));         | 
| 383 | 
+ | 
      } | 
| 384 | 
+ | 
    }   | 
| 385 | 
+ | 
 | 
| 386 | 
  | 
     | 
| 387 | 
+ | 
     | 
| 388 | 
  | 
    for (bond= mol->beginBond(bondIter); bond != NULL; bond = mol->nextBond(bondIter)) { | 
| 389 | 
  | 
      a = bond->getAtomA()->getGlobalIndex(); | 
| 390 | 
  | 
      b = bond->getAtomB()->getGlobalIndex();         | 
| 395 | 
  | 
      a = bend->getAtomA()->getGlobalIndex(); | 
| 396 | 
  | 
      b = bend->getAtomB()->getGlobalIndex();         | 
| 397 | 
  | 
      c = bend->getAtomC()->getGlobalIndex(); | 
| 398 | 
+ | 
      std::set<int> rigidSetA = getRigidSet(a, atomGroups); | 
| 399 | 
+ | 
      std::set<int> rigidSetB = getRigidSet(b, atomGroups); | 
| 400 | 
+ | 
      std::set<int> rigidSetC = getRigidSet(c, atomGroups); | 
| 401 | 
  | 
 | 
| 402 | 
< | 
      exclude_.addPair(a, b); | 
| 403 | 
< | 
      exclude_.addPair(a, c); | 
| 404 | 
< | 
      exclude_.addPair(b, c);         | 
| 402 | 
> | 
      exclude_.addPairs(rigidSetA, rigidSetB); | 
| 403 | 
> | 
      exclude_.addPairs(rigidSetA, rigidSetC); | 
| 404 | 
> | 
      exclude_.addPairs(rigidSetB, rigidSetC); | 
| 405 | 
> | 
       | 
| 406 | 
> | 
      //exclude_.addPair(a, b); | 
| 407 | 
> | 
      //exclude_.addPair(a, c); | 
| 408 | 
> | 
      //exclude_.addPair(b, c);         | 
| 409 | 
  | 
    } | 
| 410 | 
  | 
 | 
| 411 | 
  | 
    for (torsion= mol->beginTorsion(torsionIter); torsion != NULL; torsion = mol->nextTorsion(torsionIter)) { | 
| 413 | 
  | 
      b = torsion->getAtomB()->getGlobalIndex();         | 
| 414 | 
  | 
      c = torsion->getAtomC()->getGlobalIndex();         | 
| 415 | 
  | 
      d = torsion->getAtomD()->getGlobalIndex();         | 
| 416 | 
+ | 
      std::set<int> rigidSetA = getRigidSet(a, atomGroups); | 
| 417 | 
+ | 
      std::set<int> rigidSetB = getRigidSet(b, atomGroups); | 
| 418 | 
+ | 
      std::set<int> rigidSetC = getRigidSet(c, atomGroups); | 
| 419 | 
+ | 
      std::set<int> rigidSetD = getRigidSet(d, atomGroups); | 
| 420 | 
  | 
 | 
| 421 | 
+ | 
      exclude_.addPairs(rigidSetA, rigidSetB); | 
| 422 | 
+ | 
      exclude_.addPairs(rigidSetA, rigidSetC); | 
| 423 | 
+ | 
      exclude_.addPairs(rigidSetA, rigidSetD); | 
| 424 | 
+ | 
      exclude_.addPairs(rigidSetB, rigidSetC); | 
| 425 | 
+ | 
      exclude_.addPairs(rigidSetB, rigidSetD); | 
| 426 | 
+ | 
      exclude_.addPairs(rigidSetC, rigidSetD); | 
| 427 | 
+ | 
 | 
| 428 | 
+ | 
      /* | 
| 429 | 
+ | 
      exclude_.addPairs(rigidSetA.begin(), rigidSetA.end(), rigidSetB.begin(), rigidSetB.end()); | 
| 430 | 
+ | 
      exclude_.addPairs(rigidSetA.begin(), rigidSetA.end(), rigidSetC.begin(), rigidSetC.end()); | 
| 431 | 
+ | 
      exclude_.addPairs(rigidSetA.begin(), rigidSetA.end(), rigidSetD.begin(), rigidSetD.end()); | 
| 432 | 
+ | 
      exclude_.addPairs(rigidSetB.begin(), rigidSetB.end(), rigidSetC.begin(), rigidSetC.end()); | 
| 433 | 
+ | 
      exclude_.addPairs(rigidSetB.begin(), rigidSetB.end(), rigidSetD.begin(), rigidSetD.end()); | 
| 434 | 
+ | 
      exclude_.addPairs(rigidSetC.begin(), rigidSetC.end(), rigidSetD.begin(), rigidSetD.end()); | 
| 435 | 
+ | 
          | 
| 436 | 
+ | 
       | 
| 437 | 
  | 
      exclude_.addPair(a, b); | 
| 438 | 
  | 
      exclude_.addPair(a, c); | 
| 439 | 
  | 
      exclude_.addPair(a, d); | 
| 440 | 
  | 
      exclude_.addPair(b, c); | 
| 441 | 
  | 
      exclude_.addPair(b, d); | 
| 442 | 
  | 
      exclude_.addPair(c, d);         | 
| 443 | 
+ | 
      */ | 
| 444 | 
  | 
    } | 
| 445 | 
  | 
 | 
| 379 | 
– | 
    Molecule::RigidBodyIterator rbIter; | 
| 380 | 
– | 
    RigidBody* rb; | 
| 446 | 
  | 
    for (rb = mol->beginRigidBody(rbIter); rb != NULL; rb = mol->nextRigidBody(rbIter)) { | 
| 447 | 
  | 
      std::vector<Atom*> atoms = rb->getAtoms(); | 
| 448 | 
  | 
      for (int i = 0; i < atoms.size() -1 ; ++i) { | 
| 467 | 
  | 
    int b; | 
| 468 | 
  | 
    int c; | 
| 469 | 
  | 
    int d; | 
| 470 | 
+ | 
 | 
| 471 | 
+ | 
    std::map<int, std::set<int> > atomGroups; | 
| 472 | 
+ | 
 | 
| 473 | 
+ | 
    Molecule::RigidBodyIterator rbIter; | 
| 474 | 
+ | 
    RigidBody* rb; | 
| 475 | 
+ | 
    Molecule::IntegrableObjectIterator ii; | 
| 476 | 
+ | 
    StuntDouble* integrableObject; | 
| 477 | 
  | 
     | 
| 478 | 
+ | 
    for (integrableObject = mol->beginIntegrableObject(ii); integrableObject != NULL; | 
| 479 | 
+ | 
           integrableObject = mol->nextIntegrableObject(ii)) { | 
| 480 | 
+ | 
 | 
| 481 | 
+ | 
      if (integrableObject->isRigidBody()) { | 
| 482 | 
+ | 
          rb = static_cast<RigidBody*>(integrableObject); | 
| 483 | 
+ | 
          std::vector<Atom*> atoms = rb->getAtoms(); | 
| 484 | 
+ | 
          std::set<int> rigidAtoms; | 
| 485 | 
+ | 
          for (int i = 0; i < atoms.size(); ++i) { | 
| 486 | 
+ | 
            rigidAtoms.insert(atoms[i]->getGlobalIndex()); | 
| 487 | 
+ | 
          } | 
| 488 | 
+ | 
          for (int i = 0; i < atoms.size(); ++i) { | 
| 489 | 
+ | 
            atomGroups.insert(std::map<int, std::set<int> >::value_type(atoms[i]->getGlobalIndex(), rigidAtoms)); | 
| 490 | 
+ | 
          }       | 
| 491 | 
+ | 
      } else { | 
| 492 | 
+ | 
        std::set<int> oneAtomSet; | 
| 493 | 
+ | 
        oneAtomSet.insert(integrableObject->getGlobalIndex()); | 
| 494 | 
+ | 
        atomGroups.insert(std::map<int, std::set<int> >::value_type(integrableObject->getGlobalIndex(), oneAtomSet));         | 
| 495 | 
+ | 
      } | 
| 496 | 
+ | 
    }   | 
| 497 | 
+ | 
 | 
| 498 | 
+ | 
     | 
| 499 | 
  | 
    for (bond= mol->beginBond(bondIter); bond != NULL; bond = mol->nextBond(bondIter)) { | 
| 500 | 
  | 
      a = bond->getAtomA()->getGlobalIndex(); | 
| 501 | 
  | 
      b = bond->getAtomB()->getGlobalIndex();         | 
| 507 | 
  | 
      b = bend->getAtomB()->getGlobalIndex();         | 
| 508 | 
  | 
      c = bend->getAtomC()->getGlobalIndex(); | 
| 509 | 
  | 
 | 
| 510 | 
< | 
      exclude_.removePair(a, b); | 
| 511 | 
< | 
      exclude_.removePair(a, c); | 
| 512 | 
< | 
      exclude_.removePair(b, c);         | 
| 510 | 
> | 
      std::set<int> rigidSetA = getRigidSet(a, atomGroups); | 
| 511 | 
> | 
      std::set<int> rigidSetB = getRigidSet(b, atomGroups); | 
| 512 | 
> | 
      std::set<int> rigidSetC = getRigidSet(c, atomGroups); | 
| 513 | 
> | 
 | 
| 514 | 
> | 
      exclude_.removePairs(rigidSetA, rigidSetB); | 
| 515 | 
> | 
      exclude_.removePairs(rigidSetA, rigidSetC); | 
| 516 | 
> | 
      exclude_.removePairs(rigidSetB, rigidSetC); | 
| 517 | 
> | 
       | 
| 518 | 
> | 
      //exclude_.removePair(a, b); | 
| 519 | 
> | 
      //exclude_.removePair(a, c); | 
| 520 | 
> | 
      //exclude_.removePair(b, c);         | 
| 521 | 
  | 
    } | 
| 522 | 
  | 
 | 
| 523 | 
  | 
    for (torsion= mol->beginTorsion(torsionIter); torsion != NULL; torsion = mol->nextTorsion(torsionIter)) { | 
| 525 | 
  | 
      b = torsion->getAtomB()->getGlobalIndex();         | 
| 526 | 
  | 
      c = torsion->getAtomC()->getGlobalIndex();         | 
| 527 | 
  | 
      d = torsion->getAtomD()->getGlobalIndex();         | 
| 528 | 
+ | 
 | 
| 529 | 
+ | 
      std::set<int> rigidSetA = getRigidSet(a, atomGroups); | 
| 530 | 
+ | 
      std::set<int> rigidSetB = getRigidSet(b, atomGroups); | 
| 531 | 
+ | 
      std::set<int> rigidSetC = getRigidSet(c, atomGroups); | 
| 532 | 
+ | 
      std::set<int> rigidSetD = getRigidSet(d, atomGroups); | 
| 533 | 
+ | 
 | 
| 534 | 
+ | 
      exclude_.removePairs(rigidSetA, rigidSetB); | 
| 535 | 
+ | 
      exclude_.removePairs(rigidSetA, rigidSetC); | 
| 536 | 
+ | 
      exclude_.removePairs(rigidSetA, rigidSetD); | 
| 537 | 
+ | 
      exclude_.removePairs(rigidSetB, rigidSetC); | 
| 538 | 
+ | 
      exclude_.removePairs(rigidSetB, rigidSetD); | 
| 539 | 
+ | 
      exclude_.removePairs(rigidSetC, rigidSetD); | 
| 540 | 
+ | 
 | 
| 541 | 
+ | 
      /* | 
| 542 | 
+ | 
      exclude_.removePairs(rigidSetA.begin(), rigidSetA.end(), rigidSetB.begin(), rigidSetB.end()); | 
| 543 | 
+ | 
      exclude_.removePairs(rigidSetA.begin(), rigidSetA.end(), rigidSetC.begin(), rigidSetC.end()); | 
| 544 | 
+ | 
      exclude_.removePairs(rigidSetA.begin(), rigidSetA.end(), rigidSetD.begin(), rigidSetD.end()); | 
| 545 | 
+ | 
      exclude_.removePairs(rigidSetB.begin(), rigidSetB.end(), rigidSetC.begin(), rigidSetC.end()); | 
| 546 | 
+ | 
      exclude_.removePairs(rigidSetB.begin(), rigidSetB.end(), rigidSetD.begin(), rigidSetD.end()); | 
| 547 | 
+ | 
      exclude_.removePairs(rigidSetC.begin(), rigidSetC.end(), rigidSetD.begin(), rigidSetD.end()); | 
| 548 | 
  | 
 | 
| 549 | 
+ | 
       | 
| 550 | 
  | 
      exclude_.removePair(a, b); | 
| 551 | 
  | 
      exclude_.removePair(a, c); | 
| 552 | 
  | 
      exclude_.removePair(a, d); | 
| 553 | 
  | 
      exclude_.removePair(b, c); | 
| 554 | 
  | 
      exclude_.removePair(b, d); | 
| 555 | 
  | 
      exclude_.removePair(c, d);         | 
| 556 | 
+ | 
      */ | 
| 557 | 
  | 
    } | 
| 558 | 
  | 
 | 
| 436 | 
– | 
    Molecule::RigidBodyIterator rbIter; | 
| 437 | 
– | 
    RigidBody* rb; | 
| 559 | 
  | 
    for (rb = mol->beginRigidBody(rbIter); rb != NULL; rb = mol->nextRigidBody(rbIter)) { | 
| 560 | 
  | 
      std::vector<Atom*> atoms = rb->getAtoms(); | 
| 561 | 
  | 
      for (int i = 0; i < atoms.size() -1 ; ++i) { |