Combining multiple variables in CENTER_OF_MULTICOLVAR

[Jamil]

Hi everyone,

I’m trying to compute the center of excess charge (CEC) using CENTER_OF_MULTICOLVAR, but I’m running into an issue with how to combine multiple transformed coordination numbers.

Here is my input:

UNITS LENGTH=A ENERGY=eV TIME=ps
HYDROGENS: GROUP NDX_FILE=H_water_N3.ndx NDX_GROUP=H_water_N3
water_oxygens: GROUP NDX_FILE=O_water.ndx NDX_GROUP=O_water
cn1: COORDINATIONNUMBER SPECIESA=HYDROGENS SPECIESB=2641,2642,2643,2644,2645 SWITCH={COSINUS D_0=1.1 R_0=0.2} LOWMEM
cn2: COORDINATIONNUMBER SPECIESA=water_oxygens SPECIESB=HYDROGENS SWITCH={COSINUS D_0=1.1 R_0=0.2} LOWMEM
w1: MTRANSFORM_MORE DATA=cn1 SWITCH={COSINUS D_0=0 R_0=1} LOWMEM MEAN LOWEST HIGHEST
w2: MTRANSFORM_MORE DATA=cn2 SWITCH={COSINUS D_0=2 R_0=1} LOWMEM MEAN LOWEST HIGHEST
cec: CENTER_OF_MULTICOLVAR DATA=w
pcec: POSITION ATOM=cec
d: DISTANCE ATOMS=cec,2645
PRINT ARG=d,pcec.x,pcec.y,pcec.z STRIDE=1 FILE=colvars.dat

PLUMED gives the following error:

ERROR in input to action CENTER_OF_MULTICOLVAR with label cec : action labelled w does not exist or does not have vessels

What I want to do is use both w1 and w2 as weights in the CENTER_OF_MULTICOLVAR.

My question is: how should I properly combine w1 and w2 so they can be passed to CENTER_OF_MULTICOLVAR?

Thanks!

[Gareth Tribello]

Hello again

I don't understand the question you are asking.  I am going to rewrite your input using PLUMED 2.11 input syntax to explain the problem.

UNITS LENGTH=A ENERGY=eV TIME=ps

# Create a group with all the hydrogen atoms

HYDROGENS: GROUP NDX_FILE=H_water_N3.ndx NDX_GROUP=H_water_N3

# Create a group with all the oxygen atoms

water_oxygens: GROUP NDX_FILE=O_water.ndx NDX_GROUP=O_water

# Calculate the coordination of the hydrogen atoms with a particular set of atoms

# Notice that the cool kids do this by calculating the contact matrix and summing each row of the matrix by multiplying the matrix by a vector of ones 

cn1_mat: CONTACT_MATRIX GROUPA=HYDROGENS GROUPB=2641,2642,2643,2644,2645 SWITCH={COSINUS D_0=1.1 R_0=0.2}

cn1_ones: ONES SIZE=5

cn1: MATRIX_TIMES_VECTOR ARG=cn1_mat,cn1_ones

# Notice that the cool kids do this by calculating the contact matrix and summing each row of the matrix by multiplying the matrix by a vector of ones 

cn2_mat: CONTACT_MATRIX GROUPA=water_oxygens GROUPB=HYDROGENS SWITCH={COSINUS D_0=1.1 R_0=0.2}

cn2_ones: ONES SIZE=<Number of hydrogens>

cn2: MATRIX_TIMES_VECTOR ARG=cn2_mat,cn2_ones

# We will now do the transformations of the two vectors of coordination numbers by the switching functions

w1: MORE_THAN ARG=cn1 SWITCH={COSINUS D_0=0 R_0=1}

w2: MORE_THAN ARG=cn2 SWITCH={COSINUS D_0=2 R_0=1}

OK great, so we get to the end of the input above and we have two vectors above w1 and w2.  w1 contains a set of transformed coordination numbers of hydrogen atoms, while w2 contains a transformed coordination numbers of water oxygens.

At this stage, I would do either:

c: CENTER WEIGHTS=w1 ATOMS=HYDROGENS

This calculates a center of mass as:

(1/W) \sum_i w_i x_i 

where w is the sum of all the transformed coordination numbers in the vector w1 and x_i is the position of the ith atom whose transformed coordination number is w_i.

The center is thus a weighted average calculated over the positions of the hydrogens using the transformed coordination numbers as weights.

Alternatively, I would do:

c: CENTER WEIGHTS=w2 ATOMS=water_oxygens

This center is thus a weighted average over the positions of the oxygens using the transformed coordination numbers as weights.

There isn't really enough information in your question for me to work out what you mean by "properly combining" w1 and w2?  Do you want to do a weighted average over the positions of all the oxygen and hydrogen atoms?  Do you want to calculate the centers for the two types of atoms separately?  Do you want the weight ascribed to each atom to reflect the weights with the two transformed coordination numbers you have calculated?  It is not very clear.

I hope this helps

Gareth

--
You received this message because you are subscribed to the Google Groups "PLUMED users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to plumed-users...@googlegroups.com.
To view this discussion visit https://groups.google.com/d/msgid/plumed-users/63fd7172-46c3-495d-9555-3bf2e8a12a8dn%40googlegroups.com.

[Jamil]

Hi, thanks for the detailed explanation, that helps a lot.

What I want is to compute a single center that includes both hydrogens and oxygens together, where:

• hydrogens are weighted by w1

• oxygens are weighted by w2

So effectively I am trying to compute something like:

(∑i​wi(H)​xi​+∑j​wj(O)​xj​​) / (∑i​wi(H)​+∑j​wj(O) )

i.e. a weighted average over all atoms (H + O) using their respective transformed coordination numbers.

Thanks again!