Generating a 4D FES using sum_hills with multiple HILLS files (one CV per HILLS file)

[Billy Noonan]

Hi Experts,

I am running PBMETAD simulations with multiple walkers using GROMACS 2020 with plumed 2.7.0

As of GROMACS-2020, the MD code does not support the -multi flag and only has the -multidir flag available with mdrun.

Therefore, the WALKERS_MPI or MULTIPLE_WALKERS command cannot be used in the plumed input for this kind of simulation.

Here is my Plumed input:

__________________________

#GROUP DEFINITIONS

backbone_polar: GROUP NDX_FILE=bpep.ndx NDX_GROUP=backbone_polar

lipid_core: GROUP NDX_FILE=bpep.ndx NDX_GROUP=lipid_core

water: GROUP NDX_FILE=bpep.ndx NDX_GROUP=Water_and_ions

HBA_N: GROUP NDX_FILE=bpep.ndx NDX_GROUP=HBA_N

HBA_O: GROUP NDX_FILE=bpep.ndx NDX_GROUP=HBA_O

#GYRATION_CVs

exc_lipid_gyr: GYRATION TYPE=RADIUS ATOMS=lipid_core

exc_pep_gyr: GYRATION TYPE=RADIUS ATOMS=backbone_polar

#CONTACT CVs

HB:  COORDINATION GROUPA=HBA_N GROUPB=HBA_O SWITCH={TANH R_0=0.45 D_MAX=0.5} NLIST NL_CUTOFF=0.5 NL_STRIDE=50

pep_water_contact:  COORDINATION GROUPA=backbone_polar GROUPB=water SWITCH={TANH R_0=0.3 D_MAX=0.4} NLIST NL_CUTOFF=0.4 NL_STRIDE=50

lipid_water_contact:  COORDINATION GROUPA=lipid_core GROUPB=water SWITCH={TANH R_0=0.45 D_MAX=0.5} NLIST NL_CUTOFF=0.5 NL_STRIDE=50

PBMETAD ...

LABEL=metad-4CVs

ARG=exc_lipid_gyr,exc_pep_gyr,HB,lipid_water_contact

PACE=5000

HEIGHT=0.20

SIGMA=0.05,0.05,0.5,10

BIASFACTOR=1000

FILE=4-HILLS-FILE-NAMES-HERE

TEMP=310.0

GRID_MAX=50,50,50,500

GRID_MIN=0.005,0.005,0,0

GRID_BIN=200,200,50,50

WALKERS_N=3

WALKERS_ID=-insert-value-from-0-2-here-

WALKERS_DIR=./

WALKERS_RSTRIDE=1000

... PBMETAD

PRINT  ARG=exc_lipid_gyr,exc_pep_gyr,HB,lipid_water_contact,metad-4CVs.bias FILE=COLVAR-GYR-CONTACT.dat STRIDE=500

__________________________

Can someone please help me generate a 4D free energy surface from either all of the HILLS files or the COLVAR file with Plumed? Currently I have four HILLS files, one for each CV, and my understanding is that I need one or more HILLS file with all four CVs in it to use sum_hills. How do I fix this issue? I cannot use sum_hills with the COLVAR file either, as sum hills returns an error about the --sigma flag, even when I add the correct number of --sigma values:

https://www.plumed.org/doc-v2.6/user-doc/html/sum_hills.html

Trying to generate a histogram from the COLVAR file seems to be possible when using the --idw flag from sum_hills and when also specifying the correct number of --sigma values (one for each CV) but this takes too long.

I have noticed this thread here:

https://groups.google.com/g/plumed-users/c/OSQcKrGVMwg/m/UexO7I5KAgAJ

But as stated the MULTIPLE_WALKERS command in the input file is no longer supported by GROMACS after they removed the -multi command from mdrun after GROMACS-2018

How do I fix my issue?
Any help would be appreciated

Best Wishes,
Billy

Reply all

Reply to author

Forward

[Billy Noonan]

Hi All,

As per advice here:

https://groups.google.com/g/plumed-users/c/Q82ZjZ-ywxE/m/nfTdZ9pRCwAJ

And here:

https://groups.google.com/g/plumed-users/c/fhmDTvJFHYk/m/37pFmKoNAAAJ

I have modified my plumed input to this:

__________

#GROUP DEFINITIONS

backbone_polar: GROUP NDX_FILE=bpep.ndx NDX_GROUP=backbone_polar

lipid_core: GROUP NDX_FILE=bpep.ndx NDX_GROUP=lipid_core

water: GROUP NDX_FILE=bpep.ndx NDX_GROUP=Water_and_ions

HBA_N: GROUP NDX_FILE=bpep.ndx NDX_GROUP=HBA_N

HBA_O: GROUP NDX_FILE=bpep.ndx NDX_GROUP=HBA_O

#GYRATION_CVs

exc_lipid_gyr: GYRATION TYPE=RADIUS ATOMS=lipid_core

exc_pep_gyr: GYRATION TYPE=RADIUS ATOMS=backbone_polar

#CONTACT CVs

