Applying LOWER_WALLS to Distance Between a Single Atom and a Group of Atoms in OPES_META

[riyazpu]

Hi
Can anyone help me implementing a LOWER_WALLS restraint in OPES_METAD within PLUMED to control the distance between a single atom (e.g., an oxygen atom, O) and a group of atoms (e.g., 50 hydrogen atoms). When using the DISTANCE collective variable (CV), it appears to be limited to calculating the distance between two individual atoms.

[valsson]

Hello, You can define the center of mass of a group of atoms and use that to define distances; see https://www.plumed.org/doc-v2.9/user-doc/html/_c_o_m.html Best regards, Omar

…

On Wed, Apr 30, 2025 at 6:27 AM Mohd Riyaz ***@***.***> wrote: Hi Can anyone help me implementing a LOWER_WALLS restraint in OPES_METAD within PLUMED to control the distance between a single atom (e.g., an oxygen atom, O) and a group of atoms (e.g., 50 hydrogen atoms). When using the DISTANCE collective variable (CV), it appears to be limited to calculating the distance between two individual atoms. — Reply to this email directly, view it on GitHub <#1254>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/ABNCSLWYQNHVLTP2N42ZT4L24CXTRAVCNFSM6AAAAAB4FMPJCCVHI2DSMVQWIX3LMV43ERDJONRXK43TNFXW4OZYGI3DINZQGU> . You are receiving this because you are subscribed to this thread.Message ID: ***@***.***>

[riyazpu]

Thank you for your response. Actually, my system involves water, and defining the centre of mass (COM) does not guarantee that hydrogen atoms will avoid close proximity (cutoff distance 1.5 Angstrom) to the atom of interest.

The idea is that during the MD run H-atom of water can form H-bond with O of my second molecule, but should not transfer the H. that why I want to put a LOWER_WALL on the distance.

Best
Riyaz

[invemichele]

Hi Riyaz, Then I am afraid you have to define a different wall for each H atom. You should also look into coordination https://www.plumed.org/doc-v2.9/user-doc/html/_c_o_o_r_d_i_n_a_t_i_o_n.html, if you assume only one H-bond can form at the time, this might do the job. It's also probably faster since it supports neighbor lists. Best, Michele

…

On Wed, 30 Apr 2025, 14:31 Mohd Riyaz, ***@***.***> wrote: Thank you for your response. Actually, my system involves water, and defining the centre of mass (COM) does not guarantee that hydrogen atoms will avoid close proximity (cutoff distance 1.5 Angstrom) to the atom of interest. The idea is that during the MD run H-atom of water can form H-bond with O of my second molecule, but should not transfer the H. that why I want to put a LOWER_WALL on the distance. Best Riyaz — Reply to this email directly, view it on GitHub <#1254 (reply in thread)>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/ADRW265K5NT7BV2DYVS33DD24C7BPAVCNFSM6AAAAAB4FMPJCCVHI2DSMVQWIX3LMV43URDJONRXK43TNFXW4Q3PNVWWK3TUHMYTEOJZGI4DOMQ> . You are receiving this because you are subscribed to this thread.Message ID: ***@***.***>

[riyazpu]

Hi Michele
Thank you for the response. I was trying different setups like this one
"setup = [
# 1. Units (fixed syntax)
f"UNITS LENGTH=A TIME={1/(1000 * units.fs)} ENERGY=eV",
# CO atoms
"o: GROUP ATOMS=81",
"c: GROUP ATOMS=82",
"h: GROUP ATOMS=83-182,233,234,236,237,239,240,242,243",
# CO center of mass
"com_co: COM ATOMS=o,c",
# CO dist
"d_co: DISTANCE ATOMS=o,c",
# surface COM
"com_surf: COM ATOMS=56,54,35",
# Dist from surf
"d_co_surf: DISTANCE ATOMS=com_surf,c NOPBC",
# Coordination
"coord: COORDINATION GROUPA=85,86 GROUPB=1-80 R_0=2.5",
"coord_o_water: COORDINATION GROUPA=85 GROUPB=h R_0=1.5",
# Compute OH distances and extract minimum
"dOH: DISTANCES GROUPA=o GROUPB=h MIN={BETA=500.0}",
# Apply lower wall on the minimum distance
"dwallh: LOWER_WALLS ARG=dOH.min AT=1.3 KAPPA=30",
# Define upper walls
"dwall: UPPER_WALLS ARG=d_co AT=2.1 KAPPA=30 EXP=2",
"dwallb: UPPER_WALLS ARG=d_co_surf AT=3 KAPPA=30 EXP=2",# (distance CO)
# 6. METAD (traditional metadynamics)
"OPES_METAD " +
"LABEL=opes " +
"ARG=d_co,coord " +
"SIGMA=0.025,0.15 " +
"PACE=800 " +
"BARRIER=1.8 " +
"TEMP=360 " +
"FILE=KERNELS " +
"STATE_WFILE=STATE " +
"STATE_RFILE=STATE " +
"EXTRA_BIAS=dwall.bias,dwallb.bias, dwallh.bias " +
"STATE_WSTRIDE=300100",
# 7. Monitoring
"PRINT ARG=,opes.bias STRIDE=100 FILE=COLVAR""

Please check the CV and LOWER_WALLS "# Compute OH distances and extract minimum
"dOH: DISTANCES GROUPA=o GROUPB=h MIN={BETA=500.0}",
# Apply lower wall at 1.0 ?~E = 0.1 nm on the minimum distance
"dwallh: LOWER_WALLS ARG=dOH.min AT=1.3 KAPPA=30","
It's a part of the multicolvar, but for the time being it seems to be working. But I don't know if I am doing it wrong.
I will try with the coordination number also.
Thanks

