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