Simpatico
v1.10
|
An orthorhombic periodic unit cell.
Definition at line 38 of file OrthorhombicBoundary.h.
#include <OrthorhombicBoundary.h>
Public Member Functions | |
OrthorhombicBoundary () | |
Constructor. More... | |
void | setOrthorhombic (const Vector &lengths) |
Set unit cell dimensions for orthorhombic boundary. More... | |
void | setTetragonal (double a, double bc) |
Set unit cell dimensions for tetragonal boundary. More... | |
void | setCubic (double a) |
Set unit cell dimensions for a cubic boundary. More... | |
template<class Archive > | |
void | serialize (Archive &ar, const unsigned int version) |
Serialize to/from an archive. More... | |
Periodic Boundary Conditions - Primitive Cell | |
void | shift (Vector &r) const |
Shift Cartesian Vector r to its primary image. More... | |
void | shift (Vector &r, IntVector &shift) const |
Shift Cartesian Vector r to its primary image. More... | |
void | applyShift (Vector &r, int i, int t) const |
Shift Cartesian Vector r by multiple t of a Bravais lattice vector. More... | |
void | shiftGen (Vector &r) const |
Shift generalized Vector r to its primary image. More... | |
void | shiftGen (Vector &r, IntVector &shift) const |
Shift generalized Vector r to its image within the primary unit cell. More... | |
Minimum Image Pair Separations for Cartesian Vectors | |
double | distanceSq (const Vector &r1, const Vector &r2) const |
Return square distance between positions r1 and r2. More... | |
double | distanceSq (const Vector &r1, const Vector &r2, IntVector &shift) const |
Return square distance between positions r1 and r2. More... | |
double | distanceSq (const Vector &r1, const Vector &r2, Vector &dr) const |
Compute distance and separation between r1 and r2. More... | |
Minimum Image Vector Tests | |
bool | isMinImageGen (const Vector &dr) |
Is a generalized separation vector a minimimum image of itself? More... | |
bool | isMinImageCart (const Vector &dr) |
Is a Cartesian separation vector a minimimum image of itself? More... | |
Coordinate Transformations | |
void | transformCartToGen (const Vector &Rc, Vector &Rg) const |
Transform Cartesian Vector to scaled / generalized coordinates. More... | |
void | transformGenToCart (const Vector &Rg, Vector &Rc) const |
Transform Vector of generalized coordinates to Cartesian Vector. More... | |
Accessors | |
LatticeSystem | latticeSystem () |
Return actual lattice system. More... | |
const Vector & | lengths () const |
Get Vector of unit cell lengths by const reference. More... | |
double | length (int i) const |
Get length in Cartesian direction i. More... | |
double | minLength () const |
Get minimum length across primitive unit cell. More... | |
double | volume () const |
Return unit cell volume. More... | |
const Vector & | bravaisBasisVector (int i) const |
Return Bravais lattice vector i. More... | |
const Vector & | reciprocalBasisVector (int i) const |
Return reciprocal lattice basis vector i. More... | |
void | randomPosition (Random &random, Vector &r) const |
Generate random position within the primary unit cell. More... | |
bool | isValid () |
Return true if valid, or throw Exception. More... | |
Friends | |
class | ::OrthorhombicBoundaryTest |
Unit test. More... | |
std::istream & | operator>> (std::istream &in, Simp::OrthorhombicBoundary &boundary) |
istream extractor More... | |
std::ostream & | operator<< (std::ostream &out, const Simp::OrthorhombicBoundary &boundary) |
ostream inserter More... | |
Simp::OrthorhombicBoundary::OrthorhombicBoundary | ( | ) |
Constructor.
Definition at line 25 of file OrthorhombicBoundary.cpp.
References Util::Dimension, Simp::OrthoRegion::lengths_, Util::Constants::Pi, and Util::Vector::Zero.
void Simp::OrthorhombicBoundary::setOrthorhombic | ( | const Vector & | lengths | ) |
Set unit cell dimensions for orthorhombic boundary.
Also sets all related lengths and volume.
lengths | Vector of unit cell lengths |
Definition at line 47 of file OrthorhombicBoundary.cpp.
References lengths(), and Simp::OrthoRegion::maxima_.
Referenced by DdMd::NphIntegrator::integrateStep1(), DdMd::NptIntegrator::integrateStep1(), McMd::HybridNphMdMove::move(), McMd::ReplicaMove::move(), McMd::LammpsConfigIo::read(), DdMd::LammpsConfigIo::readConfig(), Tools::HoomdConfigReader::readConfig(), McMd::LammpsDumpReader::readFrame(), McMd::DCDTrajectoryReader::readFrame(), and McMd::NphIntegrator::step().
void Simp::OrthorhombicBoundary::setTetragonal | ( | double | a, |
double | bc | ||
) |
Set unit cell dimensions for tetragonal boundary.
a | unit cell dimensions in unique x-direction |
bc | unit cell length in y and z directions |
Definition at line 57 of file OrthorhombicBoundary.cpp.
References Simp::OrthoRegion::maxima_.
void Simp::OrthorhombicBoundary::setCubic | ( | double | a | ) |
Set unit cell dimensions for a cubic boundary.
a | unit cell length in x, y, and z directions. |
Definition at line 69 of file OrthorhombicBoundary.cpp.
References Util::Dimension, Simp::OrthoRegion::lengths_, Simp::OrthoRegion::maxima_, Util::Constants::Pi, and Simp::OrthoRegion::resetRegion().
void Simp::OrthorhombicBoundary::serialize | ( | Archive & | ar, |
const unsigned int | version | ||
) |
Serialize to/from an archive.
ar | saving or loading archive |
version | archive version id |
Definition at line 645 of file OrthorhombicBoundary.h.
References isValid(), Simp::OrthoRegion::serialize(), and Util::serializeEnum().
|
inline |
Shift Cartesian Vector r to its primary image.
Upon completion, each component r[i] of Vector r is shifted by a multiple of length[i] so as to lie within the allowed range minima_[i] <= r[i] < maxima_[i].
Precondition: The algorithm assumes that on input, for each i, minima_[i] - lengths_[i] < r[i] < maxima_[i] + lengths_[i]
r | Vector of Cartesian coordinates |
Definition at line 432 of file OrthorhombicBoundary.h.
References Util::Dimension, Simp::OrthoRegion::lengths_, Simp::OrthoRegion::maxima_, and Simp::OrthoRegion::minima_.
Referenced by McMd::CfbLinear::addAtom(), McMd::CfbEndBase::addEndAtom(), McMd::CfbRebridgeBase::addMiddleAtom(), McMd::Generator::attemptPlaceAtom(), McMd::McPairPotential::buildCellList(), McMd::MdPairPotential::buildPairList(), McMd::Cluster::clusterCOM(), McMd::CfbLinear::deleteAtom(), McMd::CfbEndBase::deleteEndAtom(), McMd::CfbRebridgeBase::deleteMiddleAtom(), McMd::ClusterIdentifier::identifyClusters(), McMd::McNVTChemicalPotential::load(), McMd::System::loadConfig(), McMd::AtomDisplaceMove::move(), McMd::RigidDisplaceMove::move(), McMd::SmpConfigIo::read(), McMd::MdConfigIo::readAtom(), McMd::McConfigIo::readAtom(), Tools::HoomdConfigReader::readConfig(), McMd::LammpsDumpReader::readFrame(), McMd::DCDTrajectoryReader::readFrame(), Tools::ChainMaker::readParam(), McMd::Crosslinker::sample(), McMd::G1MSD::sample(), McMd::McNVTChemicalPotential::sample(), Tools::AtomMSD::sample(), McMd::AtomMSD::sample(), McMd::ComMSD::sample(), McMd::System::saveConfig(), McMd::G1MSD::setup(), Tools::AtomMSD::setup(), McMd::AtomMSD::setup(), McMd::ComMSD::setup(), McMd::MdSystem::shiftAtoms(), and Tools::HoomdConfigWriter::writeConfig().
Shift Cartesian Vector r to its primary image.
This function maps an atomic position to its primary image, and also increments the atomic shift IntVector:
If r[i] -> r[i] - t*length_[i], then shift[i] -> shift[i] + t.
r | Vector of Cartesian coordinates |
shift | integer shifts required to obtain "true" coordinates. |
Definition at line 449 of file OrthorhombicBoundary.h.
References Util::Dimension, Simp::OrthoRegion::lengths_, Simp::OrthoRegion::maxima_, and Simp::OrthoRegion::minima_.
|
inline |
Shift Cartesian Vector r by multiple t of a Bravais lattice vector.
This function shifts the Vector r by a specified amount:
r -> r + t*a'[i]
where a[i] is Bravais lattice vector number i.
r | Cartesian position Vector |
i | direction index |
t | multiple of Bravais lattice vector i |
Definition at line 468 of file OrthorhombicBoundary.h.
References Simp::OrthoRegion::lengths_.
Referenced by DdMd::Exchanger::update().
|
inline |
Shift generalized Vector r to its primary image.
One output, each coordinate r[i] is shifted by an integer, so as to lie within the range 0 < r[i] < 1.0
Precondition: The algorithm assumes that on input, for each i=0,..,2, -1.0 < r[i] < 2.0
r | Vector of scaled / generalized coordinates |
Definition at line 474 of file OrthorhombicBoundary.h.
References Util::Dimension.
Referenced by DdMd::AtomDistributor::addAtom(), McMd::MdSpmePotential::addForces(), and McMd::MdSpmePotential::makeWaves().
Shift generalized Vector r to its image within the primary unit cell.
This function maps a vector of generalized coordinates to lie in the primary cell, and also increments the atomic shift IntVector:
If r[i] -> r[i] - t, then shift[i] -> shift[i] + t.
r | Vector of generalized coordinates (in/out) |
shift | integer shifts, modified on output (in/out) |
Definition at line 492 of file OrthorhombicBoundary.h.
References Util::Dimension.
Return square distance between positions r1 and r2.
This function returns the square distance between two Cartesian vectors, using the nearest image convention for the separation Vector.
r1 | first Cartesian position Vector |
r2 | second Cartesian position Vector |
Definition at line 540 of file OrthorhombicBoundary.h.
References Util::Dimension, Simp::OrthoRegion::halfLengths_, and Simp::OrthoRegion::lengths_.
Referenced by McMd::CfbLinear::addAtom(), McMd::CfbEndBase::addEndAtom(), McMd::BondPotentialImpl< Interaction >::addForces(), McMd::DihedralPotentialImpl< Interaction >::addForces(), McMd::AnglePotentialImpl< Interaction >::addForces(), McMd::CfbRebridgeBase::addMiddleAtom(), McMd::LinkPotentialImpl< Interaction >::atomEnergy(), McMd::BondPotentialImpl< Interaction >::atomEnergy(), McMd::DihedralPotentialImpl< Interaction >::atomEnergy(), McMd::AnglePotentialImpl< Interaction >::atomEnergy(), McMd::Generator::attemptPlaceAtom(), McMd::PairList::build(), McMd::Cluster::clusterCOM(), DdMd::BondTensorAutoCorr::computeData(), McMd::LinkPotentialImpl< Interaction >::computeEnergy(), McMd::McPairPotentialImpl< Interaction >::computeEnergy(), DdMd::AnglePotentialImpl< Interaction >::computeEnergy(), McMd::MdEwaldPairPotentialImpl< Interaction >::computeEnergy(), DdMd::BondPotentialImpl< Interaction >::computeEnergy(), McMd::DihedralPotentialImpl< Interaction >::computeEnergy(), McMd::AnglePotentialImpl< Interaction >::computeEnergy(), DdMd::DihedralPotentialImpl< Interaction >::computeEnergy(), DdMd::AnglePotentialImpl< Interaction >::computeForces(), DdMd::BondPotentialImpl< Interaction >::computeForces(), DdMd::DihedralPotentialImpl< Interaction >::computeForces(), DdMd::BondPotentialImpl< Interaction >::computeForcesAndStress(), DdMd::AnglePotentialImpl< Interaction >::computeStress(), DdMd::BondPotentialImpl< Interaction >::computeStress(), DdMd::DihedralPotentialImpl< Interaction >::computeStress(), McMd::CfbLinear::deleteAtom(), McMd::CfbEndBase::deleteEndAtom(), McMd::CfbRebridgeBase::deleteMiddleAtom(), McMd::ClusterIdentifier::initialize(), McMd::PairList::isCurrent(), McMd::McPairPotentialImpl< Interaction >::moleculeEnergy(), McMd::Cluster::momentTensor(), McMd::Sliplinker::move(), McMd::SliplinkerEnd::move(), McMd::SliplinkerAll::move(), McMd::GcSliplinkMove::move(), McMd::SliplinkMove::move(), McMd::CfbDoubleRebridgeMove::move(), McMd::GroupRebridgeBase::octaEnergy(), McMd::ComMSD::output(), McMd::Crosslinker::sample(), McMd::LinkLengthDist::sample(), McMd::G1MSD::sample(), Tools::AtomMSD::sample(), McMd::McMuExchange::sample(), Tools::PairEnergy::sample(), McMd::AtomMSD::sample(), McMd::BondLengthDist::sample(), McMd::MdPairEnergyCoefficients::sample(), McMd::ComMSD::sample(), McMd::RingRouseAutoCorr::sample(), McMd::EndtoEndXYZ::sample(), McMd::LinearRouseAutoCorr::sample(), McMd::EndtoEnd::sample(), McMd::BlockRadiusGyration::sample(), McMd::RDF::sample(), McMd::IntraBondTensorAutoCorr< SystemType >::sample(), McMd::IntraPairAutoCorr::sample(), McMd::RadiusGyration::sample(), McMd::RingOctaRebridgeMove::scanBridge(), McMd::RingTetraRebridgeMove::scanBridge(), McMd::NvtDpdVvIntegrator::setup(), and McMd::GroupRebridgeBase::tetraEnergy().
|
inline |
Return square distance between positions r1 and r2.
This function returns the square distance between two Cartesian vectors, using the nearest image convention for the separation Vector. On return, the IntVector shift contains the coefficients of the Bravais lattice vectors that were added to r1 to create a nearest image of r2.
r1 | first Cartesian position Vector |
r2 | second Cartesian position Vector |
shift | shift added to r1 to create nearest image of r2. |
Definition at line 512 of file OrthorhombicBoundary.h.
References Util::Dimension, Simp::OrthoRegion::halfLengths_, and Simp::OrthoRegion::lengths_.
|
inline |
Compute distance and separation between r1 and r2.
This function returns the square distance between two Cartesian vectors, using the nearest image convention for the separation Vector. Upon return, Vector dr contains the separation r1 - r2 computed using the nearest image convention, and the return value is the square of dr.
r1 | first Cartesian position Vector |
r2 | second Cartesian position Vector |
dr | separation Vector (upon return) |
Definition at line 564 of file OrthorhombicBoundary.h.
References Util::Dimension, Simp::OrthoRegion::halfLengths_, and Simp::OrthoRegion::lengths_.
|
inline |
Is a generalized separation vector a minimimum image of itself?
This function returns true if the scaled separation dr has a Cartesian norm less than or equal to that of any vector that can be produced by adding a Bravais vector to the corresponding Cartesian vector, or false otherwise.
dr | separation vector in generalized coordinates. |
Definition at line 615 of file OrthorhombicBoundary.h.
References Util::Dimension.
|
inline |
Is a Cartesian separation vector a minimimum image of itself?
This function returns true if the scaled separation dr has a Cartesian norm less than or equal to that of any vector that can be produced by adding a Bravais vector to dr, false otherwise.
dr | separation vector in Cartesian coordinates. |
Definition at line 629 of file OrthorhombicBoundary.h.
References Util::Dimension, and Simp::OrthoRegion::halfLengths_.
Transform Cartesian Vector to scaled / generalized coordinates.
Generalized coordinates range from 0.0 < Rg[i] < 1.0 within the primitive cell, for i=0,..,2.
Rc | Vector of Cartesian coordinates (input) |
Rg | Vector of generalized coordinates (output) |
Definition at line 591 of file OrthorhombicBoundary.h.
References Util::Dimension.
Referenced by McMd::MdEwaldPotential::addForces(), McMd::MdSpmePotential::addForces(), McMd::CellList::cellIndexFromPosition(), McMd::DeformCommand::execute(), DdMd::SerializeConfigIo::loadConfig(), McMd::MdEwaldPotential::makeWaves(), McMd::MdSpmePotential::makeWaves(), DdMd::DdMdConfigIo::readConfig(), DdMd::LammpsConfigIo::readConfig(), DdMd::DdMdOrderedConfigIo::readConfig(), McMd::McCommandManager::readStandardCommand(), DdMd::AtomStorage::transformCartToGen(), Tools::HoomdConfigWriter::writeConfig(), DdMd::DdMdTrajectoryWriter::writeFrame(), and DdMd::DdMdGroupTrajectoryWriter::writeFrame().
Transform Vector of generalized coordinates to Cartesian Vector.
Rg | Vector of generalized coordinates (input) |
Rc | Vector of Cartesian coordinates (output) |
Definition at line 603 of file OrthorhombicBoundary.h.
References Util::Dimension, and Simp::OrthoRegion::lengths_.
Referenced by McMd::DeformCommand::execute(), McMd::DdMdTrajectoryReader::readFrame(), McMd::McCommandManager::readStandardCommand(), DdMd::SerializeConfigIo::saveConfig(), DdMd::AtomStorage::transformGenToCart(), Tools::HoomdConfigWriter::writeConfig(), DdMd::DdMdConfigIo::writeConfig(), DdMd::LammpsConfigIo::writeConfig(), DdMd::DdMdOrderedConfigIo::writeConfig(), and DdMd::LammpsDumpWriter::writeFrame().
|
inline |
Return actual lattice system.
Value can be Cubic, Tetragonal, or Orthorhombic.
Definition at line 584 of file OrthorhombicBoundary.h.
|
inline |
Get Vector of unit cell lengths by const reference.
Definition at line 396 of file OrthorhombicBoundary.h.
References Simp::OrthoRegion::lengths_.
Referenced by Util::bcast< Simp::OrthorhombicBoundary >(), DdMd::Exchanger::exchange(), DdMd::NphIntegrator::integrateStep1(), DdMd::NptIntegrator::integrateStep1(), McMd::HybridNphMdMove::move(), McMd::ReplicaMove::move(), McMd::VelProf::output(), Util::recv< Simp::OrthorhombicBoundary >(), McMd::G1MSD::sample(), Tools::AtomMSD::sample(), McMd::VelProf::sample(), McMd::CompositionProfile::sample(), McMd::BoundaryAverage::sample(), McMd::AtomMSD::sample(), DdMd::OutputBoxdim::sample(), McMd::ComMSD::sample(), DdMd::OrderParamNucleation::sample(), Util::send< Simp::OrthorhombicBoundary >(), setOrthorhombic(), McMd::VelProf::setup(), Tools::PairEnergy::setup(), McMd::CellList::setup(), McMd::NphIntegrator::step(), McMd::LammpsConfigIo::write(), Tools::HoomdConfigWriter::writeConfig(), DdMd::LammpsConfigIo::writeConfig(), DdMd::LammpsDumpWriter::writeFrame(), and Tools::LammpsDumpWriter::~LammpsDumpWriter().
|
inline |
Get length in Cartesian direction i.
i | index of Cartesian direction, 0 <= i < Dimension |
Definition at line 402 of file OrthorhombicBoundary.h.
References Simp::OrthoRegion::lengths_.
Referenced by DdMd::PairPotential::buildCellList(), and DdMd::PairPotential::save().
|
inline |
Get minimum length across primitive unit cell.
Definition at line 408 of file OrthorhombicBoundary.h.
|
inline |
Return unit cell volume.
Definition at line 414 of file OrthorhombicBoundary.h.
References Simp::OrthoRegion::volume_.
Referenced by McMd::MdEwaldPotential::addForces(), McMd::MdSpmePotential::addForces(), McMd::BondPotentialImpl< Interaction >::addForces(), McMd::DihedralPotentialImpl< Interaction >::addForces(), McMd::AnglePotentialImpl< Interaction >::addForces(), McMd::MdEwaldPotential::computeEnergy(), McMd::MdSpmePotential::computeEnergy(), McMd::LinkPotentialImpl< Interaction >::computeEnergy(), McMd::McPairPotentialImpl< Interaction >::computeEnergy(), McMd::MdPairPotentialImpl< Interaction >::computeEnergy(), DdMd::BondPotentialImpl< Interaction >::computeForcesAndStress(), DdMd::PairPotentialImpl< Interaction >::computeForcesAndStress(), McMd::MdEwaldPotential::computeStress(), McMd::MdSpmePotential::computeStress(), DdMd::AnglePotentialImpl< Interaction >::computeStress(), McMd::MdEwaldPairPotentialImpl< Interaction >::computeStress(), DdMd::BondPotentialImpl< Interaction >::computeStress(), DdMd::DihedralPotentialImpl< Interaction >::computeStress(), DdMd::PairPotentialImpl< Interaction >::computeStress(), DdMd::NphIntegrator::integrateStep1(), DdMd::NptIntegrator::integrateStep1(), DdMd::NphIntegrator::integrateStep2(), DdMd::NptIntegrator::integrateStep2(), McMd::MdSystem::kineticEnergy(), McMd::HybridNphMdMove::move(), DdMd::OrderParamNucleation::output(), McMd::BoundaryAverage::sample(), DdMd::OutputBoxdim::sample(), DdMd::StressAutoCorr::sample(), McMd::RDF::sample(), McMd::StructureFactor::sample(), DdMd::StructureFactorGrid::sample(), McMd::StructureFactorGrid::sample(), DdMd::StructureFactor::sample(), McMd::StructureFactorP::sample(), McMd::IntraStructureFactor::sample(), and McMd::McSystem::unsetPotentialEnergies().
|
inline |
Return Bravais lattice vector i.
i | basis Vector index. |
Definition at line 420 of file OrthorhombicBoundary.h.
Referenced by McMd::MdEwaldPotential::makeWaves().
|
inline |
Return reciprocal lattice basis vector i.
i | basis Vector index. |
Definition at line 426 of file OrthorhombicBoundary.h.
Referenced by McMd::MdEwaldPotential::addForces(), McMd::MdEwaldPotential::computeStress(), McMd::MdEwaldPotential::makeWaves(), McMd::MdSpmePotential::makeWaves(), McMd::CompositionProfile::makeWaveVectors(), McMd::StructureFactor::makeWaveVectors(), DdMd::VanHove::makeWaveVectors(), DdMd::StructureFactor::makeWaveVectors(), McMd::StructureFactorP::makeWaveVectors(), McMd::VanHove::sample(), and McMd::IntraStructureFactor::sample().
Generate random position within the primary unit cell.
Generates Vector r[i], i=1,..,3 with minima_[i] < r[i] < maxima_[i].
random | random number generator object |
r | Vector of random coordinates (upon return) |
Definition at line 105 of file OrthorhombicBoundary.cpp.
References Util::Dimension, Simp::OrthoRegion::maxima_, Simp::OrthoRegion::minima_, and Util::Random::uniform().
Referenced by McMd::PointGenerator::attemptPlaceMolecule(), McMd::LinearGenerator::attemptPlaceMolecule(), Tools::ChainMaker::readParam(), and McMd::McNVTChemicalPotential::sample().
bool Simp::OrthorhombicBoundary::isValid | ( | ) |
Return true if valid, or throw Exception.
Definition at line 115 of file OrthorhombicBoundary.cpp.
References Util::Dimension, Util::feq(), Simp::OrthoRegion::isValid(), Simp::OrthoRegion::lengths_, Simp::OrthoRegion::minima_, Util::Constants::Pi, and UTIL_THROW.
Referenced by serialize().
|
friend |
Unit test.
Definition at line 353 of file OrthorhombicBoundary.h.
|
friend |
istream extractor
in | input stream |
boundary | OrthorhombicBoundary to be read from stream |
Definition at line 157 of file OrthorhombicBoundary.cpp.
|
friend |
ostream inserter
out | output stream |
boundary | OrthorhombicBoundary to be written to stream |
Definition at line 189 of file OrthorhombicBoundary.cpp.