Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2014-2018 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 : #include "Colvar.h"
23 : #include "ActionRegister.h"
24 : #include "tools/Pbc.h"
25 :
26 : #include <string>
27 : #include <cmath>
28 :
29 : using namespace std;
30 :
31 : namespace PLMD {
32 : namespace colvar {
33 :
34 : //+PLUMEDOC COLVAR POSITION
35 : /*
36 : Calculate the components of the position of an atom.
37 :
38 : Notice that single components will not have the proper periodicity!
39 : If you need the values to be consistent through PBC you should use SCALED_COMPONENTS,
40 : which defines values that by construction are in the -0.5,0.5 domain. This is
41 : similar to the equivalent flag for \ref DISTANCE.
42 : Also notice that by default the minimal image distance from the
43 : origin is considered (can be changed with NOPBC).
44 :
45 : \attention
46 : This variable should be used with extreme care since it allows to easily go into troubles. See comments below.
47 :
48 : This variable can be safely used only if
49 : Hamiltonian is not invariant for translation (i.e. there are other absolute positions which are biased, e.g. by position restraints)
50 : and cell size and shapes are fixed through the simulation.
51 :
52 : If you are not in this situation and still want to use the absolute position of an atom you should first fix the reference frame.
53 : This can be done e.g. using \ref FIT_TO_TEMPLATE.
54 :
55 : \par Examples
56 :
57 : \verbatim
58 : # align to a template
59 : FIT_TO_TEMPLATE REFERENCE=ref.pdb
60 : p: POSITION ATOM=3
61 : PRINT ARG=p.x,p.y,p.z
62 : \endverbatim
63 : (see also \ref FIT_TO_TEMPLATE)
64 :
65 :
66 : */
67 : //+ENDPLUMEDOC
68 :
69 36 : class Position : public Colvar {
70 : bool scaled_components;
71 : bool pbc;
72 :
73 : public:
74 : static void registerKeywords( Keywords& keys );
75 : explicit Position(const ActionOptions&);
76 : // active methods:
77 : virtual void calculate();
78 : };
79 :
80 2541 : PLUMED_REGISTER_ACTION(Position,"POSITION")
81 :
82 19 : void Position::registerKeywords( Keywords& keys ) {
83 19 : Colvar::registerKeywords( keys );
84 19 : componentsAreNotOptional(keys);
85 19 : keys.add("atoms","ATOM","the atom number");
86 19 : keys.addFlag("SCALED_COMPONENTS",false,"calculate the a, b and c scaled components of the position separately and store them as label.a, label.b and label.c");
87 19 : keys.addOutputComponent("x","default","the x-component of the atom position");
88 19 : keys.addOutputComponent("y","default","the y-component of the atom position");
89 19 : keys.addOutputComponent("z","default","the z-component of the atom position");
90 19 : keys.addOutputComponent("a","SCALED_COMPONENTS","the normalized projection on the first lattice vector of the atom position");
91 19 : keys.addOutputComponent("b","SCALED_COMPONENTS","the normalized projection on the second lattice vector of the atom position");
92 19 : keys.addOutputComponent("c","SCALED_COMPONENTS","the normalized projection on the third lattice vector of the atom position");
93 19 : }
94 :
95 18 : Position::Position(const ActionOptions&ao):
96 : PLUMED_COLVAR_INIT(ao),
97 : scaled_components(false),
98 18 : pbc(true)
99 : {
100 18 : vector<AtomNumber> atoms;
101 18 : parseAtomList("ATOM",atoms);
102 18 : if(atoms.size()!=1)
103 0 : error("Number of specified atoms should be 1");
104 18 : parseFlag("SCALED_COMPONENTS",scaled_components);
105 18 : bool nopbc=!pbc;
106 18 : parseFlag("NOPBC",nopbc);
107 18 : pbc=!nopbc;
108 18 : checkRead();
109 :
110 18 : log.printf(" for atom %d\n",atoms[0].serial());
111 18 : if(pbc) log.printf(" using periodic boundary conditions\n");
112 3 : else log.printf(" without periodic boundary conditions\n");
113 :
114 18 : if(scaled_components) {
115 4 : addComponentWithDerivatives("a"); componentIsPeriodic("a","-0.5","+0.5");
116 4 : addComponentWithDerivatives("b"); componentIsPeriodic("b","-0.5","+0.5");
117 4 : addComponentWithDerivatives("c"); componentIsPeriodic("c","-0.5","+0.5");
118 : } else {
119 14 : addComponentWithDerivatives("x"); componentIsNotPeriodic("x");
120 14 : addComponentWithDerivatives("y"); componentIsNotPeriodic("y");
121 14 : addComponentWithDerivatives("z"); componentIsNotPeriodic("z");
122 14 : log<<" WARNING: components will not have the proper periodicity - see manual\n";
123 : }
124 :
125 18 : requestAtoms(atoms);
126 18 : }
127 :
128 :
129 : // calculator
130 174 : void Position::calculate() {
131 :
132 174 : Vector distance;
133 174 : if(pbc) {
134 159 : distance=pbcDistance(Vector(0.0,0.0,0.0),getPosition(0));
135 : } else {
136 15 : distance=delta(Vector(0.0,0.0,0.0),getPosition(0));
137 : }
138 :
139 174 : if(scaled_components) {
140 41 : Value* valuea=getPntrToComponent("a");
141 41 : Value* valueb=getPntrToComponent("b");
142 41 : Value* valuec=getPntrToComponent("c");
143 41 : Vector d=getPbc().realToScaled(distance);
144 41 : setAtomsDerivatives (valuea,0,matmul(getPbc().getInvBox(),Vector(+1,0,0)));
145 41 : valuea->set(Tools::pbc(d[0]));
146 41 : setAtomsDerivatives (valueb,0,matmul(getPbc().getInvBox(),Vector(0,+1,0)));
147 41 : valueb->set(Tools::pbc(d[1]));
148 41 : setAtomsDerivatives (valuec,0,matmul(getPbc().getInvBox(),Vector(0,0,+1)));
149 41 : valuec->set(Tools::pbc(d[2]));
150 : } else {
151 133 : Value* valuex=getPntrToComponent("x");
152 133 : Value* valuey=getPntrToComponent("y");
153 133 : Value* valuez=getPntrToComponent("z");
154 :
155 133 : setAtomsDerivatives (valuex,0,Vector(+1,0,0));
156 133 : setBoxDerivatives (valuex,Tensor(distance,Vector(-1,0,0)));
157 133 : valuex->set(distance[0]);
158 :
159 133 : setAtomsDerivatives (valuey,0,Vector(0,+1,0));
160 133 : setBoxDerivatives (valuey,Tensor(distance,Vector(0,-1,0)));
161 133 : valuey->set(distance[1]);
162 :
163 133 : setAtomsDerivatives (valuez,0,Vector(0,0,+1));
164 133 : setBoxDerivatives (valuez,Tensor(distance,Vector(0,0,-1)));
165 133 : valuez->set(distance[2]);
166 : }
167 174 : }
168 :
169 : }
170 2523 : }
171 :
172 :
173 :
|
