Skip to content

Action: CONTACT_MATRIX

Module adjmat
Description Usage
Adjacency matrix in which two atoms are adjacent if they are within a certain cutoff. used in 3 tutorialsused in 8 eggs
output value type
a matrix containing the weights for the bonds between each pair of atoms matrix

Details and examples

Adjacency matrix in which two atoms are adjacent if they are within a certain cutoff.

A useful tool for developing complex collective variables is the notion of the so called adjacency matrix. An adjacency matrix can be an matrix in which the th, th element tells you whether or not the th and th atoms/molecules from a set of atoms/molecules are adjacent or not. Alternatively, you can calculate an adjacency matrix between a set of atoms and a second set of atoms. For this type of matrix the th, th element tells you whether the whether the th atom in the first group and the th atom in the second group are adjacent or not. The adjacency matrix in this case is thus .

The simplest type of adjacency matrix is a contact matrix. The elements of a contact matrix are calculated using:

where is the magnitude of the vector connecting atoms and and where is a switching function. If you want to calculate a contact matrix for one group of atoms the input would look something like this:

Click on the labels of the actions for more information on what each action computes
tested on2.11
c1: CONTACT_MATRIXAdjacency matrix in which two atoms are adjacent if they are within a certain cutoff. More details GROUPspecifies the list of atoms that should be assumed indistinguishable=1-7 SWITCHthe input for the switching function that acts upon the distance between each pair of atoms. Options for this keyword are explained in the documentation for LESS_THAN.={RATIONAL R_0=2.6 NN=6 MM=12}

Alternatively, if you want to calculate the contact matrix between two groups of atoms you would use an input like following:

Click on the labels of the actions for more information on what each action computes
tested on2.11
c2: CONTACT_MATRIXAdjacency matrix in which two atoms are adjacent if they are within a certain cutoff. More details GROUPAwhen you are calculating the adjacency matrix between two sets of atoms this keyword is used to specify the atoms along with the keyword GROUPB=1-7 GROUPBwhen you are calculating the adjacency matrix between two sets of atoms this keyword is used to specify the atoms along with the keyword GROUPA=8-20 SWITCHthe input for the switching function that acts upon the distance between each pair of atoms. Options for this keyword are explained in the documentation for LESS_THAN.={RATIONAL R_0=2.6 NN=6 MM=12}

Once you have calculated a contact matrix you can perform various matrix operations by using the tools in the matrixtools or clusters modules.

Coordination numbers

If a contact matrix, , is multiplied from the back by a vector of ones the th element of the resulting matrix tells you the number of atoms that are within of atom . In other words, the coordination numbers of the atoms can be calculated from the contact matrix by doing matrix vector multiplication.

This realisation about the relationship between the contact map and the coordination number is heavily used in PLUMED. For example, to calculate and print the coordination numbers of the first 7 atoms in the system with themselves you would use an input something like this:

Click on the labels of the actions for more information on what each action computes
tested on2.11
c1: CONTACT_MATRIXAdjacency matrix in which two atoms are adjacent if they are within a certain cutoff. More details GROUPspecifies the list of atoms that should be assumed indistinguishable=1-7 D_0 The d_0 parameter of the switching function=0.0 R_0The r_0 parameter of the switching function=2.6 NN The n parameter of the switching function =6 MM The m parameter of the switching function; 0 implies 2*NN=12
ones: ONESCreate a constant vector with all elements equal to one This action is a shortcut. More details SIZEthe number of ones that you would like to create=7
cc: MATRIX_VECTOR_PRODUCTCalculate the product of the matrix and the vector More details ARGthe label for the matrix and the vector/scalar that are being multiplied=c1,ones
PRINTPrint quantities to a file. More details ARGthe labels of the values that you would like to print to the file=cc FILEthe name of the file on which to output these quantities=colvar

This input transforms the distances using the RATIONAL switching function that is discussed in the documentation for LESS_THAN The parameters for this type of switching function can be specified using D_0, R_0, NN and MM keywords as indicated above. However, we would recommend using the following syntax instead as this allows you to use the full range of switching function types. Most importantly, if you use the following syntax you can set the D_MAX parameter, as shown below. As discussed in the optimization details section below, if you want good computational performance on large systems it is essential to set the D_MAX parameter.

