Equilibrium water uptake using MD/GCMC

Hello everyone,

I am not sure if this is the right place to post, maybe Science Talk is more appropriate.

I’m trying to predict the equilibrium water uptake of different ion exchange membranes using LAMMPS. These membranes are not rigid, they generally swell according to the generic law \rho = \phi_{\text{membrane}}\rho_{\text{membrane}} + \phi_{\text{water}}\rho_{\text{water}}, with \phi volume fractions.

I already followed the approach of this article [Prediction of equilibrium water uptake and ions diffusivities in ion-exchange membranes combining molecular dynamics and analytical models] which is purely MD but I wonder if the MD/GCMC approach would not be more elegant and rigorous.

My idea is to:

  1. Use fix/widom to obtain the \mu_{ex} of the water model I’m using.
  2. From a weakly hydrated membrane, use fix/gcmc to add water molecules to the system.
  3. Perform NPT integration to allow the membrane to swell.
  4. Loop points 2 and 3 until the number of water molecules in the system converges.

I’ve determined mu correctly, at least I think as I have double check that I obtain correct density for pure water, but I’m having trouble afterwards. I don’t see convergence but a constant increase in the number of water molecules in the system, so I am surely missing a point.

Could fluctuations in the dimensions of the box during NPT integration lead to a dilute system as I seem to be observing?
Is my approach feasible with LAMMPS?

Yes MD/GCMC are feasible in LAMMPS.
These methods are known to be extremely slow to converge, particularly when using chemical potential mu that corresponds to high relative humidity, because it leads to a dense phase of fluid in the system. May be this is what you observe here.

The combination of fix npt with fix gcmc is certainly a problematic one since the MC insertions can lead to locally very high potential energy and thus high pressure that will lead to (temporary) expansion of the system. Now with dense systems, expanding a system is much easier than shrinking it and thus you would probably better off using a different workflow, which to use fix gcmc in combination with fix nvt for a while, then run only with fix nvt to make certain the pressure has equilibrated and then run with fix npt for a bit to relax the volume of the system before you start over. In a way fix gcmc and fix npt work against each other, the former will adjust the density by inserting and removing atoms while the latter will adjust it by expansion or shrinking.