Problems with positional restraints

GROMACS version: 2021.3-MODIFIED
GROMACS modification: Yes
Here post your question

Hello,

I am currently trying to simulate water at a metal interface and to do this, I have considered using freeze groups, but opted out of it due to concerns over energy exclusions so I settled on using strong positional restraints on the metal atoms. The problem is that during the NPT equilibration, with pbc, a vacuum is made and there is a large volume both above and below the solvent. I have been trying ways to fix it, but nothing seems to work? Are there any suggestions? Attached below is the mdp file:


;PREPROCESSING
define                  = -DPOSRES
;
;RUN CONTROLS
integrator              = md            ; USE MD
dt                      = 0.001         ; TIME STEP = 0.002 ps , 2 fs
nsteps                  = 4000000       ; TIME = dt * nsteps = 
;
;OUTPUT CONTROLS
nstlog                  = 1000
nstxout                 = 100000
nstvout                 = 100000
nstfout                 = 1000
nstcalcenergy           = 1000
nstenergy               = 1000
;
; NEIGHBOR SEARCHING
cutoff-scheme           = Verlet        ; NEIGHBOR SEARCHING TYPE
nstlist                 = 20
rlist                   = 1.2 
;
; COUMLOMBIC ELECTROSTATICS
coulombtype             = pme           ; COULOMB TYPE
rcoulomb                = 1.2 
;
; VAN DER WAALS FORCES
vdwtype                 = Cut-off       ; VDW TYPE
vdw-modifier            = Force-switch
rvdw_switch             = 1.0 
rvdw                    = 1.2 
;
; TEMPERATURE COUPLING 
tcoupl                  = V-rescale     ; TEMP TYPE
tc_grps                 = system
tau_t                   = 1.0           ; TEMP TIME CONSTANT
ref_t                   = 300           ; TEMP REFERENCE
;
; PRESSURE COUPLING
pcoupl                  = C-rescale     ; P-SCHEME (BREN) = LESS ACC, FAST (Par-Nell) = MORE ACC, SLOW
pcoupltype              = isotropic     ; P-TYPE
tau_p                   = 5.0           ; P-TIME CONSTANT
compressibility         = 4.5e-5        ; ONLY ALLOWED TO FLUCTUATE IN THE Z DIRECTION
ref_p                   = 1.0
;
nstcomm                 = 100
comm_mode               = linear
;comm-grps              = SOL IPT
; VELOCITY GENERATION
gen-vel                 = yes
gen-temp                = 300           ; 200K T then work up to 300K
gen-seed                = -1            ; RANDOM SEED FOR GENERATOR
;
refcoord-scaling        = all
; PERIODIC BOUNDARY CONDITIONS
pbc                     = xyz           ; PRETTY MUCH ALWAYS USE PBC