MD simulation of hydroxyapatite [ 110 ] surface slab using GROMACS

GROMACS version: 5.1.1
GROMACS modification: Yes

Dear Users,
I’m trying to analyze the adsorption desorption affinity of Urea to hydroxyapatite [110] surface. Hydroxyapatite contains Ca, P, O (bonded to P), O (bonded to H), H. All simulations are conducted using CHARMM27 forcefield. Nose-Hoover method is employed for temperature coupling and Parinello-Rahman for pressure coupling. PME for long-range electrostatics with a cutoff distance of 14 A. Van der Waals interactions with a cutoff of 12 A with a switch radius (Verlet scheme was not applicable for this instance and i had to rely on the group cutoff). The parameters for hydroxyapatite was added to the existing forcefield files. I have two questions regarding the system.

  1. The unit cell parameters for hydroxyapatite (hexagonal P63 /m) was a=b=9.424 A c=6.879 A. The [ 110 ] surface dimensions are supposed to be a = 6.879 A, b=16.323 A, c=9.922 A. I employed the Vesta slab cutting method and also cif2cell programs --supercell tag, with a transformation matrix of [ [0,0,1],[ -1,1,0],[1,1,0]] to create the surface. But unfortunately the c vector dimension i am getting is 9.424 A (not 9.922 A), a and b dimensions are fine. Can anybody please point out the mistake i might be doing as not to be able to get the correct c vector?

  2. Although i included the parameters for Hydroxyapatite in the forcefield files, with atoms renamed as to prevent confusion with the existing ones, i am unable to produce the coordinate file using pdb2gmx using a .pdb file (the program deletes atoms as duplicates and produces a gro file with lesser number of atoms). I produced a gro file (renaming atoms) using a python script but i am still getting the message “missing atoms” in the pdb2gmx (i think because i added HAP as a residue in the aminoacids.rtp file and the program mistakes it for an aminoacid). Further on the energy minimization does not continue and i am getting the error “Segementation fault. Core dumped” at the NVT step. How can i do a proper simulation of a crystal structure using gromacs?

Can anybody please clarify my issues? any step-to-step guidance is highly appreciated. Thank you in advance