50#include "applications/staticProps/ChargeZ.hpp"
55#include "brains/Thermo.hpp"
58#include "types/FixedChargeAdapter.hpp"
59#include "types/FluctuatingChargeAdapter.hpp"
60#include "utils/simError.h"
64 ChargeZ::ChargeZ(
SimInfo* info,
const std::string& filename,
65 const std::string& sele,
int nzbins,
int axis) :
67 selectionScript_(sele), evaluator_(info), seleMan_(info), thermo_(info),
69 evaluator_.loadScriptString(sele);
70 if (!evaluator_.isDynamic()) {
71 seleMan_.setSelectionSet(evaluator_.evaluate());
76 sliceSDLists_.resize(nBins_);
77 sliceSDCount_.resize(nBins_);
78 std::fill(sliceSDCount_.begin(), sliceSDCount_.end(), 0);
80 chargeZ_.resize(nBins_);
95 setOutputName(
getPrefix(filename) +
".ChargeZ");
98 void ChargeZ::process() {
102 bool usePeriodicBoundaryConditions_ =
103 info_->getSimParams()->getUsePeriodicBoundaryConditions();
106 int nFrames = reader.getNFrames();
107 nProcessed_ = nFrames / step_;
109 for (
int istep = 0; istep < nFrames; istep += step_) {
110 reader.readFrame(istep);
113 for (
unsigned int i = 0; i < nBins_; i++) {
114 sliceSDLists_[i].clear();
118 zBox_.push_back(hmat(axis_, axis_));
120 RealType halfBoxZ_ = hmat(axis_, axis_) / 2.0;
124 area = currentSnapshot_->getYZarea();
127 area = currentSnapshot_->getXZarea();
131 area = currentSnapshot_->getXYarea();
135 areas_.push_back(area);
138 seleMan_.setSelectionSet(evaluator_.evaluate());
145 if (usePeriodicBoundaryConditions_) currentSnapshot_->
wrapVector(pos);
154 int binNo = int(nBins_ * (halfBoxZ_ + pos[axis_]) / hmat(axis_, axis_));
155 sliceSDLists_[binNo].push_back(sd);
156 sliceSDCount_[binNo]++;
160 for (
unsigned int i = 0; i < nBins_; i++) {
162 for (
unsigned int k = 0; k < sliceSDLists_[i].size(); ++k) {
164 Atom* atom =
static_cast<Atom*
>(sliceSDLists_[i][k]);
166 AtomType* atomType = atom->getAtomType();
168 if (sliceSDLists_[i][k]->isAtom()) {
170 if (fca.isFixedCharge()) { q += fca.getCharge(); }
173 if (fqa.isFluctuatingCharge()) { q += atom->getFlucQPos(); }
186 void ChargeZ::writeChargeZ() {
188 std::vector<RealType>::iterator j;
190 for (j = zBox_.begin(); j != zBox_.end(); ++j) {
193 RealType zAve = zSum / zBox_.size();
195 RealType areaSum = 0.0;
196 for (j = areas_.begin(); j != areas_.end(); ++j) {
199 RealType areaAve = areaSum / areas_.size();
201 std::ofstream rdfStream(outputFilename_.c_str());
202 if (rdfStream.is_open()) {
203 rdfStream <<
"#ChargeZ "
205 rdfStream <<
"#selection: (" << selectionScript_ <<
")\n";
206 rdfStream <<
"#" << axisLabel_ <<
"\tcharge\n";
208 for (
unsigned int i = 0; i < chargeZ_.size(); ++i) {
209 RealType z = zAve * (i + 0.5) / chargeZ_.size();
211 RealType volSlice = areaAve * zAve / zBox_.size();
213 binCharge = chargeZ_[i] / (volSlice * nProcessed_);
215 rdfStream << z <<
"\t" << binCharge <<
"\n";
219 snprintf(painCave.errMsg, MAX_SIM_ERROR_MSG_LENGTH,
220 "ChargeZ: unable to open %s\n", outputFilename_.c_str());
221 painCave.isFatal = 1;
AtomType is what OpenMD looks to for unchanging data about an atom.
bool isDynamic()
Tests if the result from evaluation of script is dynamic.
StuntDouble * nextSelected(int &i)
Finds the next selected StuntDouble in the selection.
StuntDouble * beginSelected(int &i)
Finds the first selected StuntDouble in the selection.
One of the heavy-weight classes of OpenMD, SimInfo maintains objects and variables relating to the cu...
SnapshotManager * getSnapshotManager()
Returns the snapshot manager.
Mat3x3d getHmat()
Returns the H-Matrix.
void wrapVector(Vector3d &v)
Wrapping the vector according to periodic boundary condition.
Snapshot * getCurrentSnapshot()
Returns the pointer of current snapshot.
"Don't move, or you're dead! Stand up! Captain, we've got them!"
Vector3d getPos()
Returns the current position of this stuntDouble.
void setPos(const Vector3d &pos)
Sets the current position of this stuntDouble.
This basic Periodic Table class was originally taken from the data.cpp file in OpenBabel.
std::string getPrefix(const std::string &str)