Simpatico  v1.10
List of all members | Public Member Functions | Friends
Simp::OrthorhombicBoundary Class Reference

Detailed Description

An orthorhombic periodic unit cell.

Definition at line 38 of file OrthorhombicBoundary.h.

#include <OrthorhombicBoundary.h>

Inheritance diagram for Simp::OrthorhombicBoundary:
Simp::OrthoRegion

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 Vectorlengths () 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 VectorbravaisBasisVector (int i) const
 Return Bravais lattice vector i. More...
 
const VectorreciprocalBasisVector (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...
 

Constructor & Destructor Documentation

Simp::OrthorhombicBoundary::OrthorhombicBoundary ( )

Member Function Documentation

void Simp::OrthorhombicBoundary::setOrthorhombic ( const Vector lengths)
void Simp::OrthorhombicBoundary::setTetragonal ( double  a,
double  bc 
)

Set unit cell dimensions for tetragonal boundary.

Parameters
aunit cell dimensions in unique x-direction
bcunit 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.

Parameters
aunit 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().

template<class Archive >
void Simp::OrthorhombicBoundary::serialize ( Archive &  ar,
const unsigned int  version 
)

Serialize to/from an archive.

Parameters
arsaving or loading archive
versionarchive version id

Definition at line 645 of file OrthorhombicBoundary.h.

References isValid(), Simp::OrthoRegion::serialize(), and Util::serializeEnum().

void Simp::OrthorhombicBoundary::shift ( Vector r) const
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]

Parameters
rVector 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().

void Simp::OrthorhombicBoundary::shift ( Vector r,
IntVector shift 
) const
inline

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.

See also
Atom:shift()
Parameters
rVector of Cartesian coordinates
shiftinteger 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_.

void Simp::OrthorhombicBoundary::applyShift ( Vector r,
int  i,
int  t 
) const
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.

Parameters
rCartesian position Vector
idirection index
tmultiple of Bravais lattice vector i

Definition at line 468 of file OrthorhombicBoundary.h.

References Simp::OrthoRegion::lengths_.

Referenced by DdMd::Exchanger::update().

void Simp::OrthorhombicBoundary::shiftGen ( Vector r) const
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

Parameters
rVector 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().

void Simp::OrthorhombicBoundary::shiftGen ( Vector r,
IntVector shift 
) const
inline

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.

See also
Atom:shift()
Parameters
rVector of generalized coordinates (in/out)
shiftinteger shifts, modified on output (in/out)

Definition at line 492 of file OrthorhombicBoundary.h.

References Util::Dimension.

double Simp::OrthorhombicBoundary::distanceSq ( const Vector r1,
const Vector r2 
) const
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.

Parameters
r1first Cartesian position Vector
r2second Cartesian position Vector
Returns
square of distance between r1 and r2, using nearest image.

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().

double Simp::OrthorhombicBoundary::distanceSq ( const Vector r1,
const Vector r2,
IntVector shift 
) const
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.

Parameters
r1first Cartesian position Vector
r2second Cartesian position Vector
shiftshift added to r1 to create nearest image of r2.
Returns
square of distance between r1 and r2, using nearest image.

Definition at line 512 of file OrthorhombicBoundary.h.

References Util::Dimension, Simp::OrthoRegion::halfLengths_, and Simp::OrthoRegion::lengths_.

double Simp::OrthorhombicBoundary::distanceSq ( const Vector r1,
const Vector r2,
Vector dr 
) const
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.

Parameters
r1first Cartesian position Vector
r2second Cartesian position Vector
drseparation Vector (upon return)
Returns
square of separation Vector dr

Definition at line 564 of file OrthorhombicBoundary.h.

References Util::Dimension, Simp::OrthoRegion::halfLengths_, and Simp::OrthoRegion::lengths_.

bool Simp::OrthorhombicBoundary::isMinImageGen ( const Vector dr)
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.

Parameters
drseparation vector in generalized coordinates.
Returns
true if minimum image, false otherwise

Definition at line 615 of file OrthorhombicBoundary.h.

References Util::Dimension.

bool Simp::OrthorhombicBoundary::isMinImageCart ( const Vector dr)
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.

Parameters
drseparation vector in Cartesian coordinates.
Returns
true if minimum image, false otherwise

Definition at line 629 of file OrthorhombicBoundary.h.

References Util::Dimension, and Simp::OrthoRegion::halfLengths_.

void Simp::OrthorhombicBoundary::transformCartToGen ( const Vector Rc,
Vector Rg 
) const
inline
void Simp::OrthorhombicBoundary::transformGenToCart ( const Vector Rg,
Vector Rc 
) const
inline
LatticeSystem Simp::OrthorhombicBoundary::latticeSystem ( )
inline

Return actual lattice system.

Value can be Cubic, Tetragonal, or Orthorhombic.

Definition at line 584 of file OrthorhombicBoundary.h.

const Vector & Simp::OrthorhombicBoundary::lengths ( ) const
inline
double Simp::OrthorhombicBoundary::length ( int  i) const
inline

Get length in Cartesian direction i.

Parameters
iindex 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().

double Simp::OrthorhombicBoundary::minLength ( ) const
inline

Get minimum length across primitive unit cell.

Definition at line 408 of file OrthorhombicBoundary.h.

double Simp::OrthorhombicBoundary::volume ( ) const
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().

const Vector & Simp::OrthorhombicBoundary::bravaisBasisVector ( int  i) const
inline

Return Bravais lattice vector i.

Parameters
ibasis Vector index.

Definition at line 420 of file OrthorhombicBoundary.h.

Referenced by McMd::MdEwaldPotential::makeWaves().

const Vector & Simp::OrthorhombicBoundary::reciprocalBasisVector ( int  i) const
inline
void Simp::OrthorhombicBoundary::randomPosition ( Random random,
Vector r 
) const

Generate random position within the primary unit cell.

Generates Vector r[i], i=1,..,3 with minima_[i] < r[i] < maxima_[i].

Parameters
randomrandom number generator object
rVector 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().

Friends And Related Function Documentation

friend class ::OrthorhombicBoundaryTest
friend

Unit test.

Definition at line 353 of file OrthorhombicBoundary.h.

std::istream& operator>> ( std::istream &  in,
Simp::OrthorhombicBoundary boundary 
)
friend

istream extractor

Parameters
ininput stream
boundaryOrthorhombicBoundary to be read from stream
Returns
modified input stream

Definition at line 157 of file OrthorhombicBoundary.cpp.

std::ostream& operator<< ( std::ostream &  out,
const Simp::OrthorhombicBoundary boundary 
)
friend

ostream inserter

Parameters
outoutput stream
boundaryOrthorhombicBoundary to be written to stream
Returns
modified output stream

Definition at line 189 of file OrthorhombicBoundary.cpp.


The documentation for this class was generated from the following files: