Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2011-2017 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 "MultiColvarShortcuts.h"
23 : #include "core/ActionShortcut.h"
24 : #include "core/ActionRegister.h"
25 : #include "tools/Pbc.h"
26 :
27 : #include <string>
28 : #include <cmath>
29 :
30 : //+PLUMEDOC MCOLVAR PLANES
31 : /*
32 : Calculate the components of the normals to the planes containing groups of three atoms.
33 :
34 : An example input using this shortcut is shown below:
35 :
36 : ```plumed
37 : v3: PLANES ...
38 : ATOMS1=9,10,11
39 : ATOMS2=89,90,91
40 : ATOMS3=473,474,475
41 : ATOMS4=1161,1162,1163
42 : ATOMS5=1521,1522,1523
43 : ATOMS6=1593,1594,1595
44 : ATOMS7=1601,1602,1603
45 : ATOMS8=2201,2202,2203
46 : VMEAN
47 : ...
48 : PRINT ARG=v3_vmean FILE=colvar
49 : ```
50 :
51 : As you can see if you expand the shortcut above, the input calculates the vectors perpendicular to the planes containing
52 : each group of three atoms. The mean value is then computed by adding all these three dimensional vectors together. The final
53 : value output is the norm of this mean vector. This norm is large if the orientations of the planes containing the groups of
54 : atoms are the same and small if they are not.
55 :
56 : */
57 : //+ENDPLUMEDOC
58 :
59 : namespace PLMD {
60 : namespace multicolvar {
61 :
62 : class PlaneShortcut : public ActionShortcut {
63 : private:
64 : void createVectorNormInput( const std::string& ilab, const std::string& olab, const std::string& vlab );
65 : public:
66 : static void registerKeywords( Keywords& keys );
67 : PlaneShortcut(const ActionOptions&);
68 : };
69 :
70 : PLUMED_REGISTER_ACTION(PlaneShortcut,"PLANES")
71 :
72 3 : void PlaneShortcut::registerKeywords( Keywords& keys ) {
73 3 : ActionShortcut::registerKeywords( keys );
74 3 : MultiColvarShortcuts::shortcutKeywords( keys );
75 3 : keys.add("numbered","ATOMS","the sets of atoms that you would like to calculate the planes for");
76 3 : keys.add("numbered","LOCATION","the location at which the CV is assumed to be in space");
77 6 : keys.reset_style("LOCATION","atoms");
78 3 : keys.addFlag("VMEAN",false,"calculate the norm of the mean vector.");
79 6 : keys.addOutputComponent("_vmean","VMEAN","scalar","the norm of the mean vector");
80 3 : keys.addFlag("VSUM",false,"calculate the norm of the sum of all the vectors");
81 6 : keys.addOutputComponent("_vsum","VSUM","scalar","the norm of the mean vector");
82 3 : keys.needsAction("CENTER");
83 3 : keys.needsAction("GROUP");
84 3 : keys.needsAction("PLANE");
85 3 : keys.needsAction("MEAN");
86 3 : keys.needsAction("SUM");
87 3 : keys.needsAction("COMBINE");
88 3 : keys.needsAction("CUSTOM");
89 3 : }
90 :
91 1 : PlaneShortcut::PlaneShortcut(const ActionOptions&ao):
92 : Action(ao),
93 1 : ActionShortcut(ao) {
94 : bool vmean, vsum;
95 1 : parseFlag("VMEAN",vmean);
96 2 : parseFlag("VSUM",vsum);
97 : std::string dline;
98 1 : std::string grpstr = getShortcutLabel() + "_grp: GROUP ATOMS=";
99 1 : for(unsigned i=1;; ++i) {
100 : std::string atstring;
101 18 : parseNumbered("ATOMS",i,atstring);
102 9 : if( atstring.length()==0 ) {
103 : break;
104 : }
105 : std::string locstr;
106 16 : parseNumbered("LOCATION",i,locstr);
107 8 : if( locstr.length()==0 ) {
108 : std::string num;
109 8 : Tools::convert( i, num );
110 16 : readInputLine( getShortcutLabel() + "_vatom" + num + ": CENTER ATOMS=" + atstring );
111 8 : if( i==1 ) {
112 2 : grpstr += getShortcutLabel() + "_vatom" + num;
113 : } else {
114 14 : grpstr += "," + getShortcutLabel() + "_vatom" + num;
115 : }
116 : } else {
117 0 : if( i==1 ) {
118 : grpstr += locstr;
119 : } else {
120 0 : grpstr += "," + locstr;
121 : }
122 : }
123 : std::string num;
124 8 : Tools::convert( i, num );
125 16 : dline += " ATOMS" + num + "=" + atstring;
126 8 : }
127 1 : readInputLine( grpstr );
128 2 : readInputLine( getShortcutLabel() + ": PLANE " + dline + " " + convertInputLineToString() );
129 1 : if( vmean ) {
130 2 : readInputLine( getShortcutLabel() + "_xs: MEAN ARG=" + getShortcutLabel() + ".x PERIODIC=NO");
131 2 : readInputLine( getShortcutLabel() + "_ys: MEAN ARG=" + getShortcutLabel() + ".y PERIODIC=NO");
132 2 : readInputLine( getShortcutLabel() + "_zs: MEAN ARG=" + getShortcutLabel() + ".z PERIODIC=NO");
133 : // Now calculate the total length of the vector
134 2 : createVectorNormInput( getShortcutLabel(), getShortcutLabel() + "_vmean", "s" );
135 : }
136 1 : if( vsum ) {
137 0 : readInputLine( getShortcutLabel() + "_xz: SUM ARG=" + getShortcutLabel() + ".x PERIODIC=NO");
138 0 : readInputLine( getShortcutLabel() + "_yz: SUM ARG=" + getShortcutLabel() + ".y PERIODIC=NO");
139 0 : readInputLine( getShortcutLabel() + "_zz: SUM ARG=" + getShortcutLabel() + ".z PERIODIC=NO");
140 : // Now calculate the total length of the vector
141 0 : createVectorNormInput( getShortcutLabel(), getShortcutLabel() + "_vsum", "z" );
142 : }
143 1 : }
144 :
145 1 : void PlaneShortcut::createVectorNormInput( const std::string& ilab, const std::string& olab, const std::string& vlab ) {
146 2 : readInputLine( olab + "2: COMBINE ARG=" + ilab + "_x" + vlab + "," + ilab + "_y" + vlab + "," + ilab + "_z" + vlab + " POWERS=2,2,2 PERIODIC=NO");
147 2 : readInputLine( olab + ": CUSTOM ARG=" + olab + "2 FUNC=sqrt(x) PERIODIC=NO");
148 1 : }
149 :
150 : }
151 : }
152 :
153 :
154 :
|
