query : metadynamics NPT simulation

[MamtaMohan]

Dear Colvar team,

I would like to inquire about Metadynamics simulation for a system, subjected to NPT set-up.

With Plumed there is an issue about handling virial contribution for NAMD (as per Plumed notes and verification from authors), Metadynamics simulation for NPT can have artefacts.

I was wondering if Metadynamics simulation for NPT is possible via COLVAR.

If possible is there a pressure limit?

I would appreciate your response.

Thank you. Mamta Mohan

[giacomofiorin]

Hello Mamta, NAMD includes in the virial all terms coming from the Colvars forces. You should be OK with NPT.

[MamtaMohan]

Thank you.

---

From: Giacomo Fiorin notifications@github.com Sent: Thursday, June 8, 2017 11:44:42 AM To: colvars/colvars Cc: Mohan, Mamta; Author Subject: Re: [colvars/colvars] query : metadynamics NPT simulation (#119)

Hello Mamta, NAMD includes in the virial all terms coming from the Colvars forces. You should be OK with NPT.

— You are receiving this because you authored the thread. Reply to this email directly, view it on GitHubhttps://github.com/colvars/colvars/issues/119#issuecomment-307143705, or mute the threadhttps://github.com/notifications/unsubscribe-auth/APcc_Vdw1Do29oS4Sf_lcye01FbKNd_Jks5sCBbqgaJpZM4N0P4G.

[GiovanniBussi]

Hi @giacomofiorin!

I was wondering if this means that also PLUMED would be working correctly in NPT simulations with NAMD.

1. In PLUMED you can use VOLUME as a CV. I guess there is not way to get the correct virial contribution in this case without passing it explicitly, so this will certainly not work.

2. Even in cases where your CVs depend only on atomic positions, if one uses PBCs within PLUMED (and perhaps COLVARS) to compute distances I don't think NAMD can reconstruct the correct virial.

Notice that the usual formula (r dot f) only works correctly if you do not apply PBCs in the calculation of forces. In the general case, you should use something like r_ij dot f_ij. Thus, in principle if you have a molecule that is broken through PBCs (as GROMACS does, but perhaps NAMD doesn't) or if you use e.g. water molecules in your CVs the virial estimated just from the forces might be incorrect.

Is there something that I am missing?

Thanks!

Giovanni

[giacomofiorin]

Hi @GiovanniBussi, in NAMD typically you would disable PBC wrapping to allow for easy calculations of the centers-of-mass, so that it's not an issue there.

We definitely agree that the biggest problem is when you have an explicit dependence on the unit cell parameters, like in Plumed's VOLUME, which has no equivalent in Colvars. As for the other cases, each would require ad-hoc corrections. But in the end, even though the Colvars forces are accumulated on a per-atom basis, the coordinates of each atom are unwrapped on the fly in order to compute the virial:

        Position x_orig = atom.position;
        Transform trans = atom.transform;
        Position x_atom = lattice.reverse_transform(x_orig,trans);
        extForce += f_atom;
        extVirial += outer(f_atom,x_atom);

which basically gives you r dot f again.

[GiovanniBussi]

Thanks for the clarification!

What I suspect is that if within PLUMED or COLVARS you use in any way the information of the cell size you can end up in spurious virial contributions in the end. Imagine that you compute coordination number between all water molecules. There is no single unwrapping that allows you to use the r dot f formula, since the same water molecule would need to be in different boxes depending on which molecule it is interacting with.

What do you think?

Anyway, perhaps this is paranoia. Only people doing material science applications need this kind of stuff, and they would probably not use NAMD anyway. For solvated biomolecules I think that an error of this kind would probably lead to negligible effects.

Thanks!

Giovanni

[giacomofiorin]

I think that a coordination number where the individual pair contributions are long ranged (enough to run into the wrapping problem that you mention) also carries a massive communication overhead. You could compute it on node 0 or in parallel, but either way the simulations would be also pretty slow in addition to carrying errors in the virial.

To my knowledge, PMF calculations with a collective variable that depends on distances instead of positions (i.e. coordination number vs. total dipole moment) are mostly done by people who have worked with Parrinello's lab. So I'm fairly you are taking care of supporting those. Others may use cheaper collective variables, or entirely different approaches.