| 1 |  | /* | 
| 2 | < | * Copyright (C) 2000-2004  Object Oriented Parallel Simulation Engine (OOPSE) project | 
| 3 | < | * | 
| 4 | < | * Contact: oopse@oopse.org | 
| 5 | < | * | 
| 6 | < | * This program is free software; you can redistribute it and/or | 
| 7 | < | * modify it under the terms of the GNU Lesser General Public License | 
| 8 | < | * as published by the Free Software Foundation; either version 2.1 | 
| 9 | < | * of the License, or (at your option) any later version. | 
| 10 | < | * All we ask is that proper credit is given for our work, which includes | 
| 11 | < | * - but is not limited to - adding the above copyright notice to the beginning | 
| 12 | < | * of your source code files, and to any copyright notice that you may distribute | 
| 13 | < | * with programs based on this work. | 
| 14 | < | * | 
| 15 | < | * This program is distributed in the hope that it will be useful, | 
| 16 | < | * but WITHOUT ANY WARRANTY; without even the implied warranty of | 
| 17 | < | * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the | 
| 18 | < | * GNU Lesser General Public License for more details. | 
| 19 | < | * | 
| 20 | < | * You should have received a copy of the GNU Lesser General Public License | 
| 21 | < | * along with this program; if not, write to the Free Software | 
| 22 | < | * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA. | 
| 2 | > | * Copyright (c) 2005 The University of Notre Dame. All Rights Reserved. | 
| 3 |  | * | 
| 4 | + | * The University of Notre Dame grants you ("Licensee") a | 
| 5 | + | * non-exclusive, royalty free, license to use, modify and | 
| 6 | + | * redistribute this software in source and binary code form, provided | 
| 7 | + | * that the following conditions are met: | 
| 8 | + | * | 
| 9 | + | * 1. Redistributions of source code must retain the above copyright | 
| 10 | + | *    notice, this list of conditions and the following disclaimer. | 
| 11 | + | * | 
| 12 | + | * 2. Redistributions in binary form must reproduce the above copyright | 
| 13 | + | *    notice, this list of conditions and the following disclaimer in the | 
| 14 | + | *    documentation and/or other materials provided with the | 
| 15 | + | *    distribution. | 
| 16 | + | * | 
| 17 | + | * This software is provided "AS IS," without a warranty of any | 
| 18 | + | * kind. All express or implied conditions, representations and | 
| 19 | + | * warranties, including any implied warranty of merchantability, | 
| 20 | + | * fitness for a particular purpose or non-infringement, are hereby | 
| 21 | + | * excluded.  The University of Notre Dame and its licensors shall not | 
| 22 | + | * be liable for any damages suffered by licensee as a result of | 
| 23 | + | * using, modifying or distributing the software or its | 
| 24 | + | * derivatives. In no event will the University of Notre Dame or its | 
| 25 | + | * licensors be liable for any lost revenue, profit or data, or for | 
| 26 | + | * direct, indirect, special, consequential, incidental or punitive | 
| 27 | + | * damages, however caused and regardless of the theory of liability, | 
| 28 | + | * arising out of the use of or inability to use software, even if the | 
| 29 | + | * University of Notre Dame has been advised of the possibility of | 
| 30 | + | * such damages. | 
| 31 | + | * | 
| 32 | + | * SUPPORT OPEN SCIENCE!  If you use OpenMD or its source code in your | 
| 33 | + | * research, please cite the appropriate papers when you publish your | 
| 34 | + | * work.  Good starting points are: | 
| 35 | + | * | 
| 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, 234107 (2008). | 
| 39 | + | * [4]  Kuang & Gezelter,  J. Chem. Phys. 133, 164101 (2010). | 
| 40 | + | * [5]  Vardeman, Stocker & Gezelter, J. Chem. Theory Comput. 7, 834 (2011). | 
| 41 |  | */ | 
| 42 | < |  | 
| 26 | < | /** | 
| 27 | < | * @file Snapshot.hpp | 
| 28 | < | * @author tlin | 
| 29 | < | * @date 10/20/2004 | 
| 30 | < | * @time 23:56am | 
| 31 | < | * @version 1.0 | 
| 32 | < | */ | 
| 42 | > |  | 
| 43 |  | #ifndef BRAINS_SNAPSHOT_HPP | 
| 44 |  | #define BRAINS_SNAPSHOT_HPP | 
| 45 |  |  | 
| 46 |  | #include <vector> | 
| 47 |  |  | 
| 48 | < | #include "math/Vector3.hpp" | 
| 49 | < | #include "math/SquareMatrix3.hpp" | 
| 48 | > | #include "brains/DataStorage.hpp" | 
| 49 | > | #include "nonbonded/NonBondedInteraction.hpp" | 
| 50 | > | #include "brains/Stats.hpp" | 
| 51 |  |  | 
| 52 |  | using namespace std; | 
| 53 | + | namespace OpenMD{ | 
| 54 |  |  | 
| 55 | < | namespace oopse{ | 
| 55 | > | /** | 
| 56 | > | * FrameData is a structure for holding system-wide dynamic data | 
| 57 | > | * about the simulation. | 
| 58 | > | */ | 
| 59 | > |  | 
| 60 | > | struct FrameData { | 
| 61 | > | int id;                       /**< identification number of the snapshot */ | 
| 62 | > | RealType currentTime;         /**< current time */ | 
| 63 | > | Mat3x3d  hmat;                /**< axes of the periodic box in matrix form */ | 
| 64 | > | Mat3x3d  invHmat;             /**< the inverse of the Hmat matrix */ | 
| 65 | > | Mat3x3d  bBox;                /**< axes of a bounding box in matrix form */ | 
| 66 | > | Mat3x3d  invBbox;             /**< the inverse of the bounding box */ | 
| 67 | > | bool     orthoRhombic;        /**< is this an orthorhombic periodic box? */ | 
| 68 | > | RealType totalEnergy;         /**< total energy of this frame */ | 
| 69 | > | RealType translationalKinetic; /**< translational kinetic energy of this frame */ | 
| 70 | > | RealType rotationalKinetic;   /**< rotational kinetic energy of this frame */ | 
| 71 | > | RealType kineticEnergy;       /**< kinetic energy of this frame */ | 
| 72 | > | RealType potentialEnergy;     /**< potential energy of this frame */ | 
| 73 | > | RealType shortRangePotential; /**< short-range contributions to the potential*/ | 
| 74 | > | RealType longRangePotential;  /**< long-range contributions to the potential */ | 
| 75 | > | RealType reciprocalPotential; /**< reciprocal-space contributions to the potential */ | 
| 76 | > | RealType bondPotential;       /**< bonded contribution to the potential */ | 
| 77 | > | RealType bendPotential;       /**< angle-bending contribution to the potential */ | 
| 78 | > | RealType torsionPotential;    /**< dihedral (torsion angle) contribution to the potential */ | 
| 79 | > | RealType inversionPotential;  /**< inversion (planarity) contribution to the potential */ | 
| 80 | > | potVec   lrPotentials;        /**< breakdown of long-range potentials by family */ | 
| 81 | > | potVec   excludedPotentials;  /**< breakdown of excluded potentials by family */ | 
| 82 | > | RealType restraintPotential;  /**< potential energy of restraints */ | 
| 83 | > | RealType rawPotential;        /**< unrestrained potential energy (when restraints are applied) */ | 
| 84 | > | RealType xyArea;              /**< XY area of this frame */ | 
| 85 | > | RealType volume;              /**< total volume of this frame */ | 
| 86 | > | RealType pressure;            /**< pressure of this frame */ | 
| 87 | > | RealType temperature;         /**< temperature of this frame */ | 
| 88 | > | pair<RealType, RealType> thermostat;    /**< thermostat variables */ | 
| 89 | > | RealType electronicTemperature; /**< temperature of the electronic degrees of freedom */ | 
| 90 | > | pair<RealType, RealType> electronicThermostat; /**< thermostat variables for electronic degrees of freedom */ | 
| 91 | > | Mat3x3d  barostat;            /**< barostat matrix */ | 
| 92 | > | Vector3d COM;                 /**< location of system center of mass */ | 
| 93 | > | Vector3d COMvel;              /**< system center of mass velocity */ | 
| 94 | > | Vector3d COMw;                /**< system center of mass angular velocity */ | 
| 95 | > | Mat3x3d  inertiaTensor;       /**< inertia tensor for entire system */ | 
| 96 | > | RealType gyrationalVolume;    /**< gyrational volume for entire system */ | 
| 97 | > | RealType hullVolume;          /**< hull volume for entire system */ | 
| 98 | > | Mat3x3d  stressTensor;        /**< stress tensor */ | 
| 99 | > | Mat3x3d  pressureTensor;      /**< pressure tensor */ | 
| 100 | > | Vector3d systemDipole;        /**< total system dipole moment */ | 
| 101 | > | Vector3d conductiveHeatFlux;  /**< heat flux vector (conductive only) */ | 
| 102 | > | Vector3d convectiveHeatFlux;  /**< heat flux vector (convective only) */ | 
| 103 | > | RealType conservedQuantity;   /**< anything conserved by the integrator */ | 
| 104 | > | }; | 
| 105 |  |  | 
| 45 | – | /** | 
| 46 | – | * @class Snapshot Snapshot.hpp "brains/Snapshot.hpp" | 
| 47 | – | * @brief Snapshot class is a repository class for storing dynamic data during | 
| 48 | – | *  Simulation | 
| 49 | – | * @see SimData | 
| 50 | – | */ | 
| 51 | – | class Snapshot { | 
| 52 | – | public: | 
| 106 |  |  | 
| 107 | < | Snapshot(int i) { | 
| 107 | > | /** | 
| 108 | > | * @class Snapshot | 
| 109 | > | * @brief The Snapshot class is a repository storing dynamic data during a | 
| 110 | > | * Simulation.  Every Snapshot contains FrameData (for global information) | 
| 111 | > | * as well as DataStorage (one for Atoms, one for RigidBodies, and one for | 
| 112 | > | * CutoffGroups). | 
| 113 | > | */ | 
| 114 | > | class Snapshot { | 
| 115 |  |  | 
| 116 | < | } | 
| 116 | > | public: | 
| 117 | > | Snapshot(int nAtoms, int nRigidbodies, int nCutoffGroups); | 
| 118 | > | Snapshot(int nAtoms, int nRigidbodies, int nCutoffGroups, int storageLayout); | 
| 119 | > | /** Returns the id of this Snapshot */ | 
| 120 | > | int      getID(); | 
| 121 | > | /** Sets the id of this Snapshot */ | 
| 122 | > | void     setID(int id); | 
| 123 |  |  | 
| 124 | < | Snapshot(const Snapshot& s); | 
| 124 | > | /** sets the state of the computed properties to false */ | 
| 125 | > | void     clearDerivedProperties(); | 
| 126 |  |  | 
| 127 | < | Snapshot& operator =(const Snapshot& s); | 
| 127 | > | int      getSize(); | 
| 128 | > | /** Returns the number of atoms */ | 
| 129 | > | int      getNumberOfAtoms(); | 
| 130 | > | /** Returns the number of rigid bodies */ | 
| 131 | > | int      getNumberOfRigidBodies(); | 
| 132 | > | /** Returns the number of rigid bodies */ | 
| 133 | > | int      getNumberOfCutoffGroups(); | 
| 134 | > | /** Returns the number of bytes in a FrameData structure */ | 
| 135 | > | static int getFrameDataSize(); | 
| 136 | > |  | 
| 137 | > | /** Returns the H-Matrix */ | 
| 138 | > | Mat3x3d  getHmat(); | 
| 139 | > | /** Sets the H-Matrix */ | 
| 140 | > | void     setHmat(const Mat3x3d& m); | 
| 141 | > | /** Returns the inverse H-Matrix */ | 
| 142 | > | Mat3x3d  getInvHmat(); | 
| 143 | > |  | 
| 144 | > | /** Returns the Bounding Box */ | 
| 145 | > | Mat3x3d  getBoundingBox(); | 
| 146 | > | /** Sets the Bounding Box */ | 
| 147 | > | void     setBoundingBox(const Mat3x3d& m); | 
| 148 | > | /** Returns the inverse Bounding Box*/ | 
| 149 | > | Mat3x3d  getInvBoundingBox(); | 
| 150 |  |  | 
| 151 | < | /** Returns the id of this Snapshot */ | 
| 152 | < | int getID() { | 
| 153 | < | return id_; | 
| 65 | < | } | 
| 151 | > | RealType getVolume(); | 
| 152 | > | RealType getXYarea(); | 
| 153 | > | void     setVolume(const RealType vol); | 
| 154 |  |  | 
| 155 | < | /** Sets the id of this Snapshot */ | 
| 156 | < | void setID(int id) { | 
| 69 | < | id_ = id; | 
| 70 | < | } | 
| 155 | > | /** Wrapping the vector according to periodic boundary condition*/ | 
| 156 | > | void     wrapVector(Vector3d& v); | 
| 157 |  |  | 
| 158 | < | /** */ | 
| 159 | < | Snapshot* clone() { | 
| 74 | < | return new Snapshot(*this); | 
| 75 | < | } | 
| 158 | > | /** Scaling a vector to multiples of the periodic box */ | 
| 159 | > | Vector3d scaleVector(Vector3d &v); | 
| 160 |  |  | 
| 161 | + | void     setCOM(const Vector3d &com); | 
| 162 | + | void     setCOMvel(const Vector3d &comVel); | 
| 163 | + | void     setCOMw(const Vector3d &comw); | 
| 164 |  |  | 
| 165 | < | //template<typename T> | 
| 166 | < | //static typename T::ElemPointerType getArrayPointer(vector<T>& v) { | 
| 167 | < | //    return v[0]->getArrayPointer(); | 
| 168 | < | //} | 
| 165 | > | Vector3d getCOM(); | 
| 166 | > | Vector3d getCOMvel(); | 
| 167 | > | Vector3d getCOMw(); | 
| 168 | > |  | 
| 169 | > | RealType getTime(); | 
| 170 | > | void     increaseTime(const RealType dt); | 
| 171 | > | void     setTime(const RealType time); | 
| 172 |  |  | 
| 173 | < | static double* getArrayPointer(vector<Vector3d>& v) { | 
| 174 | < | return v[0].getArrayPointer(); | 
| 175 | < | } | 
| 173 | > | void     setBondPotential(const RealType bp); | 
| 174 | > | void     setBendPotential(const RealType bp); | 
| 175 | > | void     setTorsionPotential(const RealType tp); | 
| 176 | > | void     setInversionPotential(const RealType ip); | 
| 177 | > | RealType getBondPotential(); | 
| 178 | > | RealType getBendPotential(); | 
| 179 | > | RealType getTorsionPotential(); | 
| 180 | > | RealType getInversionPotential(); | 
| 181 | > |  | 
| 182 | > | RealType getShortRangePotential(); | 
| 183 | > |  | 
| 184 | > | void     setLongRangePotential(const potVec lrPot); | 
| 185 | > | RealType getLongRangePotential(); | 
| 186 | > | potVec   getLongRangePotentials(); | 
| 187 | > |  | 
| 188 | > | void     setReciprocalPotential(const RealType rp); | 
| 189 | > | RealType getReciprocalPotential(); | 
| 190 | > |  | 
| 191 | > | void     setExcludedPotentials(const potVec exPot); | 
| 192 | > | potVec   getExcludedPotentials(); | 
| 193 | > |  | 
| 194 | > | void     setRestraintPotential(const RealType rp); | 
| 195 | > | RealType getRestraintPotential(); | 
| 196 | > |  | 
| 197 | > | void     setRawPotential(const RealType rp); | 
| 198 | > | RealType getRawPotential(); | 
| 199 | > |  | 
| 200 | > | RealType getPotentialEnergy(); | 
| 201 | > | RealType getKineticEnergy(); | 
| 202 | > | RealType getTranslationalKineticEnergy(); | 
| 203 | > | RealType getRotationalKineticEnergy(); | 
| 204 | > | void     setKineticEnergy(const RealType ke); | 
| 205 | > | void     setTranslationalKineticEnergy(const RealType tke); | 
| 206 | > | void     setRotationalKineticEnergy(const RealType rke); | 
| 207 | > | RealType getTotalEnergy(); | 
| 208 | > | void     setTotalEnergy(const RealType te); | 
| 209 | > | RealType getConservedQuantity(); | 
| 210 | > | void     setConservedQuantity(const RealType cq); | 
| 211 | > | RealType getTemperature(); | 
| 212 | > | void     setTemperature(const RealType temp); | 
| 213 | > | RealType getElectronicTemperature(); | 
| 214 | > | void     setElectronicTemperature(const RealType eTemp); | 
| 215 | > | RealType getPressure(); | 
| 216 | > | void     setPressure(const RealType pressure); | 
| 217 | > |  | 
| 218 | > | Mat3x3d  getPressureTensor(); | 
| 219 | > | void     setPressureTensor(const Mat3x3d& pressureTensor); | 
| 220 | > |  | 
| 221 | > | Mat3x3d  getStressTensor(); | 
| 222 | > | void     setStressTensor(const Mat3x3d& stressTensor); | 
| 223 | > |  | 
| 224 | > | Vector3d getConductiveHeatFlux(); | 
| 225 | > | void     setConductiveHeatFlux(const Vector3d& chf); | 
| 226 | > |  | 
| 227 | > | Vector3d getConvectiveHeatFlux(); | 
| 228 | > | void     setConvectiveHeatFlux(const Vector3d& chf); | 
| 229 | > |  | 
| 230 | > | Vector3d getHeatFlux(); | 
| 231 | > |  | 
| 232 | > | Vector3d getSystemDipole(); | 
| 233 | > | void     setSystemDipole(const Vector3d& bd); | 
| 234 | > |  | 
| 235 | > | pair<RealType, RealType> getThermostat(); | 
| 236 | > | void setThermostat(const pair<RealType, RealType>& thermostat); | 
| 237 | > |  | 
| 238 | > | pair<RealType, RealType> getElectronicThermostat(); | 
| 239 | > | void setElectronicThermostat(const pair<RealType, RealType>& eThermostat); | 
| 240 |  |  | 
| 241 | < | static double* getArrayPointer(vector<RotMat3x3d>& v) { | 
| 242 | < | return v[0].getArrayPointer(); | 
| 89 | < | } | 
| 90 | < |  | 
| 91 | < | static double* getArrayPointer(vector<double>& v) { | 
| 92 | < | assert(v.size() > 0); | 
| 93 | < | return &(v[0]); | 
| 94 | < | } | 
| 95 | < |  | 
| 96 | < | vector<Vector3d> pos; | 
| 97 | < | vector<Vector3d> vel; | 
| 98 | < | vector<Vector3d> frc; | 
| 99 | < | vector<Vector3d> trq; | 
| 100 | < | vector<RotMat3x3d> Amat; | 
| 101 | < | vector<Vector3d> mu; | 
| 102 | < | vector<Vector3d> ul; | 
| 103 | < | vector<double> zAngle; | 
| 104 | < |  | 
| 105 | < | private: | 
| 241 | > | Mat3x3d  getBarostat(); | 
| 242 | > | void     setBarostat(const Mat3x3d& barostat); | 
| 243 |  |  | 
| 244 | < | int id_; /**< identification number of the snapshot */ | 
| 245 | < | }; | 
| 244 | > | Mat3x3d  getInertiaTensor(); | 
| 245 | > | void     setInertiaTensor(const Mat3x3d& inertiaTensor); | 
| 246 |  |  | 
| 247 | + | RealType getGyrationalVolume(); | 
| 248 | + | void     setGyrationalVolume(const RealType gv); | 
| 249 | + |  | 
| 250 | + | RealType getHullVolume(); | 
| 251 | + | void     setHullVolume(const RealType hv); | 
| 252 | + |  | 
| 253 | + | void     setOrthoTolerance(RealType orthoTolerance); | 
| 254 | + |  | 
| 255 | + | DataStorage atomData; | 
| 256 | + | DataStorage rigidbodyData; | 
| 257 | + | DataStorage cgData; | 
| 258 | + | FrameData   frameData; | 
| 259 | + |  | 
| 260 | + | bool hasTotalEnergy; | 
| 261 | + | bool hasTranslationalKineticEnergy; | 
| 262 | + | bool hasRotationalKineticEnergy; | 
| 263 | + | bool hasKineticEnergy; | 
| 264 | + | bool hasShortRangePotential; | 
| 265 | + | bool hasLongRangePotential; | 
| 266 | + | bool hasPotentialEnergy; | 
| 267 | + | bool hasXYarea; | 
| 268 | + | bool hasVolume; | 
| 269 | + | bool hasPressure; | 
| 270 | + | bool hasTemperature; | 
| 271 | + | bool hasElectronicTemperature; | 
| 272 | + | bool hasCOM; | 
| 273 | + | bool hasCOMvel; | 
| 274 | + | bool hasCOMw; | 
| 275 | + | bool hasPressureTensor; | 
| 276 | + | bool hasSystemDipole; | 
| 277 | + | bool hasConvectiveHeatFlux; | 
| 278 | + | bool hasInertiaTensor; | 
| 279 | + | bool hasGyrationalVolume; | 
| 280 | + | bool hasHullVolume; | 
| 281 | + | bool hasConservedQuantity; | 
| 282 | + | bool hasBoundingBox; | 
| 283 | + |  | 
| 284 | + | private: | 
| 285 | + | RealType orthoTolerance_; | 
| 286 | + |  | 
| 287 | + | }; | 
| 288 | + |  | 
| 289 | + | typedef DataStorage (Snapshot::*DataStoragePointer); | 
| 290 |  | } | 
| 291 |  | #endif //BRAINS_SNAPSHOT_HPP |