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