Simulations using OPEP

Why you might want to use OPEP

OPEP is coarse grained force field providing a potential for proteins and RNA. It can speed up calculations significantly, particularly for large systems. However, there are some drawbacks. Firstly, there is no documentation for the OPEP force field (other than publications). Secondly, we need to create the input files through the web interface here.

OPEP input files

There are six files you need from the OPEP setup.

conf_initiale.pdb contains the starting coordinates (at least for GMIN). It is a normal .pdb file, but for the fact that there is only one bead for the entire side chain. Make sure the file has the END line included, otherwise the input reading may fail.

Protein nbevq
ATOM      1  NT  LYS     1      -3.753  -6.023   1.545
ATOM      2  CA  LYS     1      -3.633  -4.553   1.615
ATOM      3  LYS LYS     1      -3.283  -3.953   4.655
ATOM      4  C   LYS     1      -2.353  -4.113   0.925
ATOM      5  O   LYS     1      -1.593  -4.933   0.415
ATOM      6  N   PHE     2      -2.103  -2.793   0.905
ATOM      7  HN  PHE     2      -2.753  -2.163   1.335
ATOM      8  CA  PHE     2      -0.923  -2.243   0.275
ATOM      9  PHE PHE     2      -0.383  -3.513  -1.945
ATOM     10  C   PHE     2      -0.893  -0.733   0.395
ATOM     11  O   PHE     2      -1.793  -0.133   0.965
ATOM     12  N   PHE     3       0.157  -0.113  -0.145
ATOM     13  HN  PHE     3       0.877  -0.653  -0.605
ATOM     14  CA  PHE     3       0.307   1.327  -0.095
ATOM     15  PHE PHE     3      -0.463   2.267   2.215
ATOM     16  C   PHE     3       1.587   1.777  -0.795
ATOM     17  O   PHE     3       2.337   0.947  -1.305
ATOM     18  N   GLU     4       1.827   3.087  -0.805
ATOM     19  HN  GLU     4       1.177   3.727  -0.375
ATOM     20  CA  GLU     4       3.017   3.637  -1.435
ATOM     21  GLU GLU     4       3.477   2.437  -3.885
ATOM     22  C   GLU     4       3.047   5.157  -1.315
ATOM     23  OT  GLU     4       2.137   5.757  -0.745
ATOM     24  OXT GLU     4       3.977   5.797  -1.795
END

The cutoffs.dat, ichain.dat and scale.dat as well as parametres.list and parametres.top are files containing the setup and details of the force field parameters.

ichain.dat contains all fragments of the system.

1
1 24

scale.dat contains the scaling information for the potential.

cutoffs.dat contains the cut offs for the potential.

 1840.00d0       1600.0d0
 258.75d0         225.0d0
 74.60             64.0d0
 258.75d0         225.0d0

parametres.top is an AMBER style topology file. parametres.list contains details for the side chain beads.

 LYS      8      8
 PHE      6      6
 PHE      6      6
 GLU     17     17
   3 LYS   1   9 PHE   2  -0.010        6.473       31
   3 LYS   1  15 PHE   3  -0.110        6.473        5
   3 LYS   1  21 GLU   4   0.870        6.111      211
   9 PHE   2  15 PHE   3  -0.010        6.660       31
   9 PHE   2  21 GLU   4  -0.340        6.298       86
  15 PHE   3  21 GLU   4  -0.010        6.298       31
-999
   0
           4

Additional parameters can be set through OPEP_params.

use_qbug                              #turn on debug options for the force field
Potential_Scaling_Factor 1.2          #scaling of the force
Ion_Pair_Potential                    #use ion pair control
Ion_Pair_Scaling  1.2                 #scaling for ion pairs
Periodic_Boundary_Condition 100.0     #turning PBC on with box length 100.0
PDB_center_of_mass                    #centre pdb files for PBC
RANDOM_SEED      89302123             #integer for random seed generator
usextc                                #use xtc format, not pdb

Running GMIN

All that is needed for GMIN is the above files, a OPEP executable, and a data file, e.g.:

  SLOPPYCONV 1.0D-3
  TIGHTCONV  1.0D-6
  TRACKDATA
  CENTRE
  MAXERISE   1.0D-2
  MAXIT      50000 50000
  UPDATES    100
  MAXBFGS    0.1
  SAVE       500
  STEPS      250000  1.0
  STEP       0.1D
  RADIUS     100
  EDIFF 0.02
  DUMPSTRUCTURES
  RANSEED 20301057
  TEMPERATURE 1.3
  DUMPINT 10000
  GROUPROTATION 2
  OPEP PROTEIN

If you use DUMPSTRUCTURES, GMIN will save the minima as .pdb and start files. FEBH, BHPT, GROUPROTATIONS and GENRIGID should all work with OPEP.

Running OPTIM

In addition to the files outlined above and odata, OPTIM needs a start file (and potentially a finish file for DNEB) and a perm.allow file. The start and finish files are just the xyz coordinates. The perm.allow file simply contains a 0, i.e. there are no atoms to swap.

OPEP interface code

This section gives a brief overview over the changes/additions to the code for the OPEP interface.

The key file is opep_interface.F90. It contains the routine to get the number of atoms (OPEP_GET_ATOMS), which is called from getparams (OPTIM) and countatoms (GMIN), and the initialisation (OPEP_INIT). The initialisation calls DEFINITIONS(), which is in the file read_parameters.f90 that reads the OPEP_params file. We do not call the read routines for MD or ART parameters. md_initialise.F90 is called through INITIALISE. Here the potential is initialised, using protein-2006.f and RNS-2006.f. This is were all the input files are read. In md_initialise.F90 we commented out the temperature, restart and thermostat initialisation, as we do not need these. Once these routines are finished, we have internal variables for the atom masses and coordinates, which are passed on to the main body of the program.

After this we use calls to OPEP_ENERGY_AND_GRADIENT and OPEP_NUM_HESS from potential to get all information we need to run the programs. All output is handled as per usual and before we finish the programs we call OPEP_FINISH to deallocate all the memory.

Generating HiRE-RNA files

HiRE-v3

From full-atom to coarse-grained

use the script in the svn, under scripts/OPEP/FA2CG:

python FA2CG.py full_atom.pdb > coarse_grained.pdb

HiRE-v4

Generate parameters

Download and build the code at https://github.com/tc427/HiRE_generator

run the HiRE_parm executable on your pdb file.