Click on the labels of the actions for more information on what each action computes
tested on2.11
c1: CONTACT_MATRIXAdjacency matrix in which two atoms are adjacent if they are within a certain cutoff. More details GROUPspecifies the list of atoms that should be assumed indistinguishable=1-7 SWITCHthe input for the switching function that acts upon the distance between each pair of atoms. Options for this keyword are explained in the documentation for LESS_THAN.={RATIONAL R_0=2.6 NN=6 MM=12 D_MAX=5.0}
ones: ONESCreate a constant vector with all elements equal to one This action is a shortcut. More details SIZEthe number of ones that you would like to create=7
cc: MATRIX_VECTOR_PRODUCTCalculate the product of the matrix and the vector More details ARGthe label for the matrix and the vector/scalar that are being multiplied=c1,ones
PRINTPrint quantities to a file. More details ARGthe labels of the values that you would like to print to the file=cc FILEthe name of the file on which to output these quantities=colvar

Implmenting the coordination number as a product of a matrix and a vector is useful as there are many different ways to define whether two atoms/molecules and to construct a "contact" matrix based on the result. For example:

  • You could say that two molecules are connected if they are within a certain distance of each other and if they have the same orientation (see TORSIONS_MATRIX).
  • You could say that two water molecules are connected if they are hydrogen bonded to each other (see HBOND_MATRIX).
  • You could say that two atoms are connected if they are within a certain distance of each other and if they have similar values for a CV (see OUTER_PRODUCT).

When the coordination numbers is implemented in the way described above (by doing the matrix-vector multiplication) you have the flexibility to define the contact matrix that is used in the multiplication in whatever way you choose. In other words, this implementation of the coordination number is much more flexible. For example, suppose you want to calculate the number of atoms that have a coordination is greater than 3.0. You can do this with PLUMED using the following input:

Click on the labels of the actions for more information on what each action computes
tested on2.11
# Calculate the contact matrix for the first seven atoms in the system
c1: CONTACT_MATRIXAdjacency matrix in which two atoms are adjacent if they are within a certain cutoff. More details GROUPspecifies the list of atoms that should be assumed indistinguishable=1-7 SWITCHthe input for the switching function that acts upon the distance between each pair of atoms. Options for this keyword are explained in the documentation for LESS_THAN.={RATIONAL R_0=2.6 NN=6 MM=12}
# Calculate the coordination numbers for the first seven atoms in the system
ones: ONESCreate a constant vector with all elements equal to one This action is a shortcut. More details SIZEthe number of ones that you would like to create=7
cc: MATRIX_VECTOR_PRODUCTCalculate the product of the matrix and the vector More details ARGthe label for the matrix and the vector/scalar that are being multiplied=c1,ones
# Set the ith element of the vector mtc equal to one if the coordination number of atom i is greater than 3.
mtc: MORE_THANUse a switching function to determine how many of the input variables are more than a certain cutoff. More details ARGthe values input to this function=cc SWITCHThis keyword is used if you want to employ an alternative to the continuous swiching function defined above={RATIONAL D_0=3 R_0=1}
# Calculate the number of atoms with a coordination number greater than 3.
s: SUMCalculate the sum of the arguments More details ARGthe vector/matrix/grid whose elements shuld be added together=mtc PERIODICif the output of your function is periodic then you should specify the periodicity of the function=NO
PRINTPrint quantities to a file. More details ARGthe labels of the values that you would like to print to the file=s FILEthe name of the file on which to output these quantities=colvar

Alternatively, consider the CV that was used to study perovskite nucleation in this paper. The CV here measures the number of methylamonium molecules that are attached to at least 5 other methylamoniusm molecules, at least 7 lead atoms and at least 11 ionide ions. We can calculate something akin to this CV and apply a restraint on the resulting quantity by using the following input file:

Click on the labels of the actions for more information on what each action computes
tested on2.11
# Lead ions
Pb: GROUPDefine 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. More details ATOMSthe numerical indexes for the set of atoms in the group=1-64
# Iodide atoms
I: GROUPDefine 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. More details ATOMSthe numerical indexes for the set of atoms in the group=65-256
# Methylamonium "atoms" -- in the real CV these are centers of mass rather than single atoms
cn: GROUPDefine 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. More details ATOMSthe numerical indexes for the set of atoms in the group=257-320

