| 281 |
|
// do some sanity checking |
| 282 |
|
|
| 283 |
|
int selectionCount = seleMan_.getSelectionCount(); |
| 284 |
+ |
|
| 285 |
|
int nIntegrable = info->getNGlobalIntegrableObjects(); |
| 286 |
|
|
| 287 |
|
if (selectionCount > nIntegrable) { |
| 909 |
|
|
| 910 |
|
if ((c > 0.81) && (c < 1.21)) {//restrict scaling coefficients |
| 911 |
|
c = sqrt(c); |
| 912 |
< |
//std::cerr << "cold slab scaling coefficient: " << c << endl; |
| 912 |
< |
//now convert to hotBin coefficient |
| 912 |
> |
|
| 913 |
|
RealType w = 0.0; |
| 914 |
|
if (rnemdFluxType_ == rnemdFullKE) { |
| 915 |
|
x = 1.0 + px * (1.0 - c); |
| 947 |
|
} |
| 948 |
|
} |
| 949 |
|
w = sqrt(w); |
| 950 |
– |
// std::cerr << "xh= " << x << "\tyh= " << y << "\tzh= " << z |
| 951 |
– |
// << "\twh= " << w << endl; |
| 950 |
|
for (sdi = hotBin.begin(); sdi != hotBin.end(); sdi++) { |
| 951 |
|
if (rnemdFluxType_ == rnemdFullKE) { |
| 952 |
|
vel = (*sdi)->getVel(); |
| 1252 |
|
|
| 1253 |
|
if (inA) { |
| 1254 |
|
hotBin.push_back(sd); |
| 1257 |
– |
//std::cerr << "before, velocity = " << vel << endl; |
| 1255 |
|
Ph += mass * vel; |
| 1259 |
– |
//std::cerr << "after, velocity = " << vel << endl; |
| 1256 |
|
Mh += mass; |
| 1257 |
|
Kh += mass * vel.lengthSquare(); |
| 1258 |
|
if (rnemdFluxType_ == rnemdFullKE) { |
| 1300 |
|
|
| 1301 |
|
Kh *= 0.5; |
| 1302 |
|
Kc *= 0.5; |
| 1307 |
– |
|
| 1308 |
– |
// std::cerr << "Mh= " << Mh << "\tKh= " << Kh << "\tMc= " << Mc |
| 1309 |
– |
// << "\tKc= " << Kc << endl; |
| 1310 |
– |
// std::cerr << "Ph= " << Ph << "\tPc= " << Pc << endl; |
| 1303 |
|
|
| 1304 |
|
#ifdef IS_MPI |
| 1305 |
|
MPI::COMM_WORLD.Allreduce(MPI::IN_PLACE, &Ph[0], 3, MPI::REALTYPE, MPI::SUM); |
| 1331 |
|
if (hDenominator > 0.0) { |
| 1332 |
|
RealType h = sqrt(hNumerator / hDenominator); |
| 1333 |
|
if ((h > 0.9) && (h < 1.1)) { |
| 1334 |
< |
// std::cerr << "cold slab scaling coefficient: " << c << "\n"; |
| 1343 |
< |
// std::cerr << "hot slab scaling coefficient: " << h << "\n"; |
| 1334 |
> |
|
| 1335 |
|
vector<StuntDouble*>::iterator sdi; |
| 1336 |
|
Vector3d vel; |
| 1337 |
|
for (sdi = coldBin.begin(); sdi != coldBin.end(); sdi++) { |
| 1406 |
|
|
| 1407 |
|
seleMan_.setSelectionSet(evaluator_.evaluate()); |
| 1408 |
|
|
| 1409 |
< |
int selei; |
| 1409 |
> |
int selei(0); |
| 1410 |
|
StuntDouble* sd; |
| 1411 |
|
|
| 1412 |
|
vector<RealType> binMass(nBins_, 0.0); |
| 1431 |
|
sd != NULL; |
| 1432 |
|
sd = mol->nextIntegrableObject(iiter)) |
| 1433 |
|
*/ |
| 1434 |
+ |
|
| 1435 |
|
for (sd = seleMan_.beginSelected(selei); sd != NULL; |
| 1436 |
|
sd = seleMan_.nextSelected(selei)) { |
| 1437 |
< |
|
| 1437 |
> |
|
| 1438 |
|
Vector3d pos = sd->getPos(); |
| 1439 |
|
|
| 1440 |
|
// wrap the stuntdouble's position back into the box: |
| 1449 |
|
// The modulo operator is used to wrap the case when we are |
| 1450 |
|
// beyond the end of the bins back to the beginning. |
| 1451 |
|
int binNo = int(nBins_ * (pos.z() / hmat(2,2) + 0.5)) % nBins_; |
| 1452 |
< |
|
| 1452 |
> |
|
| 1453 |
|
RealType mass = sd->getMass(); |
| 1454 |
|
Vector3d vel = sd->getVel(); |
| 1455 |
|
|
| 1480 |
|
} |
| 1481 |
|
} |
| 1482 |
|
|
| 1491 |
– |
|
| 1483 |
|
#ifdef IS_MPI |
| 1484 |
|
MPI::COMM_WORLD.Allreduce(MPI::IN_PLACE, &binCount[0], |
| 1485 |
|
nBins_, MPI::INT, MPI::SUM); |
| 1510 |
|
den = binMass[i] * nBins_ * PhysicalConstants::densityConvert |
| 1511 |
|
/ currentSnap_->getVolume() ; |
| 1512 |
|
|
| 1513 |
< |
temp = 2.0 * binKE[i] / (binDOF[i] * PhysicalConstants::kb * |
| 1514 |
< |
PhysicalConstants::energyConvert); |
| 1515 |
< |
|
| 1516 |
< |
for (unsigned int j = 0; j < outputMask_.size(); ++j) { |
| 1517 |
< |
if(outputMask_[j]) { |
| 1518 |
< |
switch(j) { |
| 1519 |
< |
case Z: |
| 1520 |
< |
dynamic_cast<Accumulator *>(data_[j].accumulator[i])->add(z); |
| 1521 |
< |
break; |
| 1522 |
< |
case TEMPERATURE: |
| 1523 |
< |
dynamic_cast<Accumulator *>(data_[j].accumulator[i])->add(temp); |
| 1524 |
< |
break; |
| 1525 |
< |
case VELOCITY: |
| 1526 |
< |
dynamic_cast<VectorAccumulator *>(data_[j].accumulator[i])->add(vel); |
| 1527 |
< |
break; |
| 1528 |
< |
case DENSITY: |
| 1529 |
< |
dynamic_cast<Accumulator *>(data_[j].accumulator[i])->add(den); |
| 1530 |
< |
break; |
| 1513 |
> |
if (binCount[i] > 0) { |
| 1514 |
> |
// only add values if there are things to add |
| 1515 |
> |
temp = 2.0 * binKE[i] / (binDOF[i] * PhysicalConstants::kb * |
| 1516 |
> |
PhysicalConstants::energyConvert); |
| 1517 |
> |
|
| 1518 |
> |
for (unsigned int j = 0; j < outputMask_.size(); ++j) { |
| 1519 |
> |
if(outputMask_[j]) { |
| 1520 |
> |
switch(j) { |
| 1521 |
> |
case Z: |
| 1522 |
> |
dynamic_cast<Accumulator *>(data_[j].accumulator[i])->add(z); |
| 1523 |
> |
break; |
| 1524 |
> |
case TEMPERATURE: |
| 1525 |
> |
dynamic_cast<Accumulator *>(data_[j].accumulator[i])->add(temp); |
| 1526 |
> |
break; |
| 1527 |
> |
case VELOCITY: |
| 1528 |
> |
dynamic_cast<VectorAccumulator *>(data_[j].accumulator[i])->add(vel); |
| 1529 |
> |
break; |
| 1530 |
> |
case DENSITY: |
| 1531 |
> |
dynamic_cast<Accumulator *>(data_[j].accumulator[i])->add(den); |
| 1532 |
> |
break; |
| 1533 |
> |
} |
| 1534 |
|
} |
| 1535 |
|
} |
| 1536 |
|
} |