Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 2 : Copyright (c) 2015-2020 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/ActionShortcut.h" 23 : #include "core/ActionRegister.h" 24 : 25 : //+PLUMEDOC CONCOMP CLUSTER_DIAMETER 26 : /* 27 : Print out the diameter of one of the connected components. 28 : 29 : An example input that determines the diameter of the second largest cluster that is identified by 30 : analysing the connected components of a [CONTACT_MATRIX](CONTACT_MATRIX.md) using [DFSCLUSTERING](DFSCLUSTERING.md) is shown below: 31 : 32 : ```plumed 33 : # Calculate a contact matrix between the first 100 atoms in the configuration 34 : cm: CONTACT_MATRIX GROUP=1-100 SWITCH={CUBIC D_0=0.45 D_MAX=0.55} 35 : # Find the connected components from the contact matrix 36 : dfs: DFSCLUSTERING ARG=cm 37 : # Returns a 100-dimensional vector that is 1 if the correpsonding atom index is in the largest cluster and is zero otherwise 38 : clust: CLUSTER_WEIGHTS CLUSTERS=dfs CLUSTER=1 39 : # And determine the size of the largest cluster that was identified 40 : c1: CLUSTER_DIAMETER ARG=clust ATOMS=1-100 41 : PRINT ARG=c1 FILE=colvar 42 : ``` 43 : 44 : The diameter here is defined as the largest of all the distances between the pairs of atom in the cluster 45 : 46 : __The output from this action is NOT differentiable__ 47 : 48 : */ 49 : //+ENDPLUMEDOC 50 : 51 : namespace PLMD { 52 : namespace clusters { 53 : 54 : class ClusterDiameter : public ActionShortcut { 55 : public: 56 : static void registerKeywords(Keywords& keys); 57 : explicit ClusterDiameter(const ActionOptions&); 58 : }; 59 : 60 : PLUMED_REGISTER_ACTION(ClusterDiameter,"CLUSTER_DIAMETER") 61 : 62 4 : void ClusterDiameter::registerKeywords( Keywords& keys ) { 63 4 : ActionShortcut::registerKeywords( keys ); 64 4 : keys.add("optional","ARG","calculate ths radius of the cluster that are in this particular cluster"); 65 4 : keys.add("compulsory","ATOMS","the atoms that were used to calculate the matrix that was clustered"); 66 8 : keys.setValueDescription("scalar","the largest of all the distances between the pairs of atom in the cluster"); 67 4 : keys.needsAction("DISTANCE_MATRIX"); 68 4 : keys.needsAction("OUTER_PRODUCT"); 69 4 : keys.needsAction("CUSTOM"); 70 4 : keys.needsAction("FLATTEN"); 71 4 : keys.needsAction("HIGHEST"); 72 4 : keys.addDOI("https://doi.org/10.1021/acs.jctc.6b01073"); 73 4 : } 74 : 75 2 : ClusterDiameter::ClusterDiameter(const ActionOptions& ao): 76 : Action(ao), 77 2 : ActionShortcut(ao) { 78 : // Read in the argument 79 : std::string arg_str, atdata; 80 2 : parse("ARG",arg_str); 81 2 : parse("ATOMS",atdata); 82 : // Distance matrix 83 4 : readInputLine( getShortcutLabel() + "_dmat: DISTANCE_MATRIX GROUP=" + atdata ); 84 : // Matrix of bonds in cluster 85 4 : readInputLine( getShortcutLabel() + "_bmat: OUTER_PRODUCT FUNC=x*y ARG=" + arg_str + "," + arg_str ); 86 : // Product of matrices 87 4 : readInputLine( getShortcutLabel() + "_dcls: CUSTOM ARG=" + getShortcutLabel() + "_dmat," + getShortcutLabel() + "_bmat FUNC=x*y PERIODIC=NO"); 88 : // Convert matrix to a vector to get highest 89 4 : readInputLine( getShortcutLabel() + "_vdcls: FLATTEN ARG=" + getShortcutLabel() + "_dcls" ); 90 : // And take the highest value 91 4 : readInputLine( getShortcutLabel() + ": HIGHEST ARG=" + getShortcutLabel() + "_vdcls"); 92 2 : } 93 : 94 : } 95 : }