Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2015-2019 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 "ClusterAnalysisBase.h"
23 : #include "core/ActionRegister.h"
24 :
25 : //+PLUMEDOC CONCOMP CLUSTER_NATOMS
26 : /*
27 : Gives the number of atoms in the connected component
28 :
29 : As discussed in the section of the manual on \ref contactmatrix a useful tool for developing complex collective variables is the notion of the
30 : so called adjacency matrix. An adjacency matrix is an \f$N \times N\f$ matrix in which the \f$i\f$th, \f$j\f$th element tells you whether
31 : or not the \f$i\f$th and \f$j\f$th atoms/molecules from a set of \f$N\f$ atoms/molecules are adjacent or not. When analysing these matrix
32 : we can treat them as a graph and find connected components using some clustering algorithm. This action is used in tandem with this form of analysis
33 : to output the number of atoms that are connected together in a particular connected component. It is important to note that the quantity that is
34 : output by this action is not differentiable. As such it cannot be used as a collective variable in a biased simulation.
35 :
36 : \par Examples
37 :
38 : The following input uses PLUMED to calculate a adjacency matrix that connects a pair of atoms if they both have a coordination number that is greater
39 : than 2.0 and if they are within 6.0 nm of each other. Depth first search clustering is used to find the connected components in this matrix and then
40 : the number of atoms in the largest cluster is found. This quantity is then output to a file called colvar
41 :
42 : \plumedfile
43 : # Calculate coordination numbers
44 : c1: COORDINATIONNUMBER SPECIES=1-512 SWITCH={EXP D_0=4.0 R_0=0.5 D_MAX=6.0}
45 : # Select coordination numbers that are more than 2.0
46 : cf: MFILTER_MORE DATA=c1 SWITCH={RATIONAL D_0=2.0 R_0=0.1} LOWMEM
47 : # Build a contact matrix
48 : mat: CONTACT_MATRIX ATOMS=cf SWITCH={EXP D_0=4.0 R_0=0.5 D_MAX=6.0}
49 : # Find largest cluster
50 : dfs: DFSCLUSTERING MATRIX=mat
51 : clust1: CLUSTER_PROPERTIES CLUSTERS=dfs CLUSTER=1
52 : nat: CLUSTER_NATOMS CLUSTERS=dfs CLUSTER=1
53 : PRINT ARG=nat FILE=COLVAR
54 : \endplumedfile
55 :
56 : */
57 : //+ENDPLUMEDOC
58 :
59 : namespace PLMD {
60 : namespace adjmat {
61 :
62 48 : class ClusterSize : public ClusterAnalysisBase {
63 : private:
64 : /// The cluster we are looking for
65 : unsigned clustr;
66 : public:
67 : /// Create manual
68 : static void registerKeywords( Keywords& keys );
69 : /// Constructor
70 : explicit ClusterSize(const ActionOptions&);
71 : ///
72 : void calculate();
73 : ///
74 0 : void performTask( const unsigned& task_index, const unsigned& current, MultiValue& myvals ) const { plumed_error(); }
75 : ///
76 : void turnOnDerivatives();
77 : };
78 :
79 6476 : PLUMED_REGISTER_ACTION(ClusterSize,"CLUSTER_NATOMS")
80 :
81 25 : void ClusterSize::registerKeywords( Keywords& keys ) {
82 25 : ClusterAnalysisBase::registerKeywords( keys );
83 125 : keys.add("compulsory","CLUSTER","1","which cluster would you like to look at 1 is the largest cluster, 2 is the second largest, 3 is the the third largest and so on.");
84 25 : }
85 :
86 24 : ClusterSize::ClusterSize(const ActionOptions&ao):
87 : Action(ao),
88 24 : ClusterAnalysisBase(ao)
89 : {
90 : // Find out which cluster we want
91 48 : parse("CLUSTER",clustr);
92 :
93 24 : if( clustr<1 ) error("cannot look for a cluster larger than the largest cluster");
94 24 : if( clustr>getNumberOfNodes() ) error("cluster selected is invalid - too few atoms in system");
95 :
96 : // Create all tasks by copying those from underlying DFS object (which is actually MultiColvar)
97 24 : for(unsigned i=0; i<getNumberOfNodes(); ++i) addTaskToList(i);
98 : // And now finish the setup of everything in the base
99 24 : std::vector<AtomNumber> fake_atoms; setupMultiColvarBase( fake_atoms );
100 24 : addValue(); setNotPeriodic();
101 24 : }
102 :
103 0 : void ClusterSize::turnOnDerivatives() {
104 0 : error("cannot calculate derivatives of number of atoms in cluster. This quantity is not differentiable");
105 0 : }
106 :
107 27 : void ClusterSize::calculate() {
108 : // Retrieve the atoms in the largest cluster
109 54 : std::vector<unsigned> myatoms; retrieveAtomsInCluster( clustr, myatoms ); setValue( myatoms.size() );
110 27 : }
111 :
112 : }
113 4839 : }
|