ones64: ONESCreate a constant vector with all elements equal to one This action is a shortcut. More details SIZEthe number of ones that you would like to create=64
# Contact matrix that determines if methylamonium molecules are within 8 A of each other
cm_cncn: CONTACT_MATRIXAdjacency matrix in which two atoms are adjacent if they are within a certain cutoff. More details GROUPspecifies the list of atoms that should be assumed indistinguishable=cn SWITCHthe input for the switching function that acts upon the distance between each pair of atoms. Options for this keyword are explained in the documentation for LESS_THAN.={RATIONAL R_0=0.8}
# Coordination number of methylamounium with methylamonium
cc_cncn: MATRIX_VECTOR_PRODUCTCalculate the product of the matrix and the vector More details ARGthe label for the matrix and the vector/scalar that are being multiplied=cm_cncn,ones64
# Vector with elements that are one if coordiantion of methylamonium with methylamonium >5
mt_cncn: MORE_THANUse a switching function to determine how many of the input variables are more than a certain cutoff. More details ARGthe values input to this function=cc_cncn SWITCHThis keyword is used if you want to employ an alternative to the continuous swiching function defined above={RATIONAL R_0=5 NN=12 MM=24}

# Contact matrix that determines if methylamoinium moleulcule and Pb atom are within 7.5 A of each other
cm_cnpb: CONTACT_MATRIXAdjacency matrix in which two atoms are adjacent if they are within a certain cutoff. More details GROUPAwhen you are calculating the adjacency matrix between two sets of atoms this keyword is used to specify the atoms along with the keyword GROUPB=cn GROUPBwhen you are calculating the adjacency matrix between two sets of atoms this keyword is used to specify the atoms along with the keyword GROUPA=Pb SWITCHthe input for the switching function that acts upon the distance between each pair of atoms. Options for this keyword are explained in the documentation for LESS_THAN.={RATIONAL R_0=0.75}
# Coordination number of methylamonium with Pb
cc_cnpb: MATRIX_VECTOR_PRODUCTCalculate the product of the matrix and the vector More details ARGthe label for the matrix and the vector/scalar that are being multiplied=cm_cnpb,ones64
# Vector with elements that are one if coordination of methylamounium with lead is >7
mt_cnpb: MORE_THANUse a switching function to determine how many of the input variables are more than a certain cutoff. More details ARGthe values input to this function=cc_cnpb SWITCHThis keyword is used if you want to employ an alternative to the continuous swiching function defined above={RATIONAL R_0=7 NN=12 MM=24}

ones192: ONESCreate a constant vector with all elements equal to one This action is a shortcut. More details SIZEthe number of ones that you would like to create=192
# Contact matrix that determines if methylamoinium moleulcule and I atom are within 6.5 A of each other
cm_cnI: CONTACT_MATRIXAdjacency matrix in which two atoms are adjacent if they are within a certain cutoff. More details GROUPAwhen you are calculating the adjacency matrix between two sets of atoms this keyword is used to specify the atoms along with the keyword GROUPB=cn GROUPBwhen you are calculating the adjacency matrix between two sets of atoms this keyword is used to specify the atoms along with the keyword GROUPA=I SWITCHthe input for the switching function that acts upon the distance between each pair of atoms. Options for this keyword are explained in the documentation for LESS_THAN.={RATIONAL R_0=0.65}
# Coordination number of methylamonium with I
cc_cnI: MATRIX_VECTOR_PRODUCTCalculate the product of the matrix and the vector More details ARGthe label for the matrix and the vector/scalar that are being multiplied=cm_cnI,ones192
# Vector with elements that are one if coordination of methylamounium with lead is >11
mt_cnI: MORE_THANUse a switching function to determine how many of the input variables are more than a certain cutoff. More details ARGthe values input to this function=cc_cnI SWITCHThis keyword is used if you want to employ an alternative to the continuous swiching function defined above={RATIONAL R_0=11 NN=12 MM=24}

