1 /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 Copyright (c) 2014-2021 The plumed team
3 (see the PEOPLE file at the root of the distribution for a list of names)
4
5 See http://www.plumed.org for more information.
6
7 This file is part of plumed, version 2.
8
9 plumed is free software: you can redistribute it and/or modify
10 it under the terms of the GNU Lesser General Public License as published by
11 the Free Software Foundation, either version 3 of the License, or
12 (at your option) any later version.
13
14 plumed is distributed in the hope that it will be useful,
15 but WITHOUT ANY WARRANTY; without even the implied warranty of
16 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17 GNU Lesser General Public License for more details.
18
19 You should have received a copy of the GNU Lesser General Public License
20 along with plumed. If not, see <http://www.gnu.org/licenses/>.
21 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ */
22 #ifndef __PLUMED_multicolvar_VolumeGradientBase_h
23 #define __PLUMED_multicolvar_VolumeGradientBase_h
24
25 #include "BridgedMultiColvarFunction.h"
26
27 namespace PLMD {
28 namespace multicolvar {
29
30 class VolumeGradientBase : public BridgedMultiColvarFunction {
31 friend class MultiColvarBase;
32 private:
33 /// This is used to store forces temporarily in apply
34 std::vector<double> tmpforces;
35 protected:
36 /// Get the cell box
37 const Tensor & getBox() const;
38 /// Get reference to Pbc
39 const Pbc & getPbc() const;
40 /// Calculate distance between two points
41 Vector pbcDistance( const Vector& v1, const Vector& v2) const;
42 /// Get position of atom
43 Vector getPosition( int iatom ) const ;
44 /// Request the atoms
45 void requestAtoms( const std::vector<AtomNumber>& atoms );
46 /// Set the number in the volume
47 void setNumberInVolume( const unsigned&, const unsigned&, const double&, const Vector&, const Tensor&, const std::vector<Vector>&, MultiValue& ) const ;
48 public:
49 static void registerKeywords( Keywords& keys );
50 explicit VolumeGradientBase(const ActionOptions&);
51 /// Do jobs required before tasks are undertaken
52 void doJobsRequiredBeforeTaskList() override;
53 /// Actually do what we are asked
54 void completeTask( const unsigned& curr, MultiValue& invals, MultiValue& outvals ) const override;
55 /// Calculate what is in the volumes
56 virtual void calculateAllVolumes( const unsigned& curr, MultiValue& outvals ) const=0;
57 /// Setup the regions that this is based on
58 virtual void setupRegions()=0;
59 /// Forces here are applied through the bridge
60 void addBridgeForces( const std::vector<double>& bb );
61 };
62
63 inline
getBox()64 const Tensor & VolumeGradientBase::getBox()const {
65 return getPntrToMultiColvar()->getBox();
66 }
67
68 inline
getPbc()69 const Pbc & VolumeGradientBase::getPbc() const {
70 return getPntrToMultiColvar()->getPbc();
71 }
72
73 inline
pbcDistance(const Vector & v1,const Vector & v2)74 Vector VolumeGradientBase::pbcDistance( const Vector& v1, const Vector& v2) const {
75 return getPntrToMultiColvar()->pbcDistance(v1,v2);
76 }
77
78 inline
getPosition(int iatom)79 Vector VolumeGradientBase::getPosition( int iatom ) const {
80 if( !checkNumericalDerivatives() ) return ActionAtomistic::getPosition(iatom);
81 // This is for numerical derivatives of quantity wrt to the local atoms
82 Vector tmp_p = ActionAtomistic::getPosition(iatom);
83 if( bridgeVariable<3*getNumberOfAtoms() ) {
84 if( static_cast<int>(bridgeVariable)>=3*iatom && static_cast<int>(bridgeVariable)<(iatom+1)*3 ) tmp_p[bridgeVariable%3]+=sqrt(epsilon);
85 }
86 // This makes sure that numerical derivatives of virial are calculated correctly
87 tmp_p = ActionAtomistic::getPbc().realToScaled( tmp_p );
88 tmp_p = getPbc().scaledToReal( tmp_p );
89 return tmp_p;
90 }
91
92 }
93 }
94 #endif
95