Best
Riyaz

[invemichele]

Hi, Just a small suggestion. I would remove the EXTRA_BIAS line. OPES will not push hard on the walls anyway, since it's limited by the BARRIER. If the repulsive walls are steep, reweighing for them might introduce some large noise. In general I would not use EXTRA_BIAS except in rare cases, like if the wall is on the same CV pushed by OPES and it's really critical to stay below a certain value. I will add something about this to the documentation. Also, consider using OPES_METAD_EXPLORE if you are having troubles going back and forth. Best, Michele

…

On Thu, May 1, 2025 at 8:24 PM Mohd Riyaz ***@***.***> wrote: Hi Michele Thank you for the response. I was trying different setups like this one "setup = [ # 1. Units (fixed syntax) f"UNITS LENGTH=A TIME={1/(1000 * units.fs)} ENERGY=eV", # CO atoms "o: GROUP ATOMS=81", "c: GROUP ATOMS=82", "h: GROUP ATOMS=83-182,233,234,236,237,239,240,242,243", # CO center of mass "com_co: COM ATOMS=o,c", # CO dist "d_co: DISTANCE ATOMS=o,c", # surface COM "com_surf: COM ATOMS=56,54,35", # Dist from surf "d_co_surf: DISTANCE ATOMS=com_surf,c NOPBC", # Coordination "coord: COORDINATION GROUPA=85,86 GROUPB=1-80 R_0=2.5", "coord_o_water: COORDINATION GROUPA=85 GROUPB=h R_0=1.5", # Compute OH distances and extract minimum "dOH: DISTANCES GROUPA=o GROUPB=h MIN={BETA=500.0}", # Apply lower wall on the minimum distance "dwallh: LOWER_WALLS ARG=dOH.min AT=1.3 KAPPA=30", # Define upper walls "dwall: UPPER_WALLS ARG=d_co AT=2.1 KAPPA=30 EXP=2", "dwallb: UPPER_WALLS ARG=d_co_surf AT=3 KAPPA=30 EXP=2",# (distance CO) # 6. METAD (traditional metadynamics) "OPES_METAD " + \ "LABEL=opes " + \ "ARG=d_co,coord " + \ "SIGMA=0.025,0.15 " + \ "PACE=800 " + \ "BARRIER=1.8 " + \ "TEMP=360 " + \ "FILE=KERNELS " + \ "STATE_WFILE=STATE " + \ "STATE_RFILE=STATE " + \ "EXTRA_BIAS=dwall.bias,dwallb.bias, dwallh.bias " + \ "STATE_WSTRIDE=300*100", # 7. Monitoring "PRINT ARG=*,opes.bias STRIDE=100 FILE=COLVAR"" Please check the CV and LOWER_WALLS "# Compute OH distances and extract minimum "dOH: DISTANCES GROUPA=o GROUPB=h MIN={BETA=500.0}", # Apply lower wall at 1.0 ?~E = 0.1 nm on the minimum distance "dwallh: LOWER_WALLS ARG=dOH.min AT=1.3 KAPPA=30"," It's a part of the multicolvar, but for the time being it seems to be working. But I don't know if I am doing it wrong. I will try with the coordination number also. Thanks Best Riyaz — Reply to this email directly, view it on GitHub <#1254 (reply in thread)>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/ADRW26YOZU66EI453UW7HA324JREBAVCNFSM6AAAAAB4FMPJCCVHI2DSMVQWIX3LMV43URDJONRXK43TNFXW4Q3PNVWWK3TUHMYTGMBQGY2DSMQ> . You are receiving this because you commented.Message ID: ***@***.***>

[riyazpu]

Thanks for the suggestions.
I have another question, I am running a molecular dynamics simulation for CO dissociation on a metal surface in vacuum using OPES_META. The free energy barrier obtained from COLVAR by reweighting is consistently higher than the barrier derived from the STATE file, even when I run the simulation for more than 1 ns. Could you suggest possible reasons for this discrepancy?

[invemichele]

Hi, The barrier is the hardest part of the FES to estimate, because it's where one collects fewer samples, thus it has the largest uncertainty. It could just be an artifact of the kernel density estimation. you should try using several different bandwidth values, and see if there is a range where the FES estimate looks similar. Another possible thing you can try is to try to remove an initial transient of the simulation and see if the FES estimate changes. If you are interested specifically in the transition state you could also try some more tailored method, like https://arxiv.org/abs/2410.17029 Best, Michele

…

On Tue, 6 May 2025, 17:34 Mohd Riyaz, ***@***.***> wrote: Thanks for the suggestions. I have another question, I am running a molecular dynamics simulation for CO dissociation on a metal surface in vacuum using OPES_META. The free energy barrier obtained from COLVAR by reweighting is consistently higher than the barrier derived from the STATE file, even when I run the simulation for more than 1 ns. Could you suggest possible reasons for this discrepancy? — Reply to this email directly, view it on GitHub <#1254 (comment)>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/ADRW263U6RLWYY6FDN7TBNL25DJA3AVCNFSM6AAAAAB4FMPJCCVHI2DSMVQWIX3LMV43URDJONRXK43TNFXW4Q3PNVWWK3TUHMYTGMBVGE4DMNI> . You are receiving this because you commented.Message ID: ***@***.***>

[riyazpu]

Thankyou!!