Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2016-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 "core/ActionRegister.h"
23 : #include "core/ActionShortcut.h"
24 : #include "core/PlumedMain.h"
25 : #include "core/ActionSet.h"
26 : #include "core/ActionWithValue.h"
27 :
28 : //+PLUMEDOC MCOLVAR EUCLIDEAN_DISTANCE
29 : /*
30 : Calculate the euclidean distance between two vectors of arguments
31 :
32 : If we have two $n$-dimensional vectors $u$ and $v$ we can calculate the
33 : [Euclidean distance](https://en.wikipedia.org/wiki/Euclidean_distance) between the two points as
34 :
35 : $$
36 : d = \sqrt{ \sum_{i=1}^n (u_i - v_i)^2 }
37 : $$
38 :
39 : which can be expressed in matrix form as:
40 :
41 : $$
42 : d^2 = (u-v)^T (u-v)
43 : $$
44 :
45 : The inputs below shows an example where this is used to calculate the Euclidean distance
46 : between the instaneous values of some torsional angles and some reference values
47 : for these torsion. In this first example the input values are vectors:
48 :
49 : ```plumed
50 : c: CONSTANT VALUES=1,2,3
51 : d: DISTANCE ATOMS1=1,2 ATOMS2=3,4 ATOMS3=5,6
52 : dd: EUCLIDEAN_DISTANCE ARG1=c ARG2=d
53 : PRINT ARG=dd FILE=colvar
54 : ```
55 :
56 : while this second example does the same thing but uses scalars in input.
57 :
58 : ```plumed
59 : c1: CONSTANT VALUE=1
60 : d1: DISTANCE ATOMS=1,2
61 : c2: CONSTANT VALUE=2
62 : d2: DISTANCE ATOMS=3,4
63 : c3: CONSTANT VALUE=3
64 : d3: DISTANCE ATOMS=5,6
65 : dd: EUCLIDEAN_DISTANCE ARG1=c1,c2,c3 ARG2=d1,d2,d3
66 : PRINT ARG=dd FILE=colvar
67 : ```
68 :
69 : ## Calculating multiple distances
70 :
71 : Suppose that we now have $m$ reference configurations we can define the following $m$ distances
72 : from these reference configurations:
73 :
74 : $$
75 : d_j^2 = (u-v_j)^T (u-v_j)
76 : $$
77 :
78 : Lets suppose that we put the $m$, $n$-dimensional $(u-v_j)$ vectors in this expression into a
79 : $n\times m$ matrix, $A$, by using the [DISPLACEMENT](DISPLACEMENT.md) command. It is then
80 : straightforward to show that the $d_j^2$ values in the above expression are then the diagonal
81 : elements of the matrix product $A^T A$.
82 :
83 : We can use this idea to calculate multiple EUCLIDEAN_DISTANCE values in the following inputs.
84 : This first example calculates the three distances between the instaneoues values of two torsions
85 : and three reference configurations.
86 :
87 : ```plumed
88 : ref_psi: CONSTANT VALUES=2.25,1.3,-1.5
89 : ref_phi: CONSTANT VALUES=-1.91,-0.6,2.4
90 :
91 : psi: TORSION ATOMS=1,2,3,4
92 : phi: TORSION ATOMS=13,14,15,16
93 :
94 : dd: EUCLIDEAN_DISTANCE ARG1=psi,phi ARG2=ref_psi,ref_phi
95 : PRINT ARG=dd FILE=colvar
96 : ```
97 :
98 : This section example calculates the three distances between a single reference value for the two
99 : torsions and three instances of this pair of torsions.
100 :
101 : ```plumed
102 : ref_psi: CONSTANT VALUES=2.25
103 : ref_phi: CONSTANT VALUES=-1.91
104 :
105 : psi: TORSION ATOMS1=1,2,3,4 ATOMS2=5,6,7,8 ATOMS3=9,10,11,12
106 : phi: TORSION ATOMS1=13,14,15,16 ATOMS2=17,18,19,20 ATOMS3=21,22,23,24
107 :
108 : dd: EUCLIDEAN_DISTANCE ARG1=psi,phi ARG2=ref_psi,ref_phi
109 : PRINT ARG=dd FILE=colvar
110 : ```
111 :
112 : This final example then computes three distances between three pairs of torsional angles and threee
113 : reference values for these three values.
114 :
115 : ```plumed
116 : ref_psi: CONSTANT VALUES=2.25,1.3,-1.5
117 : ref_phi: CONSTANT VALUES=-1.91,-0.6,2.4
118 :
119 : psi: TORSION ATOMS1=1,2,3,4 ATOMS2=5,6,7,8 ATOMS3=9,10,11,12
120 : phi: TORSION ATOMS1=13,14,15,16 ATOMS2=17,18,19,20 ATOMS3=21,22,23,24
121 :
122 : dd: EUCLIDEAN_DISTANCE ARG1=psi,phi ARG2=ref_psi,ref_phi
123 : PRINT ARG=dd FILE=colvar
124 : ```
125 :
126 : */
127 : //+ENDPLUMEDOC
128 :
129 : namespace PLMD {
130 : namespace refdist {
131 :
132 : class EuclideanDistance : public ActionShortcut {
133 : public:
134 : static void registerKeywords( Keywords& keys );
135 : explicit EuclideanDistance(const ActionOptions&ao);
136 : };
137 :
138 : PLUMED_REGISTER_ACTION(EuclideanDistance,"EUCLIDEAN_DISTANCE")
139 :
140 35 : void EuclideanDistance::registerKeywords( Keywords& keys ) {
141 35 : ActionShortcut::registerKeywords(keys);
142 35 : keys.add("compulsory","ARG1","The poin that we are calculating the distance from");
143 35 : keys.add("compulsory","ARG2","The point that we are calculating the distance to");
144 35 : keys.addFlag("SQUARED",false,"The squared distance should be calculated");
145 70 : keys.setValueDescription("scalar/vector","the euclidean distances between the input vectors");
146 35 : keys.needsAction("DISPLACEMENT");
147 35 : keys.needsAction("CUSTOM");
148 35 : keys.needsAction("TRANSPOSE");
149 35 : keys.needsAction("MATRIX_PRODUCT_DIAGONAL");
150 35 : }
151 :
152 28 : EuclideanDistance::EuclideanDistance( const ActionOptions& ao):
153 : Action(ao),
154 28 : ActionShortcut(ao) {
155 : std::string arg1, arg2;
156 28 : parse("ARG1",arg1);
157 28 : parse("ARG2",arg2);
158 : // Vectors are in rows here
159 56 : readInputLine( getShortcutLabel() + "_diff: DISPLACEMENT ARG1=" + arg1 + " ARG2=" + arg2 );
160 : // Get the action that computes the differences
161 28 : ActionWithValue* av = plumed.getActionSet().selectWithLabel<ActionWithValue*>( getShortcutLabel() + "_diff");
162 28 : plumed_assert( av );
163 : // Check if squared
164 : bool squared;
165 28 : parseFlag("SQUARED",squared);
166 28 : std::string olab = getShortcutLabel();
167 28 : if( !squared ) {
168 : olab += "_2";
169 : }
170 : // Deal with an annoying corner case when displacement has a single argument
171 28 : if( av->copyOutput(0)->getRank()==0 ) {
172 0 : readInputLine( olab + ": CUSTOM ARG=" + getShortcutLabel() + "_diff FUNC=x*x PERIODIC=NO");
173 : } else {
174 : // Notice that the vectors are in the columns here
175 56 : readInputLine( getShortcutLabel() + "_diffT: TRANSPOSE ARG=" + getShortcutLabel() + "_diff");
176 56 : readInputLine( olab + ": MATRIX_PRODUCT_DIAGONAL ARG=" + getShortcutLabel() + "_diff," + getShortcutLabel() + "_diffT");
177 : }
178 28 : if( !squared ) {
179 46 : readInputLine( getShortcutLabel() + ": CUSTOM ARG=" + getShortcutLabel() + "_2 FUNC=sqrt(x) PERIODIC=NO");
180 : }
181 28 : }
182 :
183 : }
184 : }
|