HB:  COORDINATION GROUPA=HBA_N GROUPB=HBA_O SWITCH={TANH R_0=0.45 D_MAX=0.5} NLIST NL_CUTOFF=0.5 NL_STRIDE=50

pep_water_contact:  COORDINATION GROUPA=backbone_polar GROUPB=water SWITCH={TANH R_0=0.3 D_MAX=0.4} NLIST NL_CUTOFF=0.4 NL_STRIDE=50

lipid_water_contact:  COORDINATION GROUPA=lipid_core GROUPB=water SWITCH={TANH R_0=0.45 D_MAX=0.5} NLIST NL_CUTOFF=0.5 NL_STRIDE=50

PBMETAD ...

LABEL=metad-4CVs

ARG=exc_lipid_gyr,exc_pep_gyr,HB,lipid_water_contact

PACE=5000

HEIGHT=0.20

SIGMA=0.05,0.05,0.5,10

BIASFACTOR=1000

FILE=../HILLS_4CVs

TEMP=310.0

GRID_MAX=50,50,50,500

GRID_MIN=0.005,0.005,0,0

GRID_BIN=200,200,50,50

WALKERS_MPI

... PBMETAD

PRINT  ARG=exc_lipid_gyr,exc_pep_gyr,HB,lipid_water_contact,metad-4CVs.bias FILE=COLVAR-GYR-CONTACT.dat STRIDE=500

_______

With the following command:

mpirun -np 12 gmx_mpi mdrun -v -deffnm meta-contact -plumed plumed-GYR-CONTACT.dat -multidir WALKER0 WALKER1 WALKER2

Plumed is complaining that I need an equal number of HILLS and CVs

From multiple HILLS, how would I generate an FES in 2+ dimensions, when each HILLS file only has one CV biased?

Plumed sum_hills complains when I try to supply the multiple HILLS files.

Best Wishes,
Billy

[Billy Noonan]

Hi Everyone,

By carefully re-reading the quick replies already on the Plumed mailing list, I have found a script that does what I would expect, and am sharing it here to help others who might have the same question.
If you use PB-MetaD you will need to replace "METAD" with "PB-METAD" in each script, and state one HILLS file for each CV, using the same syntax. With METAD, you can generate an FES in 2+ dimensions by using Plumed sum_hills. By PB-MetaD, you have to reweight each 1D HILLS file into a multi-dimensional COLVAR file as described on the Plumed website. The directory in which you state where the HILLS file(s) are in the Plumed input depends on whether you are using communicating between WALKERS using the file system or whether you are doing this with MPI. You will obviously need one directory in your working directory in your working directory for each WALKER, named as per your gmx mdrun command. All necessary files (.tpr files, plumed inputs and .ndx files) need to be in the WALKER directories. I hope this clarifies everything for people experiencing similar problems to me.

Best Wishes,
Billy

____________________________________________

With contact terms as CVs, and when using MPI to communicate between WALKERS, this should be used as plumed input:

#GROUP DEFINITIONS

backbone_polar: GROUP NDX_FILE=bpep.ndx NDX_GROUP=backbone_polar

lipid_core: GROUP NDX_FILE=bpep.ndx NDX_GROUP=lipid_core

water: GROUP NDX_FILE=bpep.ndx NDX_GROUP=Water_and_ions

HBA_N: GROUP NDX_FILE=bpep.ndx NDX_GROUP=HBA_N

HBA_O: GROUP NDX_FILE=bpep.ndx NDX_GROUP=HBA_O

charged: GROUP NDX_FILE=bpep.ndx NDX_GROUP=charged_groups

#GYRATION_CVs

exc_lipid_gyr: GYRATION TYPE=RADIUS ATOMS=lipid_core

exc_pep_gyr: GYRATION TYPE=RADIUS ATOMS=backbone_polar

#CONTACT CVs

HB:  COORDINATION GROUPA=HBA_N GROUPB=HBA_O SWITCH={TANH R_0=0.4 D_MAX=0.6} NLIST NL_CUTOFF=0.6 NL_STRIDE=50

c1: COORDINATION GROUPA=charged SWITCH={TANH R_0=0.4 D_MAX=0.6} NLIST NL_CUTOFF=0.6 NL_STRIDE=50

charged_c:  COMBINE ARG=c1 COEFFICIENTS=2 PERIODIC=NO

c2: COORDINATION GROUPA=lipid_core SWITCH={TANH R_0=0.4 D_MAX=0.8} NLIST NL_CUTOFF=0.8 NL_STRIDE=50

lipid_c:  COMBINE ARG=c2 COEFFICIENTS=2 PERIODIC=NO

pep_water_contact:  COORDINATION GROUPA=backbone_polar GROUPB=water SWITCH={TANH R_0=0.3 D_MAX=0.4} NLIST NL_CUTOFF=0.4 NL_STRIDE=50

lipid_water_contact:  COORDINATION GROUPA=lipid_core GROUPB=water SWITCH={TANH R_0=0.45 D_MAX=0.5} NLIST NL_CUTOFF=0.5 NL_STRIDE=50

