Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2014-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 "AdjacencyMatrixVessel.h"
24 : #include "core/ActionRegister.h"
25 :
26 : //+PLUMEDOC CONCOMP CLUSTER_PROPERTIES
27 : /*
28 : Calculate properties of the distribution of some quantities that are part of a connected component
29 :
30 : This collective variable was developed for looking at nucleation phenomena, where you are
31 : interested in using studying the behavior of atoms in small aggregates or nuclei. In these sorts of
32 : problems you might be interested in the degree the atoms in a nucleus have adopted their crystalline
33 : structure or (in the case of heterogenous nucleation of a solute from a solvent) you might be
34 : interested in how many atoms are present in the largest cluster \cite tribello-clustering.
35 :
36 : \par Examples
37 :
38 : The input below calculates the coordination numbers of atoms 1-100 and then computes the an adjacency
39 : matrix whose elements measures whether atoms \f$i\f$ and \f$j\f$ are within 0.55 nm of each other. The action
40 : labelled dfs then treats the elements of this matrix as zero or ones and thus thinks of the matrix as defining
41 : a graph. This dfs action then finds the largest connected component in this graph. The sum of the coordination
42 : numbers for the atoms in this largest connected component are then computed and this quantity is output to a colvar
43 : file. The way this input can be used is described in detail in \cite tribello-clustering.
44 :
45 : \plumedfile
46 : lq: COORDINATIONNUMBER SPECIES=1-100 SWITCH={CUBIC D_0=0.45 D_MAX=0.55} LOWMEM
47 : cm: CONTACT_MATRIX ATOMS=lq SWITCH={CUBIC D_0=0.45 D_MAX=0.55}
48 : dfs: DFSCLUSTERING MATRIX=cm
49 : clust1: CLUSTER_PROPERTIES CLUSTERS=dfs CLUSTER=1 SUM
50 : PRINT ARG=clust1.* FILE=colvar
51 : \endplumedfile
52 :
53 : */
54 : //+ENDPLUMEDOC
55 :
56 : namespace PLMD {
57 : namespace adjmat {
58 :
59 54 : class ClusterProperties : public ClusterAnalysisBase {
60 : private:
61 : /// The cluster we are looking for
62 : unsigned clustr;
63 : public:
64 : /// Create manual
65 : static void registerKeywords( Keywords& keys );
66 : /// Constructor
67 : explicit ClusterProperties(const ActionOptions&);
68 : /// Do the calculation
69 : void calculate();
70 : /// We can use ActionWithVessel to run all the calculation
71 : void performTask( const unsigned&, const unsigned&, MultiValue& ) const ;
72 : };
73 :
74 6479 : PLUMED_REGISTER_ACTION(ClusterProperties,"CLUSTER_PROPERTIES")
75 :
76 28 : void ClusterProperties::registerKeywords( Keywords& keys ) {
77 28 : ClusterAnalysisBase::registerKeywords( keys );
78 140 : 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.");
79 112 : keys.use("MEAN"); keys.use("MORE_THAN"); keys.use("LESS_THAN");
80 84 : if( keys.reserved("VMEAN") ) keys.use("VMEAN");
81 84 : if( keys.reserved("VSUM") ) keys.use("VSUM");
82 140 : keys.use("BETWEEN"); keys.use("HISTOGRAM"); keys.use("MOMENTS"); keys.use("ALT_MIN");
83 168 : keys.use("MIN"); keys.use("MAX"); keys.use("SUM"); keys.use("LOWEST"); keys.use("HIGHEST");
84 28 : }
85 :
86 27 : ClusterProperties::ClusterProperties(const ActionOptions&ao):
87 : Action(ao),
88 27 : ClusterAnalysisBase(ao)
89 : {
90 : // Find out which cluster we want
91 54 : parse("CLUSTER",clustr);
92 :
93 27 : if( clustr<1 ) error("cannot look for a cluster larger than the largest cluster");
94 27 : 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 27 : for(unsigned i=0; i<getNumberOfNodes(); ++i) addTaskToList(i);
98 :
99 : // And now finish the setup of everything in the base
100 27 : std::vector<AtomNumber> fake_atoms; setupMultiColvarBase( fake_atoms );
101 27 : }
102 :
103 34 : void ClusterProperties::calculate() {
104 : // Retrieve the atoms in the largest cluster
105 34 : std::vector<unsigned> myatoms; retrieveAtomsInCluster( clustr, myatoms );
106 : // Activate the relevant tasks
107 34 : deactivateAllTasks();
108 5424 : for(unsigned i=0; i<myatoms.size(); ++i) taskFlags[myatoms[i]]=1;
109 34 : lockContributors();
110 : // Now do the calculation
111 34 : runAllTasks();
112 34 : }
113 :
114 1795 : void ClusterProperties::performTask( const unsigned& task_index, const unsigned& current, MultiValue& myvals ) const {
115 1795 : std::vector<double> vals( myvals.getNumberOfValues() ); getPropertiesOfNode( current, vals );
116 1795 : if( !doNotCalculateDerivatives() ) getNodePropertyDerivatives( current, myvals );
117 14360 : for(unsigned k=0; k<vals.size(); ++k) myvals.setValue( k, vals[k] );
118 1795 : }
119 :
120 : }
121 4839 : }
|