# Element wise product of these three input vectors.
# mm[i]==1 if coordination number of corrsponding methylamounium with methylamonium is >5
# and if coordination of methylamounium with Pb is >7 and if coordination of methylamounium with I > 11
mm: CUSTOMCalculate a combination of variables using a custom expression. More details ARGthe values input to this function=mt_cncn,mt_cnpb,mt_cnI FUNCthe function you wish to evaluate=x*y*z PERIODICif the output of your function is periodic then you should specify the periodicity of the function=NO

# Sum of coordination numbers and thus equal to number of methylamoniums with desired coordination numbers
ff: SUMCalculate the sum of the arguments More details ARGthe vector/matrix/grid whose elements shuld be added together=mm PERIODICif the output of your function is periodic then you should specify the periodicity of the function=NO

rr: RESTRAINTAdds harmonic and/or linear restraints on one or more variables. This action has hidden defaults. More details ARGthe values the harmonic restraint acts upon=ff ATthe position of the restraint=62 KAPPA specifies that the restraint is harmonic and what the values of the force constants on each of the variables are=10

COMPONENTS flag

If you add the flag COMPONENTS to the input as shown below:

Click on the labels of the actions for more information on what each action computes
tested on2.11
c4: CONTACT_MATRIXAdjacency matrix in which two atoms are adjacent if they are within a certain cutoff. More details GROUPspecifies the list of atoms that should be assumed indistinguishable=1-7 COMPONENTS also calculate the components of the vector connecting the atoms in the contact matrix SWITCHthe input for the switching function that acts upon the distance between each pair of atoms. Options for this keyword are explained in the documentation for LESS_THAN.={RATIONAL R_0=2.6 NN=6 MM=12}

then four matrices with the labels c4.w, c4.x, c4.y and c4.z are output by the action. The matrix with the label c4.w is the adjacency matrix that would be output if you had not added the COMPONENTS flag. The component of the matrices c4.x, c4.y and c4.z contain the , and components of the vector connecting atoms and . Importantly, however, the components of these vectors are only stored in c4.x, c4.y and c4.z if the elements of c4.w are non-zero.

NOPBC flag

By default PLUMED calculates the distances that are transformed by the switching function in the CONTACT_MATRIX action in a way that accounts for the periodic boundary conditions. If for any reason you do not want PLUMED to account for the periodic boundary conditions you can use the NOPBC flag as shown below:

Click on the labels of the actions for more information on what each action computes
tested on2.11
c: CONTACT_MATRIXAdjacency matrix in which two atoms are adjacent if they are within a certain cutoff. More details GROUPspecifies the list of atoms that should be assumed indistinguishable=1-7 NOPBC don't use pbc SWITCHthe input for the switching function that acts upon the distance between each pair of atoms. Options for this keyword are explained in the documentation for LESS_THAN.={RATIONAL R_0=2.6 NN=6 MM=12}

If you also add the COMPONENTS flag in this input then the vector connecting each pairs of atoms will also be calculated without taking the periodic boundary conditions into account.

Optimisation details

Adjacency matrices are sparse. Each atom is only be connected to a small number of neighbours and the vast majority of the elements of the contact matrix are thus zero. To reduce the amount of memory that PLUMED requires PLUMED uses sparse matrix storage. Consequently, whenever you calculate and store a contact matrix only the elements of the matrix that are non-zero are stored. The same thing holds for the additional values that are created when you use the COMPONENTS flag. The components of the vectors connecting atoms are only stored when the elements of c4.w are non-zero.

We can also use the sparsity of the adjacency matrix to make the time required to compute a contact matrix scale linearly rather than quadratically with the number of atoms. Element of the contact matrix is only non-zero if two atoms are within a cutoff, . We can determine that many pairs of atoms are further appart than without computing the distance between these atoms by using divide and conquer strategies such as linked lists and neighbour lists. To turn on these features you need to set the D_MAX parameter in the switching functions. The value you pass to the D_MAX keyword is used as the cutoff in the link cell algorithm.

In theory we could further optimize the implementation of the CONTACT_MATRIX action by exploiting neighbor lists. If we were to do this we would likely add two further keywords as shown below:

