Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2013-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 "Steinhardt.h"
23 : #include "LocalSteinhardt.h"
24 : #include "core/ActionRegister.h"
25 :
26 : //+PLUMEDOC MCOLVAR Q6
27 : /*
28 : Calculate 6th order Steinhardt parameters.
29 :
30 : The 6th order Steinhardt parameters allow us to measure the degree to which the first coordination shell
31 : around an atom is ordered. The Steinhardt parameter for atom, \f$i\f$ is complex vector whose components are
32 : calculated using the following formula:
33 :
34 : \f[
35 : q_{6m}(i) = \frac{\sum_j \sigma( r_{ij} ) Y_{6m}(\mathbf{r}_{ij}) }{\sum_j \sigma( r_{ij} ) }
36 : \f]
37 :
38 : where \f$Y_{6m}\f$ is one of the 6th order spherical harmonics so \f$m\f$ is a number that runs from \f$-6\f$ to
39 : \f$+6\f$. The function \f$\sigma( r_{ij} )\f$ is a \ref switchingfunction that acts on the distance between
40 : atoms \f$i\f$ and \f$j\f$. The parameters of this function should be set so that it the function is equal to one
41 : when atom \f$j\f$ is in the first coordination sphere of atom \f$i\f$ and is zero otherwise.
42 :
43 : The Steinhardt parameters can be used to measure the degree of order in the system in a variety of different ways. The
44 : simplest way of measuring whether or not the coordination sphere is ordered is to simply take the norm of the above vector i.e.
45 :
46 : \f[
47 : Q_6(i) = \sqrt{ \sum_{m=-6}^6 q_{6m}(i)^{*} q_{6m}(i) }
48 : \f]
49 :
50 : This norm is small when the coordination shell is disordered and larger when the coordination shell is ordered. Furthermore, when
51 : the keywords LESS_THAN, MIN, MAX, HISTOGRAM, MEAN and so on are used with this colvar it is the distribution of these normed quantities
52 : that is investigated.
53 :
54 : Other measures of order can be taken by averaging the components of the individual \f$q_6\f$ vectors individually or by taking dot products of
55 : the \f$q_{6}\f$ vectors on adjacent atoms. More information on these variables can be found in the documentation for \ref LOCAL_Q6,
56 : \ref LOCAL_AVERAGE and \ref NLINKS.
57 :
58 : \par Examples
59 :
60 : The following command calculates the average Q6 parameter for the 64 atoms in a box of Lennard Jones and prints this
61 : quantity to a file called colvar:
62 :
63 : \plumedfile
64 : Q6 SPECIES=1-64 D_0=1.3 R_0=0.2 MEAN LABEL=q6
65 : PRINT ARG=q6.mean FILE=colvar
66 : \endplumedfile
67 :
68 : The following command calculates the histogram of Q6 parameters for the 64 atoms in a box of Lennard Jones and prints these
69 : quantities to a file called colvar:
70 :
71 : \plumedfile
72 : Q6 SPECIES=1-64 D_0=1.3 R_0=0.2 HISTOGRAM={GAUSSIAN LOWER=0.0 UPPER=1.0 NBINS=20 SMEAR=0.1} LABEL=q6
73 : PRINT ARG=q6.* FILE=colvar
74 : \endplumedfile
75 :
76 : The following command could be used to measure the Q6 paramters that describe the arrangement of chlorine ions around the
77 : sodium atoms in NaCl. The imagined system here is composed of 64 NaCl formula units and the atoms are arranged in the input
78 : with the 64 Na\f$^+\f$ ions followed by the 64 Cl\f$-\f$ ions. Once again the average Q6 paramter is calculated and output to a
79 : file called colvar
80 :
81 : \plumedfile
82 : Q6 SPECIESA=1-64 SPECIESB=65-128 D_0=1.3 R_0=0.2 MEAN LABEL=q6
83 : PRINT ARG=q6.mean FILE=colvar
84 : \endplumedfile
85 :
86 : */
87 : //+ENDPLUMEDOC
88 :
89 : //+PLUMEDOC MCOLVARF LOCAL_Q6
90 : /*
91 : Calculate the local degree of order around an atoms by taking the average dot product between the \f$q_6\f$ vector on the central atom and the \f$q_6\f$ vector
92 : on the atoms in the first coordination sphere.
93 :
94 : The \ref Q6 command allows one to calculate one complex vectors for each of the atoms in your system that describe the degree of order in the coordination sphere
95 : around a particular atom. The difficulty with these vectors comes when combining the order parameters from all of the individual atoms/molecules so as to get a
96 : measure of the global degree of order for the system. The simplest way of doing this - calculating the average Steinhardt parameter - can be problematic. If one is
97 : examining nucleation say only the order parameters for those atoms in the nucleus will change significantly when the nucleus forms. The order parameters for the
98 : atoms in the surrounding liquid will remain pretty much the same. As such if one models a small nucleus embedded in a very large amount of solution/melt any
99 : change in the average order parameter will be negligible. Substantial changes in the value of this average can be observed in simulations of nucleation but only
100 : because the number of atoms is relatively small.
101 :
102 : When the average \ref Q6 parameter is used to bias the dynamics a problems
103 : can occur. These averaged coordinates cannot distinguish between the correct,
104 : single-nucleus pathway and a concerted pathway in which all the atoms rearrange
105 : themselves into their solid-like configuration simultaneously. This second type
106 : of pathway would be impossible in reality because there is a large entropic
107 : barrier that prevents concerted processes like this from happening. However,
108 : in the finite sized systems that are commonly simulated this barrier is reduced
109 : substantially. As a result in simulations where average Steinhardt parameters
110 : are biased there are often quite dramatic system size effects
111 :
112 : If one wants to simulate nucleation using some form on
113 : biased dynamics what is really required is an order parameter that measures:
114 :
115 : - Whether or not the coordination spheres around atoms are ordered
116 : - Whether or not the atoms that are ordered are clustered together in a crystalline nucleus
117 :
118 : \ref LOCAL_AVERAGE and \ref NLINKS are variables that can be combined with the Steinhardt parameteters allow to calculate variables that satisfy these requirements.
119 : LOCAL_Q6 is another variable that can be used in these sorts of calculations. The LOCAL_Q6 parameter for a particular atom is a number that measures the extent to
120 : which the orientation of the atoms in the first coordination sphere of an atom match the orientation of the central atom. It does this by calculating the following
121 : quantity for each of the atoms in the system:
122 :
123 : \f[
124 : s_i = \frac{ \sum_j \sigma( r_{ij} ) \sum_{m=-6}^6 q_{6m}^{*}(i)q_{6m}(j) }{ \sum_j \sigma( r_{ij} ) }
125 : \f]
126 :
127 : where \f$q_{6m}(i)\f$ and \f$q_{6m}(j)\f$ are the 6th order Steinhardt vectors calculated for atom \f$i\f$ and atom \f$j\f$ respectively and the asterix denotes
128 : complex conjugation. The function
129 : \f$\sigma( r_{ij} )\f$ is a \ref switchingfunction that acts on the distance between atoms \f$i\f$ and \f$j\f$. The parameters of this function should be set
130 : so that it the function is equal to one when atom \f$j\f$ is in the first coordination sphere of atom \f$i\f$ and is zero otherwise. The sum in the numerator
131 : of this expression is the dot product of the Steinhardt parameters for atoms \f$i\f$ and \f$j\f$ and thus measures the degree to which the orientations of these
132 : adjacent atoms is correlated.
133 :
134 : \par Examples
135 :
136 : The following command calculates the average value of the LOCAL_Q6 parameter for the 64 Lennard Jones atoms in the system under study and prints this
137 : quantity to a file called colvar.
138 :
139 : \plumedfile
140 : Q6 SPECIES=1-64 D_0=1.3 R_0=0.2 LABEL=q6
141 : LOCAL_Q6 SPECIES=q6 SWITCH={RATIONAL D_0=1.3 R_0=0.2} MEAN LABEL=lq6
142 : PRINT ARG=lq6.mean FILE=colvar
143 : \endplumedfile
144 :
145 : The following input calculates the distribution of LOCAL_Q6 parameters at any given time and outputs this information to a file.
146 :
147 : \plumedfile
148 : Q6 SPECIES=1-64 D_0=1.3 R_0=0.2 LABEL=q6
149 : LOCAL_Q6 SPECIES=q6 SWITCH={RATIONAL D_0=1.3 R_0=0.2} HISTOGRAM={GAUSSIAN LOWER=0.0 UPPER=1.0 NBINS=20 SMEAR=0.1} LABEL=lq6
150 : PRINT ARG=lq6.* FILE=colvar
151 : \endplumedfile
152 :
153 : The following calculates the LOCAL_Q6 parameters for atoms 1-5 only. For each of these atoms comparisons of the geometry of the coordination sphere
154 : are done with those of all the other atoms in the system. The final quantity is the average and is outputted to a file
155 :
156 : \plumedfile
157 : Q6 SPECIESA=1-5 SPECIESB=1-64 D_0=1.3 R_0=0.2 LABEL=q6a
158 : Q6 SPECIESA=6-64 SPECIESB=1-64 D_0=1.3 R_0=0.2 LABEL=q6b
159 :
160 : LOCAL_Q6 SPECIES=q4a,q4b SWITCH={RATIONAL D_0=1.3 R_0=0.2} MEAN LOWMEM LABEL=w4
161 : PRINT ARG=w6.* FILE=colvar
162 : \endplumedfile
163 :
164 : */
165 : //+ENDPLUMEDOC
166 :
167 : namespace PLMD {
168 : namespace crystallization {
169 :
170 30 : class Q6 : public Steinhardt {
171 : public:
172 : static void registerKeywords( Keywords& keys );
173 : explicit Q6( const ActionOptions& ao );
174 : };
175 :
176 6467 : PLUMED_REGISTER_ACTION(Q6,"Q6")
177 : typedef LocalSteinhardt<Q6> LOCAL_Q6;
178 6457 : PLUMED_REGISTER_ACTION(LOCAL_Q6,"LOCAL_Q6")
179 :
180 16 : void Q6::registerKeywords( Keywords& keys ) {
181 16 : Steinhardt::registerKeywords( keys );
182 16 : }
183 :
184 15 : Q6::Q6(const ActionOptions& ao ):
185 : Action(ao),
186 15 : Steinhardt(ao)
187 : {
188 15 : setAngularMomentum(6);
189 :
190 15 : normaliz.resize( 7 );
191 15 : normaliz[0] = sqrt( ( 13.0*720.0 ) / (4.0*pi*720.0) );
192 15 : normaliz[1] = -sqrt( ( 13.0*120.0 ) / (4.0*pi*5040) );
193 15 : normaliz[2] = sqrt( ( 13.0*24) / (4.0*pi*40320) );
194 15 : normaliz[3] = -sqrt( ( 13.0*6) / (4.0*pi*362880) );
195 15 : normaliz[4] = sqrt( (13.0*2) / (4.0*pi*3628800) );
196 15 : normaliz[5] = -sqrt( (13.0*1) / (4.0*pi*39916800) );
197 15 : normaliz[6] = sqrt( (13.0*1) / (4.0*pi*479001600) );
198 :
199 15 : coeff_poly.resize( 7 );
200 30 : coeff_poly[0]=-0.3125; coeff_poly[1]=0.0;
201 30 : coeff_poly[2]=6.5625; coeff_poly[3]=0.0;
202 30 : coeff_poly[4]=-19.6875; coeff_poly[5]=0.0;
203 15 : coeff_poly[6]=14.4375;
204 15 : }
205 :
206 : }
207 4839 : }
208 :
|