Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2011-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/Angle.h"
25 :
26 : #include <string>
27 : #include <cmath>
28 :
29 : using namespace std;
30 :
31 : namespace PLMD {
32 : namespace colvar {
33 :
34 : //+PLUMEDOC COLVAR ANGLE
35 : /*
36 : Calculate an angle.
37 :
38 : This command can be used to compute the angle between three atoms. Alternatively
39 : if four atoms appear in the atom
40 : specification it calculates the angle between
41 : two vectors identified by two pairs of atoms.
42 :
43 : If _three_ atoms are given, the angle is defined as:
44 : \f[
45 : \theta=\arccos\left(\frac{ {\bf r}_{21}\cdot {\bf r}_{23}}{
46 : |{\bf r}_{21}| |{\bf r}_{23}|}\right)
47 : \f]
48 : Here \f$ {\bf r}_{ij}\f$ is the distance vector among the
49 : i-th and the j-th listed atom.
50 :
51 : If _four_ atoms are given, the angle is defined as:
52 : \f[
53 : \theta=\arccos\left(\frac{ {\bf r}_{21}\cdot {\bf r}_{34}}{
54 : |{\bf r}_{21}| |{\bf r}_{34}|}\right)
55 : \f]
56 :
57 : Notice that angles defined in this way are non-periodic variables and
58 : their value is limited by definition between 0 and \f$\pi\f$.
59 :
60 : The vectors \f$ {\bf r}_{ij}\f$ are by default evaluated taking
61 : periodic boundary conditions into account.
62 : This behavior can be changed with the NOPBC flag.
63 :
64 : \par Examples
65 :
66 : This command tells plumed to calculate the angle between the vector connecting atom 1 to atom 2 and
67 : the vector connecting atom 2 to atom 3 and to print it on file COLVAR1. At the same time,
68 : the angle between vector connecting atom 1 to atom 2 and the vector connecting atom 3 to atom 4 is printed
69 : on file COLVAR2.
70 : \verbatim
71 :
72 : a: ANGLE ATOMS=1,2,3
73 : # equivalently one could state:
74 : # a: ANGLE ATOMS=1,2,2,3
75 :
76 : b: ANGLE ATOMS=1,2,3,4
77 :
78 : PRINT ARG=a FILE=COLVAR1
79 : PRINT ARG=b FILE=COLVAR2
80 : \endverbatim
81 : (see also \ref PRINT)
82 :
83 :
84 : */
85 : //+ENDPLUMEDOC
86 :
87 22 : class Angle : public Colvar {
88 : bool pbc;
89 :
90 : public:
91 : explicit Angle(const ActionOptions&);
92 : // active methods:
93 : virtual void calculate();
94 : static void registerKeywords( Keywords& keys );
95 : };
96 :
97 2534 : PLUMED_REGISTER_ACTION(Angle,"ANGLE")
98 :
99 12 : void Angle::registerKeywords( Keywords& keys ) {
100 12 : Colvar::registerKeywords(keys);
101 12 : keys.add("atoms","ATOMS","the list of atoms involved in this collective variable (either 3 or 4 atoms)");
102 12 : }
103 :
104 11 : Angle::Angle(const ActionOptions&ao):
105 : PLUMED_COLVAR_INIT(ao),
106 11 : pbc(true)
107 : {
108 11 : vector<AtomNumber> atoms;
109 11 : parseAtomList("ATOMS",atoms);
110 11 : bool nopbc=!pbc;
111 11 : parseFlag("NOPBC",nopbc);
112 11 : pbc=!nopbc;
113 :
114 11 : if(atoms.size()==3) {
115 7 : log.printf(" between atoms %d %d %d\n",atoms[0].serial(),atoms[1].serial(),atoms[2].serial());
116 7 : atoms.resize(4);
117 7 : atoms[3]=atoms[2];
118 7 : atoms[2]=atoms[1];
119 4 : } else if(atoms.size()==4) {
120 4 : log.printf(" between lines %d-%d and %d-%d\n",atoms[0].serial(),atoms[1].serial(),atoms[2].serial(),atoms[3].serial());
121 0 : } else error("Number of specified atoms should be either 3 or 4");
122 :
123 11 : if(pbc) log.printf(" using periodic boundary conditions\n");
124 0 : else log.printf(" without periodic boundary conditions\n");
125 :
126 11 : addValueWithDerivatives(); setNotPeriodic();
127 11 : requestAtoms(atoms);
128 11 : checkRead();
129 11 : }
130 :
131 : // calculator
132 219 : void Angle::calculate() {
133 :
134 219 : if(pbc) makeWhole();
135 :
136 219 : Vector dij,dik;
137 219 : dij=delta(getPosition(2),getPosition(3));
138 219 : dik=delta(getPosition(1),getPosition(0));
139 219 : Vector ddij,ddik;
140 : PLMD::Angle a;
141 219 : double angle=a.compute(dij,dik,ddij,ddik);
142 219 : setAtomsDerivatives(0,ddik);
143 219 : setAtomsDerivatives(1,-ddik);
144 219 : setAtomsDerivatives(2,-ddij);
145 219 : setAtomsDerivatives(3,ddij);
146 219 : setValue (angle);
147 219 : setBoxDerivativesNoPbc();
148 219 : }
149 :
150 : }
151 2523 : }
152 :
153 :
154 :
|