Click on the labels of the actions for more information on what each action computes
tested on2.11
c4: CONTACT_MATRIXAdjacency matrix in which two atoms are adjacent if they are within a certain cutoff. More details GROUPspecifies the list of atoms that should be assumed indistinguishable=1-10000 SWITCHthe input for the switching function that acts upon the distance between each pair of atoms. Options for this keyword are explained in the documentation for LESS_THAN.={RATIONAL R_0=0.1 NN=6 MM=12 D_MAX=0.5} NL_CUTOFF The cutoff for the neighbor list=0.7 NL_STRIDE The frequency with which we are updating the atoms in the neighbor list=5

The NL_CUTOFF keyword would be used to specify the cutoff (in nm) to use when constructing neighbor lists. This value would need to be slightly larger than the D_MAX parameter for the switching function. The NL_STRIDE keyword would then be used to specify how frequently the neighbour list should be updated. Thus far we have not found it necessary to implement this algorithm. We have been happy with the performance even if we use the linked list algorithm to update the neighbors on every step. If you feel that you need this CV to perform better please get in touch as adding a neighbor list for this action should be relatively straightforward.

The MASK keyword

Suppose that you want to calculate the average coordination number for the atoms that are within a sphere in the center of your simulation box. You can do so by exploiting an input similar to the one shown below:

Click on the labels of the actions for more information on what each action computes
tested on2.11
# The atoms that are of interest
ow: GROUPDefine 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. More details ATOMSthe numerical indexes for the set of atoms in the group=1-16500
# Fixed virtual atom which serves as the probe volume's center (pos. in nm)
center: FIXEDATOMAdd a virtual atom in a fixed position. This action has hidden defaults. More details ATcoordinates of the virtual atom=2.5,2.5,2.5
# Vector in which element i is one if atom i is in sphere of interest and zero otherwise
sphere: INSPHEREThis quantity can be used to calculate functions of the distribution of collective variables for the atoms that lie in a particular, user-specified part of of the cell. More details ATOMSthe group of atoms that you would like to investigate=ow CENTERthe atom whose vicinity we are interested in examining=center RADIUSthe switching function that tells us the extent of the sphereical region of interest. Options for this keyword are explained in the documentation for LESS_THAN.={GAUSSIAN D_0=0.5 R_0=0.01 D_MAX=0.52}
# Calculates cooordination numbers
cmap: CONTACT_MATRIXAdjacency matrix in which two atoms are adjacent if they are within a certain cutoff. More details GROUPspecifies the list of atoms that should be assumed indistinguishable=ow SWITCHthe input for the switching function that acts upon the distance between each pair of atoms. Options for this keyword are explained in the documentation for LESS_THAN.={GAUSSIAN D_0=0.32 R_0=0.01 D_MAX=0.34}
ones: ONESCreate a constant vector with all elements equal to one This action is a shortcut. More details SIZEthe number of ones that you would like to create=16500
cc: MATRIX_VECTOR_PRODUCTCalculate the product of the matrix and the vector More details ARGthe label for the matrix and the vector/scalar that are being multiplied=cmap,ones
# Multiply coordination numbers by sphere vector
prod: CUSTOMCalculate a combination of variables using a custom expression. More details ARGthe values input to this function=cc,sphere FUNCthe function you wish to evaluate=x*y PERIODICif the output of your function is periodic then you should specify the periodicity of the function=NO
# Sum of coordination numbers for atoms that are in the sphere of interest
numer: SUMCalculate the sum of the arguments More details ARGthe vector/matrix/grid whose elements shuld be added together=prod PERIODICif the output of your function is periodic then you should specify the periodicity of the function=NO
# Number of atoms that are in sphere of interest
denom: SUMCalculate the sum of the arguments More details ARGthe vector/matrix/grid whose elements shuld be added together=sphere PERIODICif the output of your function is periodic then you should specify the periodicity of the function=NO
# Average coordination number for atoms in sphere of interest
av: CUSTOMCalculate a combination of variables using a custom expression. More details ARGthe values input to this function=prod,sphere FUNCthe function you wish to evaluate=x/y PERIODICif the output of your function is periodic then you should specify the periodicity of the function=NO
# And print out final CV to a file
PRINTPrint quantities to a file. More details ARGthe labels of the values that you would like to print to the file=av FILEthe name of the file on which to output these quantities=colvar STRIDE the frequency with which the quantities of interest should be output=1

