| 175 |
|
cutoffGroup = createCutoffGroup(mol, *fai, localIndexMan); |
| 176 |
|
mol->addCutoffGroup(cutoffGroup); |
| 177 |
|
} |
| 178 |
< |
//create constraints |
| 178 |
> |
|
| 179 |
> |
//create bonded constraintPairs: |
| 180 |
|
createConstraintPair(mol); |
| 181 |
+ |
|
| 182 |
+ |
//create non-bonded constraintPairs |
| 183 |
+ |
for (int i = 0; i < molStamp->getNConstraints(); ++i) { |
| 184 |
+ |
ConstraintStamp* cStamp = molStamp->getConstraintStamp(i); |
| 185 |
+ |
Atom* atomA; |
| 186 |
+ |
Atom* atomB; |
| 187 |
+ |
|
| 188 |
+ |
atomA = mol->getAtomAt(cStamp->getA()); |
| 189 |
+ |
atomB = mol->getAtomAt(cStamp->getB()); |
| 190 |
+ |
assert( atomA && atomB ); |
| 191 |
+ |
|
| 192 |
+ |
RealType distance; |
| 193 |
+ |
bool printConstraintForce = false; |
| 194 |
+ |
|
| 195 |
+ |
if (!cStamp->haveConstrainedDistance()) { |
| 196 |
+ |
sprintf(painCave.errMsg, |
| 197 |
+ |
"Constraint Error: A non-bond constraint was specified\n" |
| 198 |
+ |
"\twithout providing a value for the constrainedDistance.\n"); |
| 199 |
+ |
painCave.isFatal = 1; |
| 200 |
+ |
simError(); |
| 201 |
+ |
} else { |
| 202 |
+ |
distance = cStamp->getConstrainedDistance(); |
| 203 |
+ |
} |
| 204 |
+ |
|
| 205 |
+ |
if (cStamp->havePrintConstraintForce()) { |
| 206 |
+ |
printConstraintForce = cStamp->getPrintConstraintForce(); |
| 207 |
+ |
} |
| 208 |
+ |
|
| 209 |
+ |
ConstraintElem* consElemA = new ConstraintElem(atomA); |
| 210 |
+ |
ConstraintElem* consElemB = new ConstraintElem(atomB); |
| 211 |
+ |
ConstraintPair* cPair = new ConstraintPair(consElemA, consElemB, distance, |
| 212 |
+ |
printConstraintForce); |
| 213 |
+ |
mol->addConstraintPair(cPair); |
| 214 |
+ |
} |
| 215 |
+ |
|
| 216 |
+ |
// now create the constraint elements: |
| 217 |
+ |
|
| 218 |
|
createConstraintElem(mol); |
| 219 |
|
|
| 220 |
|
// Does this molecule stamp define a total constrained charge value? |
| 582 |
|
//add bond constraints |
| 583 |
|
Molecule::BondIterator bi; |
| 584 |
|
Bond* bond; |
| 585 |
+ |
ConstraintPair* cPair; |
| 586 |
+ |
|
| 587 |
|
for (bond = mol->beginBond(bi); bond != NULL; bond = mol->nextBond(bi)) { |
| 588 |
|
|
| 589 |
|
BondType* bt = bond->getBondType(); |
| 593 |
|
|
| 594 |
|
ConstraintElem* consElemA = new ConstraintElem(bond->getAtomA()); |
| 595 |
|
ConstraintElem* consElemB = new ConstraintElem(bond->getAtomB()); |
| 596 |
< |
ConstraintPair* consPair = new ConstraintPair(consElemA, consElemB, fbt->getEquilibriumBondLength()); |
| 597 |
< |
mol->addConstraintPair(consPair); |
| 596 |
> |
cPair = new ConstraintPair(consElemA, consElemB, |
| 597 |
> |
fbt->getEquilibriumBondLength(), false); |
| 598 |
> |
mol->addConstraintPair(cPair); |
| 599 |
|
} |
| 600 |
|
} |
| 601 |
|
|
| 607 |
|
ConstraintPair* consPair; |
| 608 |
|
Molecule::ConstraintPairIterator cpi; |
| 609 |
|
std::set<StuntDouble*> sdSet; |
| 610 |
< |
for (consPair = mol->beginConstraintPair(cpi); consPair != NULL; consPair = mol->nextConstraintPair(cpi)) { |
| 610 |
> |
for (consPair = mol->beginConstraintPair(cpi); consPair != NULL; |
| 611 |
> |
consPair = mol->nextConstraintPair(cpi)) { |
| 612 |
|
|
| 613 |
|
StuntDouble* sdA = consPair->getConsElem1()->getStuntDouble(); |
| 614 |
|
if (sdSet.find(sdA) == sdSet.end()){ |
| 615 |
|
sdSet.insert(sdA); |
| 616 |
|
mol->addConstraintElem(new ConstraintElem(sdA)); |
| 617 |
|
} |
| 618 |
< |
|
| 618 |
> |
|
| 619 |
|
StuntDouble* sdB = consPair->getConsElem2()->getStuntDouble(); |
| 620 |
|
if (sdSet.find(sdB) == sdSet.end()){ |
| 621 |
|
sdSet.insert(sdB); |
| 622 |
|
mol->addConstraintElem(new ConstraintElem(sdB)); |
| 623 |
< |
} |
| 582 |
< |
|
| 623 |
> |
} |
| 624 |
|
} |
| 584 |
– |
|
| 625 |
|
} |
| 586 |
– |
|
| 626 |
|
} |
