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 10592 : PLUMED_REGISTER_ACTION(Position,"POSITION")
85 :
86 88 : void Position::registerKeywords( Keywords& keys ) {
87 88 : Colvar::registerKeywords( keys );
88 88 : componentsAreNotOptional(keys);
89 176 : keys.add("atoms","ATOM","the atom number");
90 176 : 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 176 : keys.addOutputComponent("x","default","the x-component of the atom position");
92 176 : keys.addOutputComponent("y","default","the y-component of the atom position");
93 176 : keys.addOutputComponent("z","default","the z-component of the atom position");
94 264 : keys.addOutputComponent("a","SCALED_COMPONENTS","the normalized projection on the first lattice vector of the atom position");
95 264 : keys.addOutputComponent("b","SCALED_COMPONENTS","the normalized projection on the second lattice vector of the atom position");
96 264 : keys.addOutputComponent("c","SCALED_COMPONENTS","the normalized projection on the third lattice vector of the atom position");
97 88 : }
98 :
99 87 : Position::Position(const ActionOptions&ao):
100 : PLUMED_COLVAR_INIT(ao),
101 87 : scaled_components(false),
102 87 : pbc(true)
103 : {
104 : std::vector<AtomNumber> atoms;
105 174 : parseAtomList("ATOM",atoms);
106 87 : if(atoms.size()!=1)
107 1 : error("Number of specified atoms should be 1");
108 86 : parseFlag("SCALED_COMPONENTS",scaled_components);
109 86 : bool nopbc=!pbc;
110 86 : parseFlag("NOPBC",nopbc);
111 86 : pbc=!nopbc;
112 86 : checkRead();
113 :
114 86 : log.printf(" for atom %d\n",atoms[0].serial());
115 86 : if(pbc) log.printf(" using periodic boundary conditions\n");
116 4 : else log.printf(" without periodic boundary conditions\n");
117 :
118 86 : if(scaled_components) {
119 12 : addComponentWithDerivatives("a"); componentIsPeriodic("a","-0.5","+0.5");
120 12 : addComponentWithDerivatives("b"); componentIsPeriodic("b","-0.5","+0.5");
121 12 : addComponentWithDerivatives("c"); componentIsPeriodic("c","-0.5","+0.5");
122 : } else {
123 164 : addComponentWithDerivatives("x"); componentIsNotPeriodic("x");
124 164 : addComponentWithDerivatives("y"); componentIsNotPeriodic("y");
125 165 : addComponentWithDerivatives("z"); componentIsNotPeriodic("z");
126 82 : log<<" WARNING: components will not have the proper periodicity - see manual\n";
127 : }
128 :
129 86 : requestAtoms(atoms);
130 88 : }
131 :
132 :
133 : // calculator
134 8037 : void Position::calculate() {
135 :
136 8037 : Vector distance;
137 8037 : if(pbc) {
138 7992 : distance=pbcDistance(Vector(0.0,0.0,0.0),getPosition(0));
139 : } else {
140 45 : distance=delta(Vector(0.0,0.0,0.0),getPosition(0));
141 : }
142 :
143 8037 : if(scaled_components) {
144 41 : Value* valuea=getPntrToComponent("a");
145 41 : Value* valueb=getPntrToComponent("b");
146 82 : Value* valuec=getPntrToComponent("c");
147 41 : Vector d=getPbc().realToScaled(distance);
148 41 : setAtomsDerivatives (valuea,0,matmul(getPbc().getInvBox(),Vector(+1,0,0)));
149 41 : valuea->set(Tools::pbc(d[0]));
150 41 : setAtomsDerivatives (valueb,0,matmul(getPbc().getInvBox(),Vector(0,+1,0)));
151 41 : valueb->set(Tools::pbc(d[1]));
152 41 : setAtomsDerivatives (valuec,0,matmul(getPbc().getInvBox(),Vector(0,0,+1)));
153 41 : valuec->set(Tools::pbc(d[2]));
154 : } else {
155 7996 : Value* valuex=getPntrToComponent("x");
156 7996 : Value* valuey=getPntrToComponent("y");
157 7996 : Value* valuez=getPntrToComponent("z");
158 :
159 7996 : setAtomsDerivatives (valuex,0,Vector(+1,0,0));
160 7996 : setBoxDerivatives (valuex,Tensor(distance,Vector(-1,0,0)));
161 7996 : valuex->set(distance[0]);
162 :
163 7996 : setAtomsDerivatives (valuey,0,Vector(0,+1,0));
164 7996 : setBoxDerivatives (valuey,Tensor(distance,Vector(0,-1,0)));
165 7996 : valuey->set(distance[1]);
166 :
167 7996 : setAtomsDerivatives (valuez,0,Vector(0,0,+1));
168 7996 : setBoxDerivatives (valuez,Tensor(distance,Vector(0,0,-1)));
169 7996 : valuez->set(distance[2]);
170 : }
171 8037 : }
172 :
173 : }
174 : }
175 :
176 :
177 :
|