This calculation is slow because you have to calculate the coordination numbers of all the atoms even though only a small subset of these quanitties are required to compute the average coordination number in the sphere. To avoid all these unecessary calculations you use the MASK keyword as shown below:

Click on the labels of the actions for more information on what each action computes
tested on2.11
# The atoms that are of interest
ow: GROUPDefine 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. More details ATOMSthe numerical indexes for the set of atoms in the group=1-16500
# Fixed virtual atom which serves as the probe volume's center (pos. in nm)
center: FIXEDATOMAdd a virtual atom in a fixed position. This action has hidden defaults. More details ATcoordinates of the virtual atom=2.5,2.5,2.5
# Vector in which element i is one if atom i is in sphere of interest and zero otherwise
sphere: INSPHEREThis quantity can be used to calculate functions of the distribution of collective variables for the atoms that lie in a particular, user-specified part of of the cell. More details ATOMSthe group of atoms that you would like to investigate=ow CENTERthe atom whose vicinity we are interested in examining=center RADIUSthe switching function that tells us the extent of the sphereical region of interest. Options for this keyword are explained in the documentation for LESS_THAN.={GAUSSIAN D_0=0.5 R_0=0.01 D_MAX=0.52}
# Calculates cooordination numbers
cmap: CONTACT_MATRIXAdjacency matrix in which two atoms are adjacent if they are within a certain cutoff. More details GROUPspecifies the list of atoms that should be assumed indistinguishable=ow SWITCHthe input for the switching function that acts upon the distance between each pair of atoms. Options for this keyword are explained in the documentation for LESS_THAN.={GAUSSIAN D_0=0.32 R_0=0.01 D_MAX=0.34} MASKa vector that is used to used to determine which rows of the adjancency matrix to compute=sphere
ones: ONESCreate a constant vector with all elements equal to one This action is a shortcut. More details SIZEthe number of ones that you would like to create=16500
cc: MATRIX_VECTOR_PRODUCTCalculate the product of the matrix and the vector More details ARGthe label for the matrix and the vector/scalar that are being multiplied=cmap,ones
# Multiply coordination numbers by sphere vector
prod: CUSTOMCalculate a combination of variables using a custom expression. More details ARGthe values input to this function=cc,sphere FUNCthe function you wish to evaluate=x*y PERIODICif the output of your function is periodic then you should specify the periodicity of the function=NO
# Sum of coordination numbers for atoms that are in the sphere of interest
numer: SUMCalculate the sum of the arguments More details ARGthe vector/matrix/grid whose elements shuld be added together=prod PERIODICif the output of your function is periodic then you should specify the periodicity of the function=NO
# Number of atoms that are in sphere of interest
denom: SUMCalculate the sum of the arguments More details ARGthe vector/matrix/grid whose elements shuld be added together=sphere PERIODICif the output of your function is periodic then you should specify the periodicity of the function=NO
# Average coordination number for atoms in sphere of interest
av: CUSTOMCalculate a combination of variables using a custom expression. More details ARGthe values input to this function=prod,sphere FUNCthe function you wish to evaluate=x/y PERIODICif the output of your function is periodic then you should specify the periodicity of the function=NO
# And print out final CV to a file
PRINTPrint quantities to a file. More details ARGthe labels of the values that you would like to print to the file=av FILEthe name of the file on which to output these quantities=colvar STRIDE the frequency with which the quantities of interest should be output=1

Adding the instruction MASK=sphere to the CONTACT_MATRIX line in this input tells PLUMED to only calculate the th row in the adjacency matrix if the th element of the vector sphere is non-zero. In other words, by adding this command we have ensured that we are not calculating coordination numbers for atoms that are not in the sphere that is of interest. In this way we can thus reduce the computational expense of the calculation enormously.

Notice, that there are other places where we can use this same trick. For example, we could have used MASK as shown below for our calculation of the CV for perovskite nucleation as shown below:

