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). |
39 |
> |
* [4] Kuang & Gezelter, J. Chem. Phys. 133, 164101 (2010). |
40 |
> |
* [5] Vardeman, Stocker & Gezelter, J. Chem. Theory Comput. 7, 834 (2011). |
41 |
|
*/ |
42 |
|
|
43 |
|
#define _LARGEFILE_SOURCE64 |
58 |
|
#include "utils/simError.h" |
59 |
|
#include "utils/MemoryUtils.hpp" |
60 |
|
#include "utils/StringTokenizer.hpp" |
61 |
+ |
#include "brains/Thermo.hpp" |
62 |
|
|
63 |
|
#ifdef IS_MPI |
62 |
– |
|
64 |
|
#include <mpi.h> |
65 |
< |
#define TAKE_THIS_TAG_CHAR 0 |
65 |
< |
#define TAKE_THIS_TAG_INT 1 |
65 |
> |
#endif |
66 |
|
|
67 |
– |
#endif // is_mpi |
67 |
|
|
69 |
– |
|
68 |
|
namespace OpenMD { |
69 |
|
|
70 |
|
DumpReader::DumpReader(SimInfo* info, const std::string& filename) |
140 |
|
prevPos = currPos; |
141 |
|
bool foundOpenSnapshotTag = false; |
142 |
|
bool foundClosedSnapshotTag = false; |
143 |
+ |
bool foundOpenSiteDataTag = false; |
144 |
|
while(inFile_->getline(buffer, bufferSize)) { |
145 |
|
++lineNo; |
146 |
|
|
242 |
|
|
243 |
|
if (needCOMprops_) { |
244 |
|
Snapshot* s = info_->getSnapshotManager()->getCurrentSnapshot(); |
245 |
+ |
Thermo thermo(info_); |
246 |
|
Vector3d com; |
247 |
– |
Vector3d comvel; |
248 |
– |
Vector3d comw; |
249 |
– |
if (needPos_ && needVel_){ |
250 |
– |
info_->getComAll(com, comvel); |
251 |
– |
comw = info_->getAngularMomentum(); |
252 |
– |
}else{ |
253 |
– |
com = info_->getCom(); |
254 |
– |
comvel = 0.0; |
255 |
– |
comw = 0.0; |
256 |
– |
} |
257 |
– |
s->setCOMprops(com, comvel, comw); |
258 |
– |
} |
247 |
|
|
248 |
+ |
if (needPos_ && needVel_) { |
249 |
+ |
Vector3d comvel; |
250 |
+ |
Vector3d comw; |
251 |
+ |
thermo.getComAll(com, comvel); |
252 |
+ |
comw = thermo.getAngularMomentum(); |
253 |
+ |
} else { |
254 |
+ |
com = thermo.getCom(); |
255 |
+ |
} |
256 |
+ |
} |
257 |
|
} |
258 |
|
|
259 |
|
void DumpReader::readSet(int whichFrame) { |
321 |
|
|
322 |
|
inputStream.getline(buffer, bufferSize); |
323 |
|
line = buffer; |
324 |
< |
if (line.find("</Snapshot>") == std::string::npos) { |
325 |
< |
sprintf(painCave.errMsg, |
326 |
< |
"DumpReader Error: can not find </Snapshot>\n"); |
327 |
< |
painCave.isFatal = 1; |
328 |
< |
simError(); |
329 |
< |
} |
330 |
< |
|
324 |
> |
|
325 |
> |
if (line.find("<SiteData>") != std::string::npos) { |
326 |
> |
//read SiteData |
327 |
> |
readSiteData(inputStream); |
328 |
> |
} else { |
329 |
> |
if (line.find("</Snapshot>") == std::string::npos) { |
330 |
> |
sprintf(painCave.errMsg, |
331 |
> |
"DumpReader Error: can not find </Snapshot>\n"); |
332 |
> |
painCave.isFatal = 1; |
333 |
> |
simError(); |
334 |
> |
} |
335 |
> |
} |
336 |
|
} |
337 |
|
|
338 |
|
void DumpReader::parseDumpLine(const std::string& line) { |
468 |
|
integrableObject->setTrq(torque); |
469 |
|
break; |
470 |
|
} |
471 |
+ |
case 'u' : { |
472 |
+ |
|
473 |
+ |
RealType particlePot; |
474 |
+ |
particlePot = tokenizer.nextTokenAsDouble(); |
475 |
+ |
integrableObject->setParticlePot(particlePot); |
476 |
+ |
break; |
477 |
+ |
} |
478 |
+ |
case 'c' : { |
479 |
+ |
|
480 |
+ |
RealType flucQPos; |
481 |
+ |
flucQPos = tokenizer.nextTokenAsDouble(); |
482 |
+ |
integrableObject->setFlucQPos(flucQPos); |
483 |
+ |
break; |
484 |
+ |
} |
485 |
+ |
case 'w' : { |
486 |
+ |
|
487 |
+ |
RealType flucQVel; |
488 |
+ |
flucQVel = tokenizer.nextTokenAsDouble(); |
489 |
+ |
integrableObject->setFlucQVel(flucQVel); |
490 |
+ |
break; |
491 |
+ |
} |
492 |
+ |
case 'g' : { |
493 |
+ |
|
494 |
+ |
RealType flucQFrc; |
495 |
+ |
flucQFrc = tokenizer.nextTokenAsDouble(); |
496 |
+ |
integrableObject->setFlucQFrc(flucQFrc); |
497 |
+ |
break; |
498 |
+ |
} |
499 |
+ |
case 'e' : { |
500 |
+ |
|
501 |
+ |
Vector3d eField; |
502 |
+ |
eField[0] = tokenizer.nextTokenAsDouble(); |
503 |
+ |
eField[1] = tokenizer.nextTokenAsDouble(); |
504 |
+ |
eField[2] = tokenizer.nextTokenAsDouble(); |
505 |
+ |
integrableObject->setElectricField(eField); |
506 |
+ |
break; |
507 |
+ |
} |
508 |
|
default: { |
509 |
|
sprintf(painCave.errMsg, |
510 |
|
"DumpReader Error: %s is an unrecognized type\n", type.c_str()); |
519 |
|
} |
520 |
|
|
521 |
|
|
522 |
< |
void DumpReader::readStuntDoubles(std::istream& inputStream) { |
522 |
> |
void DumpReader::parseSiteLine(const std::string& line) { |
523 |
> |
|
524 |
> |
StringTokenizer tokenizer(line); |
525 |
> |
int nTokens; |
526 |
> |
|
527 |
> |
nTokens = tokenizer.countTokens(); |
528 |
> |
|
529 |
> |
if (nTokens < 2) { |
530 |
> |
sprintf(painCave.errMsg, |
531 |
> |
"DumpReader Error: Not enough Tokens.\n%s\n", line.c_str()); |
532 |
> |
painCave.isFatal = 1; |
533 |
> |
simError(); |
534 |
> |
} |
535 |
> |
|
536 |
> |
/** |
537 |
> |
* The first token is the global integrable object index. |
538 |
> |
*/ |
539 |
> |
|
540 |
> |
int index = tokenizer.nextTokenAsInt(); |
541 |
> |
StuntDouble* integrableObject = info_->getIOIndexToIntegrableObject(index); |
542 |
> |
if (integrableObject == NULL) { |
543 |
> |
return; |
544 |
> |
} |
545 |
> |
StuntDouble* sd = integrableObject; |
546 |
|
|
547 |
+ |
/** |
548 |
+ |
* Test to see if the next token is an integer or not. If not, |
549 |
+ |
* we've got data on the integrable object itself. If there is an |
550 |
+ |
* integer, we're parsing data for a site on a rigid body. |
551 |
+ |
*/ |
552 |
+ |
|
553 |
+ |
std::string indexTest = tokenizer.peekNextToken(); |
554 |
+ |
std::istringstream i(indexTest); |
555 |
+ |
int siteIndex; |
556 |
+ |
if (i >> siteIndex) { |
557 |
+ |
// chew up this token and parse as an int: |
558 |
+ |
siteIndex = tokenizer.nextTokenAsInt(); |
559 |
+ |
RigidBody* rb = static_cast<RigidBody*>(integrableObject); |
560 |
+ |
sd = rb->getAtoms()[siteIndex]; |
561 |
+ |
} |
562 |
+ |
|
563 |
+ |
/** |
564 |
+ |
* The next token contains information on what follows. |
565 |
+ |
*/ |
566 |
+ |
std::string type = tokenizer.nextToken(); |
567 |
+ |
int size = type.size(); |
568 |
+ |
|
569 |
+ |
for(int i = 0; i < size; ++i) { |
570 |
+ |
switch(type[i]) { |
571 |
+ |
|
572 |
+ |
case 'u' : { |
573 |
+ |
|
574 |
+ |
RealType particlePot; |
575 |
+ |
particlePot = tokenizer.nextTokenAsDouble(); |
576 |
+ |
sd->setParticlePot(particlePot); |
577 |
+ |
break; |
578 |
+ |
} |
579 |
+ |
case 'c' : { |
580 |
+ |
|
581 |
+ |
RealType flucQPos; |
582 |
+ |
flucQPos = tokenizer.nextTokenAsDouble(); |
583 |
+ |
sd->setFlucQPos(flucQPos); |
584 |
+ |
break; |
585 |
+ |
} |
586 |
+ |
case 'w' : { |
587 |
+ |
|
588 |
+ |
RealType flucQVel; |
589 |
+ |
flucQVel = tokenizer.nextTokenAsDouble(); |
590 |
+ |
sd->setFlucQVel(flucQVel); |
591 |
+ |
break; |
592 |
+ |
} |
593 |
+ |
case 'g' : { |
594 |
+ |
|
595 |
+ |
RealType flucQFrc; |
596 |
+ |
flucQFrc = tokenizer.nextTokenAsDouble(); |
597 |
+ |
sd->setFlucQFrc(flucQFrc); |
598 |
+ |
break; |
599 |
+ |
} |
600 |
+ |
case 'e' : { |
601 |
+ |
|
602 |
+ |
Vector3d eField; |
603 |
+ |
eField[0] = tokenizer.nextTokenAsDouble(); |
604 |
+ |
eField[1] = tokenizer.nextTokenAsDouble(); |
605 |
+ |
eField[2] = tokenizer.nextTokenAsDouble(); |
606 |
+ |
sd->setElectricField(eField); |
607 |
+ |
break; |
608 |
+ |
} |
609 |
+ |
default: { |
610 |
+ |
sprintf(painCave.errMsg, |
611 |
+ |
"DumpReader Error: %s is an unrecognized type\n", type.c_str()); |
612 |
+ |
painCave.isFatal = 1; |
613 |
+ |
simError(); |
614 |
+ |
break; |
615 |
+ |
} |
616 |
+ |
} |
617 |
+ |
} |
618 |
+ |
} |
619 |
+ |
|
620 |
+ |
|
621 |
+ |
void DumpReader::readStuntDoubles(std::istream& inputStream) { |
622 |
+ |
|
623 |
|
inputStream.getline(buffer, bufferSize); |
624 |
|
std::string line(buffer); |
625 |
|
|
642 |
|
|
643 |
|
} |
644 |
|
|
645 |
+ |
void DumpReader::readSiteData(std::istream& inputStream) { |
646 |
+ |
|
647 |
+ |
inputStream.getline(buffer, bufferSize); |
648 |
+ |
std::string line(buffer); |
649 |
+ |
|
650 |
+ |
if (line.find("<SiteData>") == std::string::npos) { |
651 |
+ |
// site data isn't required for a simulation, so skip |
652 |
+ |
return; |
653 |
+ |
} |
654 |
+ |
|
655 |
+ |
while(inputStream.getline(buffer, bufferSize)) { |
656 |
+ |
line = buffer; |
657 |
+ |
|
658 |
+ |
if(line.find("</SiteData>") != std::string::npos) { |
659 |
+ |
break; |
660 |
+ |
} |
661 |
+ |
|
662 |
+ |
parseSiteLine(line); |
663 |
+ |
} |
664 |
+ |
|
665 |
+ |
} |
666 |
+ |
|
667 |
|
void DumpReader::readFrameProperties(std::istream& inputStream) { |
668 |
|
|
669 |
|
Snapshot* s = info_->getSnapshotManager()->getCurrentSnapshot(); |
709 |
|
hmat(2, 2) = tokenizer.nextTokenAsDouble(); |
710 |
|
s->setHmat(hmat); |
711 |
|
} else if (propertyName == "Thermostat") { |
712 |
< |
RealType chi = tokenizer.nextTokenAsDouble(); |
713 |
< |
RealType integralOfChiDt = tokenizer.nextTokenAsDouble(); |
714 |
< |
s->setChi(chi); |
715 |
< |
s->setIntegralOfChiDt(integralOfChiDt); |
712 |
> |
pair<RealType, RealType> thermostat; |
713 |
> |
thermostat.first = tokenizer.nextTokenAsDouble(); |
714 |
> |
thermostat.second = tokenizer.nextTokenAsDouble(); |
715 |
> |
s->setThermostat(thermostat); |
716 |
|
} else if (propertyName == "Barostat") { |
717 |
|
Mat3x3d eta; |
718 |
|
eta(0, 0) = tokenizer.nextTokenAsDouble(); |
724 |
|
eta(2, 0) = tokenizer.nextTokenAsDouble(); |
725 |
|
eta(2, 1) = tokenizer.nextTokenAsDouble(); |
726 |
|
eta(2, 2) = tokenizer.nextTokenAsDouble(); |
727 |
< |
s->setEta(eta); |
727 |
> |
s->setBarostat(eta); |
728 |
|
} else { |
729 |
|
sprintf(painCave.errMsg, |
730 |
|
"DumpReader Error: %s is an invalid property in <FrameData>\n", propertyName.c_str()); |