SPICA: Surface Property fItting Coarse grAined model

TutorialpSPICA

pSPICA

pSPICA is a CG FF designed with a polar CG water model. The CG modeling method of pSPICA is the same as the SPICA FF. The LJ12-5 function is adopted to reproduce the experimental density, surface tension, and dielectric permittivity of water. A brief summary of the model is given here, and the publication with full details here .
Note that to use the LJ12-5 function form, LAMMPS patch 3 AUG 2022 or later version is required.

Example for lipid membrane system:

Here we will show you how to use the pSPICA to simulate a lipid membrane system composed of 128 DMPC lipids and 1152 polar CG water particles with LAMMPS, and start the simulation from an initial configuration obtained through the SPICA tutorial for a lipid membrane system.
The equilibrated configuration, spica_conf.pdb, has the standard SPICA CG water, representing three water molecules as a single CG particle. So, we replace the SPICA CG water in the configuration with the pSPICA polar CG water using the wat2polar command of spica-tools. The command line will be:

$ cg_spica wat2polar spica_conf.pdb pspica_conf.pdb

where the first and second argument following wat2polar are input (using the SPICA water) and output (using the pSPICA water, pspica_conf.pdb) configuration files, respectively.
Next, to generate LAMMPS input files, we run the setup_lmp command for the initial configuration. The script requires parameter and topology files of the pSPICA force field. For this system, we prepared PDMPC.top and PWAT.top as the topology files, and pSPICA_lipid.prm as the parameter file of pSPICA lipids and water. We put the prefix of pSPICA topology files, "P", to distinguish those from SPICA ones. Please be careful the magnitude of partial charge in pSPICA is different from that in SPICA. For example, the CG segment PH, corresponding to phosphate group, has a partial charge of -0.5990 in pSPICA topology files. On the other hand, the segment has a partial charge of +0.1118 in SPICA topology files. These values are taken to consider the relative permittivity implicitly (pSPICA:εr = 3.2 and SPICA:εr = 80).
setup_lmp will provide the LAMMPS readable files from topology, parameter, and configuration files. The number of molecules contained in the target system must follow topology file names in the command. In addition, the topology file names must be in sequence indentical to the system configuration file.

$ cg_spica setup_lmp PDMPC.top 128 PWAT.top 1152 pSPICA_lipid.prm pspica_conf.pdb

Check the generated two LAMMPS input files, DATA.FILE and PARM.FILE. The remaining file to perform the simulation with pSPICA is a LAMMPS setup file.
We use in.npt as the setup file to optimize the configuration and simulate the system with the NPT ensemble. Please make sure that the file includes a mandatory command, fix shake. pSPICA CG water represented by two particles having oposite charge each other has a rigid bond to connect the particles.
The two commands in the in.npt should look like:

...
fix Fix0 all shake 1.0e-6 25 0 b 9
...

The b of the "fix shake" command lists the bond types that will be constrained. We have to be careful about the last argument of the command because it corresponds to the bond index that is fixed. The bond index and parameter in the PARM.FILE look like:

...
bond_coeff 7 8.4000 3.4800 # GL EST2
bond_coeff 8 8.4000 3.4800 # EST2 CM
bond_coeff 9 8.4000 3.4800 # WH WO

angle_coeff 1 sdk 3.1000 112.0000 lj9_6 0.5550 4.6200 # NC PH GL
...

We can see that the index of bond between WH-WO, composing the CG water, is 9 there. This is why the last argument of "fix shake" is set to 9 in the in.npt. Whenever you build new system configuration files and prepare input files to carry out simulations with pSPICA, you have to check the bond index of WO-WH (and CG ions), and change the last argument of "fix shake" to the bond index number. In addition, please be careful to write the "fix shake" command between the "minimize" and "fix npt" command. Otherwise, some errors in shake determinant or energy calculation will occur in the begining of the simulations.
The LAMMPS simulation can now be run in the directory where the DATA.FILE, PARM.FILE, and the in.npt are, with the following command:

$ lmp_mpi -in in.npt -log output_npt.log

Example for lipid-protein system:

We can also simulate lipid-protein systems using the pSPICA FF. After gaining the initial configuration of a target lipid-protein system by following tutorial for protein systems, we will apply the "wat2polar" command to replace the uncharged water of the SPICA FF in the initial configuration, system.cg.pdb, with the polar water of pSPICA:

$ cg_spica wat2polar system.cg.pdb p_system.cg.pdb

Then we will excute "setup_lmp" using topology files for the pSPICA FF. Although topology files for protein can be generated by following tutorial for protein systems, one must be careful about values of the charge in the pSPICA FF. The ENM command used in the SPICA protein tutorial has an option, -pspica, to assign partial charges for pSPICA. So with protein.cg.pdb, you can generate the protein topology file where values of partial charges are 0.5990 using the following command:

$ cg_spica ENM protein.cg.pdb p_protein.cg.top -pspica

p_protein.cg.top is a protein topology file for the pSPICA FF.
Then we will use "setup_lmp" command to generate LAMMPS input files, similarly to tutorial for protein systems. Download the pSPICA protein parameter set (pSPICA_protein.prm) and topology files of DOPC(PDOPC.top) and CLA (PCLA.top).

$ cg_spica setup_lmp p_protein.cg.top 1 PDOPC.top 128 PWAT.top 2134 PCLA.top 4 pSPICA_protein.prm p_system.cg.pdb

A LAMMPS input in.npt for the lipid-protein system is almost same as that for the lipid membrane system. Then we can perform LAMMPS simulations of the lipid-protein system using DATA.FILE, PARM.FILE, and the LAMMPS setup file in.npt. Make sure the command "fix shake" are properly set in the input, mentioned above.

SPICA Force Field

Okayama University