Click on the labels of the actions for more information on what each action computes
tested on2.11
# Lead ions
Pb: GROUPDefine 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. More details ATOMSthe numerical indexes for the set of atoms in the group=1-64
# Iodide atoms
I: GROUPDefine 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. More details ATOMSthe numerical indexes for the set of atoms in the group=65-256
# Methylamonium "atoms" -- in the real CV these are centers of mass rather than single atoms
cn: GROUPDefine 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. More details ATOMSthe numerical indexes for the set of atoms in the group=257-320

ones192: ONESCreate a constant vector with all elements equal to one This action is a shortcut. More details SIZEthe number of ones that you would like to create=192
# Contact matrix that determines if methylamoinium moleulcule and I atom are within 6.5 A of each other
cm_cnI: CONTACT_MATRIXAdjacency matrix in which two atoms are adjacent if they are within a certain cutoff. More details GROUPAwhen you are calculating the adjacency matrix between two sets of atoms this keyword is used to specify the atoms along with the keyword GROUPB=cn GROUPBwhen you are calculating the adjacency matrix between two sets of atoms this keyword is used to specify the atoms along with the keyword GROUPA=I SWITCHthe input for the switching function that acts upon the distance between each pair of atoms. Options for this keyword are explained in the documentation for LESS_THAN.={RATIONAL R_0=0.65 D_MAX=1.0}
# Coordination number of methylamonium with I
cc_cnI: MATRIX_VECTOR_PRODUCTCalculate the product of the matrix and the vector More details ARGthe label for the matrix and the vector/scalar that are being multiplied=cm_cnI,ones192
# Vector with elements that are one if coordination of methylamounium with lead is >11
mt_cnI: MORE_THANUse a switching function to determine how many of the input variables are more than a certain cutoff. More details ARGthe values input to this function=cc_cnI SWITCHThis keyword is used if you want to employ an alternative to the continuous swiching function defined above={RATIONAL R_0=11 NN=12 MM=24}

ones64: ONESCreate a constant vector with all elements equal to one This action is a shortcut. More details SIZEthe number of ones that you would like to create=64
# Contact matrix that determines if methylamoinium moleulcule and Pb atom are within 7.5 A of each other
cm_cnpb: CONTACT_MATRIXAdjacency matrix in which two atoms are adjacent if they are within a certain cutoff. More details GROUPAwhen you are calculating the adjacency matrix between two sets of atoms this keyword is used to specify the atoms along with the keyword GROUPB=cn GROUPBwhen you are calculating the adjacency matrix between two sets of atoms this keyword is used to specify the atoms along with the keyword GROUPA=Pb SWITCHthe input for the switching function that acts upon the distance between each pair of atoms. Options for this keyword are explained in the documentation for LESS_THAN.={RATIONAL R_0=0.75 D_MAX=1.5} MASKa vector that is used to used to determine which rows of the adjancency matrix to compute=mt_cnI
# Coordination number of methylamonium with Pb
cc_cnpb: MATRIX_VECTOR_PRODUCTCalculate the product of the matrix and the vector More details ARGthe label for the matrix and the vector/scalar that are being multiplied=cm_cnpb,ones64
# Vector with elements that are one if coordination of methylamounium with lead is >7
mt_cnpb: MORE_THANUse a switching function to determine how many of the input variables are more than a certain cutoff. More details ARGthe values input to this function=cc_cnpb SWITCHThis keyword is used if you want to employ an alternative to the continuous swiching function defined above={RATIONAL R_0=7 NN=12 MM=24}

# Contact matrix that determines if methylamonium molecules are within 8 A of each other
cm_cncn: CONTACT_MATRIXAdjacency matrix in which two atoms are adjacent if they are within a certain cutoff. More details GROUPspecifies the list of atoms that should be assumed indistinguishable=cn SWITCHthe input for the switching function that acts upon the distance between each pair of atoms. Options for this keyword are explained in the documentation for LESS_THAN.={RATIONAL R_0=0.8 D_MAX=1.75} MASKa vector that is used to used to determine which rows of the adjancency matrix to compute=mt_cnI
# Coordination number of methylamounium with methylamonium
cc_cncn: MATRIX_VECTOR_PRODUCTCalculate the product of the matrix and the vector More details ARGthe label for the matrix and the vector/scalar that are being multiplied=cm_cncn,ones64
# Vector with elements that are one if coordiantion of methylamonium with methylamonium >5
mt_cncn: MORE_THANUse a switching function to determine how many of the input variables are more than a certain cutoff. More details ARGthe values input to this function=cc_cncn SWITCHThis keyword is used if you want to employ an alternative to the continuous swiching function defined above={RATIONAL R_0=5 NN=12 MM=24}

