Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2011-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 :
23 : #include "core/ActionRegister.h"
24 : #include "core/ActionAtomistic.h"
25 : #include "core/Atoms.h"
26 : #include "tools/IFile.h"
27 : #include "tools/Tools.h"
28 : #include <string>
29 : #include <vector>
30 : #include <algorithm>
31 :
32 : namespace PLMD {
33 : namespace generic {
34 :
35 : //+PLUMEDOC GENERIC GROUP
36 : /*
37 : Define a group of atoms so that a particular list of atoms can be referenced with a single label in definitions of CVs or virtual atoms.
38 :
39 : Atoms can be listed as comma separated numbers (i.e. `1,2,3,10,45,7,9`) , simple positive ranges
40 : (i.e. `20-40`), ranges with a stride either positive or negative (i.e. `20-40:2` or `80-50:-2`) or as
41 : comma separated combinations of all the former methods (`1,2,4,5,10-20,21-40:2,80-50:-2`).
42 :
43 : Moreover, lists can be imported from ndx files (GROMACS format). Use `NDX_FILE` to set the name of
44 : the index file and `NDX_GROUP` to set the name of the group to be imported (default is first one).
45 :
46 : It is also possible to remove atoms from a list and or sort them using keywords `REMOVE`, `SORT`, and `UNIQUE`.
47 : The flow is the following:
48 : - If `ATOMS` is present, then take the ordered list of atoms from the `ATOMS` keyword as a starting list.
49 : - Alternatively, if `NDX_FILE` is present, use the list obtained from the gromacs group.
50 : - If `REMOVE` is present, then remove the first occurrence of each of these atoms from the list.
51 : If one tries to remove an atom that was not listed plumed adds a notice in the output.
52 : An atom that is present twice in the original list might be removed twice.
53 : - If `SORT` is present, then the resulting list is sorted by increasing serial number.
54 : - If `UNIQUE` is present, then the resulting list is sorted by increasing serial number _and_ duplicate elements are removed.
55 :
56 : Notice that this command just creates a shortcut, and does not imply any real calculation.
57 : So, having a huge group defined does not slow down your calculation in any way.
58 : It is just convenient to better organize input files. Might be used in combination with
59 : the \ref INCLUDE command so as to store long group definitions in a separate file.
60 :
61 :
62 : \par Examples
63 :
64 : This command create a group of atoms containing atoms 1, 4, 7, 11 and 14 (labeled 'o'), and another containing
65 : atoms 2, 3, 5, 6, 8, 9, 12, and 13 (labeled 'h'):
66 : \plumedfile
67 : o: GROUP ATOMS=1,4,7,11,14
68 : h: GROUP ATOMS=2,3,5,6,8,9,12,13
69 : # compute the coordination among the two groups
70 : c: COORDINATION GROUPA=o GROUPB=h R_0=0.3
71 : # same could have been obtained without GROUP, just writing:
72 : # c: COORDINATION GROUPA=1,4,7,11,14 GROUPB=2,3,5,6,8,9,12,13
73 :
74 : # print the coordination on file 'colvar'
75 : PRINT ARG=c FILE=colvar
76 : \endplumedfile
77 :
78 : Groups can be conveniently stored in a separate file.
79 : E.g. one could create a file named `groups.dat` which reads
80 : \plumedfile
81 : #SETTINGS FILENAME=groups.dat
82 : # this is groups.dat
83 : o: GROUP ATOMS=1,4,7,11,14
84 : h: GROUP ATOMS=2,3,5,6,8,9,12,13
85 : \endplumedfile
86 : and then include it in the main 'plumed.dat' file
87 : \plumedfile
88 : INCLUDE FILE=groups.dat
89 : # compute the coordination among the two groups
90 : c: COORDINATION GROUPA=o GROUPB=h R_0=0.3
91 : # print the coordination on file 'colvar'
92 : PRINT ARG=c FILE=colvar
93 : \endplumedfile
94 : The `groups.dat` file could be very long and include lists of thousand atoms without cluttering the main plumed.dat file.
95 :
96 : A GROMACS index file such as the one shown below:
97 :
98 : \auxfile{index.ndx}
99 : [ Protein ]
100 : 1 3 5 7 9
101 : 2 4 6 8 10
102 : [ Group2 ]
103 : 30 31 32 33 34 35 36 37 38 39 40
104 : 5
105 : \endauxfile
106 :
107 : can also be imported by using the GROUP keyword as shown below
108 : \plumedfile
109 : # import group named 'Protein' from file index.ndx
110 : pro: GROUP NDX_FILE=index.ndx NDX_GROUP=Protein
111 : # dump all the atoms of the protein on a trajectory file
112 : DUMPATOMS ATOMS=pro FILE=traj.gro
113 : \endplumedfile
114 :
115 : A list can be edited with `REMOVE`. For instance, if you
116 : are using a water model with three atoms per molecule, you can
117 : easily construct the list of hydrogen atoms in this manner
118 : \plumedfile
119 : # take one atom every three, that is oxygens
120 : ox: GROUP ATOMS=1-90:3
121 : # take the remaining atoms, that is hydrogens
122 : hy: GROUP ATOMS=1-90 REMOVE=ox
123 : DUMPATOMS ATOMS=ox FILE=ox.gro
124 : DUMPATOMS ATOMS=hy FILE=hy.gro
125 : \endplumedfile
126 :
127 :
128 : */
129 : //+ENDPLUMEDOC
130 :
131 : class Group:
132 : public ActionAtomistic {
133 :
134 : public:
135 : explicit Group(const ActionOptions&ao);
136 : ~Group();
137 : static void registerKeywords( Keywords& keys );
138 0 : void calculate() override {}
139 0 : void apply() override {}
140 : };
141 :
142 14059 : PLUMED_REGISTER_ACTION(Group,"GROUP")
143 :
144 137 : Group::Group(const ActionOptions&ao):
145 : Action(ao),
146 137 : ActionAtomistic(ao) {
147 : std::vector<AtomNumber> atoms;
148 274 : parseAtomList("ATOMS",atoms);
149 : std::string ndxfile,ndxgroup;
150 137 : parse("NDX_FILE",ndxfile);
151 274 : parse("NDX_GROUP",ndxgroup);
152 137 : if(ndxfile.length()>0 && atoms.size()>0) {
153 0 : error("either use explicit atom list or import from index file");
154 : }
155 137 : if(ndxfile.length()==0 && ndxgroup.size()>0) {
156 0 : error("NDX_GROUP can be only used is NDX_FILE is also used");
157 : }
158 :
159 137 : if(ndxfile.length()>0) {
160 23 : if(ndxgroup.size()>0) {
161 44 : log<<" importing group '"+ndxgroup+"'";
162 : } else {
163 1 : log<<" importing first group";
164 : }
165 23 : log<<" from index file "<<ndxfile<<"\n";
166 :
167 23 : IFile ifile;
168 23 : ifile.open(ndxfile);
169 : std::string line;
170 : std::string groupname;
171 : bool firstgroup=true;
172 : bool groupfound=false;
173 5785 : while(ifile.getline(line)) {
174 5762 : std::vector<std::string> words=Tools::getWords(line);
175 11654 : if(words.size()>=3 && words[0]=="[" && words[2]=="]") {
176 214 : if(groupname.length()>0) {
177 : firstgroup=false;
178 : }
179 : groupname=words[1];
180 214 : if(groupname==ndxgroup || ndxgroup.length()==0) {
181 : groupfound=true;
182 : }
183 5548 : } else if(groupname==ndxgroup || (firstgroup && ndxgroup.length()==0)) {
184 9375 : for(unsigned i=0; i<words.size(); i++) {
185 : AtomNumber at;
186 8770 : Tools::convert(words[i],at);
187 8770 : atoms.push_back(at);
188 : }
189 : }
190 5762 : }
191 23 : if(!groupfound) {
192 0 : error("group has not been found in index file");
193 : }
194 23 : }
195 :
196 : std::vector<AtomNumber> remove;
197 274 : parseAtomList("REMOVE",remove);
198 137 : if(remove.size()>0) {
199 : std::vector<AtomNumber> notfound;
200 : unsigned k=0;
201 2 : log<<" removing these atoms from the list:";
202 114 : for(unsigned i=0; i<remove.size(); i++) {
203 112 : const auto it = find(atoms.begin(),atoms.end(),remove[i]);
204 112 : if(it!=atoms.end()) {
205 111 : if(k%25==0) {
206 6 : log<<"\n";
207 : }
208 111 : log<<" "<<(*it).serial();
209 111 : k++;
210 : atoms.erase(it);
211 : } else {
212 1 : notfound.push_back(remove[i]);
213 : }
214 : }
215 2 : log<<"\n";
216 2 : if(notfound.size()>0) {
217 1 : log<<" the following atoms were not found:";
218 2 : for(unsigned i=0; i<notfound.size(); i++) {
219 1 : log<<" "<<notfound[i].serial();
220 : }
221 1 : log<<"\n";
222 : }
223 : }
224 :
225 137 : bool sortme=false;
226 137 : parseFlag("SORT",sortme);
227 137 : if(sortme) {
228 1 : log<<" atoms are sorted\n";
229 1 : sort(atoms.begin(),atoms.end());
230 : }
231 137 : bool unique=false;
232 137 : parseFlag("UNIQUE",unique);
233 137 : if(unique) {
234 1 : log<<" sorting atoms and removing duplicates\n";
235 1 : Tools::removeDuplicates(atoms);
236 : }
237 :
238 137 : this->atoms.insertGroup(getLabel(),atoms);
239 137 : log.printf(" list of atoms:");
240 186022 : for(unsigned i=0; i<atoms.size(); i++) {
241 185885 : if(i%25==0) {
242 7517 : log<<"\n";
243 : }
244 185885 : log<<" "<<atoms[i].serial();
245 : }
246 137 : log.printf("\n");
247 137 : }
248 :
249 141 : void Group::registerKeywords( Keywords& keys ) {
250 141 : Action::registerKeywords( keys );
251 141 : ActionAtomistic::registerKeywords( keys );
252 282 : keys.add("atoms", "ATOMS", "the numerical indexes for the set of atoms in the group");
253 282 : keys.add("atoms", "REMOVE","remove these atoms from the list");
254 282 : keys.addFlag("SORT",false,"sort the resulting list");
255 282 : keys.addFlag("UNIQUE",false,"sort atoms and remove duplicated ones");
256 282 : keys.add("optional", "NDX_FILE", "the name of index file (gromacs syntax)");
257 282 : keys.add("optional", "NDX_GROUP", "the name of the group to be imported (gromacs syntax) - first group found is used by default");
258 141 : }
259 :
260 274 : Group::~Group() {
261 137 : atoms.removeGroup(getLabel());
262 274 : }
263 :
264 : }
265 : }
266 :
|