charged_water_contact: COORDINATION GROUPA=charged GROUPB=water SWITCH={TANH R_0=0.3 D_MAX=0.4} NLIST NL_CUTOFF=0.4 NL_STRIDE=50

METAD ...

LABEL=metad-3CVs

ARG=HB,charged_c,lipid_c

PACE=5000

HEIGHT=0.02

SIGMA=0.01

BIASFACTOR=1000

ADAPTIVE=GEOM

FILE=HILLS_CONTACT

TEMP=310.0

GRID_MAX=100,350,900

GRID_MIN=0,0,0

GRID_BIN=50,200,200

WALKERS_MPI

WALKERS_DIR=../

WALKERS_RSTRIDE=1000

... METAD

PRINT  ARG=exc_lipid_gyr,exc_pep_gyr,HB,charged_c,lipid_c,lipid_water_contact,charged_water_contact,metad-3CVs.bias FILE=../COLVAR-GYR-CONTACT.dat STRIDE=500

With this gmx command:

mpirun -np 12 gmx_mpi mdrun -v -deffnm meta-contact -plumed plumed-GYR-CONTACT.dat -multidir WALKER0 WALKER1 WALKER2 -cpi

____________________________________________

When using the file system to communicate between files, this should be used:

#GROUP DEFINITIONS

backbone_polar: GROUP NDX_FILE=bpep.ndx NDX_GROUP=backbone_polar

lipid_core: GROUP NDX_FILE=bpep.ndx NDX_GROUP=lipid_core

water: GROUP NDX_FILE=bpep.ndx NDX_GROUP=Water_and_ions

HBA_N: GROUP NDX_FILE=bpep.ndx NDX_GROUP=HBA_N

HBA_O: GROUP NDX_FILE=bpep.ndx NDX_GROUP=HBA_O

charged: GROUP NDX_FILE=bpep.ndx NDX_GROUP=charged_groups

#GYRATION_CVs

exc_lipid_gyr: GYRATION TYPE=RADIUS ATOMS=lipid_core

exc_pep_gyr: GYRATION TYPE=RADIUS ATOMS=backbone_polar

#CONTACT CVs

HB:  COORDINATION GROUPA=HBA_N GROUPB=HBA_O SWITCH={TANH R_0=0.4 D_MAX=0.6} NLIST NL_CUTOFF=0.6 NL_STRIDE=50

c1: COORDINATION GROUPA=charged SWITCH={TANH R_0=0.4 D_MAX=0.6} NLIST NL_CUTOFF=0.6 NL_STRIDE=50

charged_c:  COMBINE ARG=c1 COEFFICIENTS=2 PERIODIC=NO

c2: COORDINATION GROUPA=lipid_core SWITCH={TANH R_0=0.4 D_MAX=0.8} NLIST NL_CUTOFF=0.8 NL_STRIDE=50

lipid_c:  COMBINE ARG=c2 COEFFICIENTS=2 PERIODIC=NO

pep_water_contact:  COORDINATION GROUPA=backbone_polar GROUPB=water SWITCH={TANH R_0=0.3 D_MAX=0.4} NLIST NL_CUTOFF=0.4 NL_STRIDE=50

lipid_water_contact:  COORDINATION GROUPA=lipid_core GROUPB=water SWITCH={TANH R_0=0.45 D_MAX=0.5} NLIST NL_CUTOFF=0.5 NL_STRIDE=50

charged_water_contact: COORDINATION GROUPA=charged GROUPB=water SWITCH={TANH R_0=0.3 D_MAX=0.4} NLIST NL_CUTOFF=0.4 NL_STRIDE=50

METAD ...

LABEL=metad-3CVs

ARG=HB,charged_c,lipid_c

PACE=5000

HEIGHT=0.02

SIGMA=0.01

BIASFACTOR=1000

ADAPTIVE=GEOM

FILE=../HILLS_CONTACT

TEMP=310.0

GRID_MAX=100,350,900

GRID_MIN=0,0,0

GRID_BIN=50,200,200

WALKERS_ID=[insert a number from 0 to 2 here for three WALKERS]

WALKERS_N=3 [change this number to number of WALKERS, indexed at 1]

WALKERS_DIR=../

WALKERS_RSTRIDE=1000

... METAD

PRINT  ARG=exc_lipid_gyr,exc_pep_gyr,HB,charged_c,lipid_c,lipid_water_contact,charged_water_contact,metad-3CVs.bias FILE=../COLVAR-GYR-CONTACT.dat STRIDE=500

With this gmx command:

mpirun -np 12 gmx_mpi mdrun -v -deffnm meta-contact -plumed plumed-GYR-CONTACT.dat -multidir WALKER0 WALKER1 WALKER2 -cpi

Best WIshes,
Billy

[rose rahmani]

Hi Billy,

I did multiple walker metadynamics and I need to get the FES from it. I have some questions about that. How do you beave with different HILLS? Do you concatenate them somehow(and what about their timing) and use sum_hills?

Best regards,

Roja