# Element wise product of these three input vectors.
# mm[i]==1 if coordination number of corrsponding methylamounium with methylamonium is >5
# and if coordination of methylamounium with Pb is >7 and if coordination of methylamounium with I > 11
mm: CUSTOMCalculate a combination of variables using a custom expression. More details ARGthe values input to this function=mt_cncn,mt_cnpb,mt_cnI FUNCthe function you wish to evaluate=x*y*z PERIODICif the output of your function is periodic then you should specify the periodicity of the function=NO

# Sum of coordination numbers and thus equal to number of methylamoniums with desired coordination numbers
ff: SUMCalculate the sum of the arguments More details ARGthe vector/matrix/grid whose elements shuld be added together=mm PERIODICif the output of your function is periodic then you should specify the periodicity of the function=NO

rr: RESTRAINTAdds harmonic and/or linear restraints on one or more variables. This action has hidden defaults. More details ARGthe values the harmonic restraint acts upon=ff ATthe position of the restraint=62 KAPPA specifies that the restraint is harmonic and what the values of the force constants on each of the variables are=10

This trick works here because when we find that there are methylamoinium moleulcules with fewer than 11 lead atoms in there first coordination sphere we know that there is no point calculating the second two coordination numbers.

Input

The arguments and atoms that serve as the input for this action are specified using one or more of the keywords in the following table.

Keyword Type Description
MASK vector a vector that is used to used to determine which rows of the adjancency matrix to compute
GROUPA atoms when you are calculating the adjacency matrix between two sets of atoms this keyword is used to specify the atoms along with the keyword GROUPB
GROUPB atoms when you are calculating the adjacency matrix between two sets of atoms this keyword is used to specify the atoms along with the keyword GROUPA

Output components

This action can calculate the values in the following table when the associated keyword is included in the input for the action. These values can be referenced elsewhere in the input by using this Action's label followed by a dot and the name of the value required from the list below.

Name Type Keyword Description
w matrix COMPONENTS a matrix containing the weights for the bonds between each pair of atoms
x matrix COMPONENTS the projection of the bond on the x axis
y matrix COMPONENTS the projection of the bond on the y axis
z matrix COMPONENTS the projection of the bond on the z axis

Full list of keywords

The following table describes the keywords and options that can be used with this action

Keyword Type Default Description
MASK input none a vector that is used to used to determine which rows of the adjancency matrix to compute
GROUPA input none when you are calculating the adjacency matrix between two sets of atoms this keyword is used to specify the atoms along with the keyword GROUPB
GROUPB input none when you are calculating the adjacency matrix between two sets of atoms this keyword is used to specify the atoms along with the keyword GROUPA
NL_CUTOFF compulsory 0.0 The cutoff for the neighbor list
NL_STRIDE compulsory 1 The frequency with which we are updating the atoms in the neighbor list
NN compulsory 6 The n parameter of the switching function
MM compulsory 0 The m parameter of the switching function; 0 implies 2*NN
D_0 compulsory 0.0 The d_0 parameter of the switching function
R_0 compulsory none The r_0 parameter of the switching function
COMPONENTS optional false also calculate the components of the vector connecting the atoms in the contact matrix
NOPBC optional false don't use pbc
SERIAL optional false do the calculation in serial. Further information about this flag can be found here.
USEGPU optional false run this calculation on the GPU. Further information about this flag can be found here.
GROUP optional not used specifies the list of atoms that should be assumed indistinguishable
SWITCH optional not used the input for the switching function that acts upon the distance between each pair of atoms. Options for this keyword are explained in the documentation for LESS_THAN.

deprecated keywords

The keywords in the following table can still be used with this action but have been deprecated

Keyword Description
ATOMS You should use GROUP instead of this keyword which was used in older versions of PLUMED and is provided for back compatibility only

References

More information about how this action can be used is available in the following articles: