Gromacs simulation class

class gromacs_py.gromacs.gmx5.GmxSys(name=None, coor_file=None, top_file=None, tpr=None)

Gromacs system encapsulation class.

This class can be used to launch most of gromacs commands (pdb2gmx, grompp, mdrun, trjconv, editconf, genconf, …). After each steps, outputs file paths of gromacs commands are store in the class variable, like corr_file, top_file, tpr … Most function will need the corr_file and/or top_file variable to be defined.

The GmxSys object can be considered as a md simulation system. Each operation on the object will affect the object variables.

The variables nt, ntmpi and gpu_id are only used by functions which run simulations run_simulation() like em() or production().

Parameters:
  • name (str) – generic name of the system
  • sim_name (str, optional) – name of the simulation (used to create .tpr .log .edr …)
  • coor_file (str, optional) – path of the coordinate file (.pdb, .gro)
  • top_file (str, optional) – path of the .top file
  • tpr (str, optional) – path of the .tpr file
  • mdp (str, optional) – path of the .mdp file
  • xtc (str, optional) – path of the .xtc file
  • edr (str, optional) – path of the .edr file
  • log (str, optional) – path of the .log file
  • ndx (str, optional) – path of the .ndx file
  • nt (int, default=0) – Total number of threads to start
  • ntmpi (int, default=0) – Number of thread-MPI threads to start
  • gpu_id (str, default=None) – List of GPU device id-s to use, specifies the per-node PP rank to GPU mapping
  • sys_history (list of GmxSys()) – List of previous GmxSys() states
Example: 
>>> TEST_OUT = str(getfixture('tmpdir'))
>>> prot = GmxSys(name='1y0m', coor_file=TEST_PATH+'/1y0m.pdb')
>>> ###################################
>>> ####   Create the topologie:   ###
>>> ###################################
>>> prot.prepare_top(out_folder=os.path.join(TEST_OUT, 'top_SH3')) #doctest: +ELLIPSIS
Succeed to read file .../test/input/1y0m.pdb ,  648 atoms found
Succeed to save file tmp_pdb2pqr.pdb
pdb2pqr.py --ff CHARMM --ffout CHARMM --chain tmp_pdb2pqr.pdb 00_1y0m.pqr
Succeed to read file 00_1y0m.pqr ,  996 atoms found
Chain: A  Residue: 0 to 60
Succeed to save file 01_1y0m_good_his.pdb
-Create topologie
gmx pdb2gmx -f 01_1y0m_good_his.pdb -o 1y0m_pdb2gmx.pdb -p 1y0m_pdb2gmx.top -i 1y0m_posre.itp -water tip3p -ff charmm36-jul2017 -ignh -vsite hydrogens
Molecule topologie present in 1y0m_pdb2gmx.top , extract the topologie in a separate file: 1y0m_pdb2gmx.itp
Protein_chain_A
-ITP file: 1y0m_pdb2gmx.itp
-molecules defined in the itp file:
* Protein_chain_A
Rewrite topologie: 1y0m_pdb2gmx.top
>>> ###################################
>>> ####    Add water and ions:     ###
>>> ###################################
>>> prot.solvate_add_ions(out_folder=os.path.join(TEST_OUT, 'top_sys')) #doctest: +ELLIPSIS
-Create pbc box
gmx editconf -f .../top_SH3/1y0m_pdb2gmx.pdb -o .../top_SH3/1y0m_pdb2gmx_box.pdb -bt dodecahedron -d 1.1
-Solvate the pbc box
Copy topologie file and dependancies
Copy topologie file and dependancies
-Create the tpr file  genion_1y0m_water_ion.tpr
gmx grompp -f .../gromacs/template/mini.mdp -c 1y0m_water.pdb -r 1y0m_water.pdb -p 1y0m_water_ion.top -po out_mini.mdp -o genion_1y0m_water_ion.tpr -maxwarn 1
-Add ions to the system with an ionic concentration of 0.15 M , sytem charge = 0.0 water num= 4775
Add ions : NA : 12   CL : 12
gmx genion -s genion_1y0m_water_ion.tpr -p 1y0m_water_ion.top -o 1y0m_water_ion.gro -np 12 -pname NA -nn 12 -nname CL
>>> ###################################
>>> ####    Minimize the system     ###
>>> ###################################
>>> prot.em(out_folder=os.path.join(TEST_OUT, 'em_SH3'), nsteps=100, constraints='none')
-Create the tpr file  1y0m.tpr
gmx grompp -f 1y0m.mdp -c ../top_sys/1y0m_water_ion.gro -r ../top_sys/1y0m_water_ion.gro -p ../top_sys/1y0m_water_ion.top -po out_1y0m.mdp -o 1y0m.tpr -maxwarn 1
-Launch the simulation 1y0m.tpr
gmx mdrun -s 1y0m.tpr -deffnm 1y0m -nt 0 -ntmpi 0 -nsteps -2 -nocopyright
>>> ###################################
>>> ####    Create a D peptide      ###
>>> ###################################
>>> pep = GmxSys(name='D')
>>> pep.create_peptide(sequence='D', out_folder=os.path.join(TEST_OUT, 'top_D'), em_nsteps=100, equi_nsteps=0) #doctest: +ELLIPSIS
-Make peptide: D
residue name:X
residue name:D
Succeed to save file .../top_D/D.pdb
-Create topologie
gmx pdb2gmx -f ../D.pdb -o D_pdb2gmx.pdb -p D_pdb2gmx.top -i D_posre.itp -water tip3p -ff charmm36-jul2017 -ignh -ter -vsite hydrogens
Molecule topologie present in D_pdb2gmx.top , extract the topologie in a separate file: D_pdb2gmx.itp
Protein_chain_P
-ITP file: D_pdb2gmx.itp
-molecules defined in the itp file:
* Protein_chain_P
Rewrite topologie: D_pdb2gmx.top
-Create pbc box
gmx editconf -f .../top_D/00_top/D_pdb2gmx.pdb -o .../top_D/00_top/D_pdb2gmx_box.pdb -bt dodecahedron -d 1.0
-Create the tpr file  D.tpr
gmx grompp -f D.mdp -c ../00_top/D_pdb2gmx_box.pdb -r ../00_top/D_pdb2gmx_box.pdb -p ../00_top/D_pdb2gmx.top -po out_D.mdp -o D.tpr -maxwarn 1
-Launch the simulation D.tpr
gmx mdrun -s D.tpr -deffnm D -nt 0 -ntmpi 0 -nsteps -2 -nocopyright
>>> #######################################################
>>> ### Insert 4 copy of the peptide in the SH3 system: ###
>>> #######################################################
>>> prot.insert_mol_sys(mol_gromacs=pep, mol_num=4, new_name='SH3_D', out_folder=os.path.join(TEST_OUT, 'top_D_SH3')) #doctest: +ELLIPSIS
-Copy pbc box using genconf
Succeed to read file ../top_D/01_mini/D_copy_box.pdb ,  88 atoms found
Succeed to save file ../top_D/01_mini/D_copy_box.pdb
AA num: 1
-Convert trj/coor
gmx trjconv -f ../em_SH3/1y0m.gro -o ../em_SH3/1y0m_compact.pdb -s ../em_SH3/1y0m.tpr -ur compact -pbc mol
Concat files: ['../em_SH3/1y0m_compact.pdb', '../top_D/01_mini/D_copy_box.pdb']
Succeed to save concat file:  SH3_D_pre_mix.pdb
Succeed to read file SH3_D_pre_mix.pdb ,  15429 atoms found
Insert mol in system
Insert 4 mol of 2 residues each
insert mol   1, water mol   ..., time=0...
insert mol   2, water mol   ..., time=0...
insert mol   3, water mol   ..., time=0...
insert mol   4, water mol   ..., time=0...
Delete ... overlapping water atoms
Succeed to save file SH3_D.pdb
Peptide
Add 4 mol D_pdb2gmx.itp
Succeed to read file SH3_D.pdb ,  15... atoms found
Water num: 47...
CHARGE: -4.0
Should neutralize the system
Copy topologie file and dependancies
-Create the tpr file  genion_SH3_D_neutral.tpr
gmx grompp -f .../template/mini.mdp -c SH3_D.pdb -r SH3_D.pdb -p SH3_D_neutral.top -po out_mini.mdp -o genion_SH3_D_neutral.tpr -maxwarn 1
-Add ions to the system with an ionic concentration of 0 M , sytem charge = -4.0 water num= 47...
Add ions : NA : 4   CL : 0
gmx genion -s genion_SH3_D_neutral.tpr -p SH3_D_neutral.top -o SH3_D_neutral.gro -np 4 -pname NA -nn 0 -nname CL
>>> ################################
>>> ####   Minimize the system   ###
>>> ################################
>>> prot.em_2_steps(out_folder=os.path.join(TEST_OUT, 'top_D_SH3'), no_constr_nsteps=100, constr_nsteps=100)
-Create the tpr file  Init_em_1y0m.tpr
gmx grompp -f Init_em_1y0m.mdp -c SH3_D_neutral.gro -r SH3_D_neutral.gro -p SH3_D_neutral.top -po out_Init_em_1y0m.mdp -o Init_em_1y0m.tpr -maxwarn 1
-Launch the simulation Init_em_1y0m.tpr
gmx mdrun -s Init_em_1y0m.tpr -deffnm Init_em_1y0m -nt 0 -ntmpi 0 -nsteps -2 -nocopyright
-Create the tpr file  1y0m.tpr
gmx grompp -f 1y0m.mdp -c Init_em_1y0m.gro -r Init_em_1y0m.gro -p SH3_D_neutral.top -po out_1y0m.mdp -o 1y0m.tpr -maxwarn 1
-Launch the simulation 1y0m.tpr
gmx mdrun -s 1y0m.tpr -deffnm 1y0m -nt 0 -ntmpi 0 -nsteps -2 -nocopyright
>>> ##################################
>>> ####    Show system history    ###
>>> ##################################
>>> prot.display_history() #doctest: +ELLIPSIS
State -3:
<BLANKLINE>
name         : 1y0m
sim_name     : genion_1y0m_water_ion
coor_file    : .../top_sys/1y0m_water_ion.gro
top_file     : .../top_sys/1y0m_water_ion.top
tpr          : .../top_sys/genion_1y0m_water_ion.tpr
mdp          : gromacs_py/gromacs/template/mini.mdp
nt           : 0
ntmpi        : 0
sys_history  : 0
<BLANKLINE>
State -2:
<BLANKLINE>
name         : 1y0m
sim_name     : genion_SH3_D_neutral
coor_file    : .../top_D_SH3/SH3_D_neutral.gro
top_file     : .../top_D_SH3/SH3_D_neutral.top
tpr          : .../top_D_SH3/genion_SH3_D_neutral.tpr
mdp          : gromacs_py/gromacs/template/mini.mdp
xtc          : .../em_SH3/1y0m.trr
edr          : .../em_SH3/1y0m.edr
log          : .../em_SH3/1y0m.log
nt           : 0
ntmpi        : 0
sys_history  : 0
<BLANKLINE>
State -1:
<BLANKLINE>
name         : 1y0m
sim_name     : Init_em_1y0m
coor_file    : .../top_D_SH3/Init_em_1y0m.gro
top_file     : .../top_D_SH3/SH3_D_neutral.top
tpr          : .../top_D_SH3/Init_em_1y0m.tpr
mdp          : .../top_D_SH3/Init_em_1y0m.mdp
xtc          : .../top_D_SH3/Init_em_1y0m.trr
edr          : .../top_D_SH3/Init_em_1y0m.edr
log          : .../top_D_SH3/Init_em_1y0m.log
nt           : 0
ntmpi        : 0
sys_history  : 0
<BLANKLINE>
>>> ###################################
>>> ####   Equilibrate the system   ###
>>> ###################################
>>> equi_template_mdp = os.path.join(GROMACS_MOD_DIRNAME, "template/equi.mdp")
>>> mdp_options = {'nsteps': 100, 'define': '-DPOSRES', 'dt': 0.001}
>>> prot.run_md_sim(out_folder=os.path.join(TEST_OUT, 'equi_HA_D_SH3'), name="equi_HA_D_SH3",                        mdp_template=equi_template_mdp,                        mdp_options=mdp_options)
-Create the tpr file  equi_HA_D_SH3.tpr
gmx grompp -f equi_HA_D_SH3.mdp -c ../top_D_SH3/1y0m.gro -r ../top_D_SH3/1y0m.gro -p ../top_D_SH3/SH3_D_neutral.top -po out_equi_HA_D_SH3.mdp -o equi_HA_D_SH3.tpr -maxwarn 0
-Launch the simulation equi_HA_D_SH3.tpr
gmx mdrun -s equi_HA_D_SH3.tpr -deffnm equi_HA_D_SH3 -nt 0 -ntmpi 0 -nsteps -2 -nocopyright
>>> prot.get_simulation_time() #doctest: +ELLIPSIS
-Get simulation time from : .../equi_HA_D_SH3/equi_HA_D_SH3.cpt
gmx check -f .../equi_HA_D_SH3/equi_HA_D_SH3.cpt
0.1
>>> prot.convert_trj(traj=False) #doctest: +ELLIPSIS
-Convert trj/coor
gmx trjconv -f .../equi_HA_D_SH3/equi_HA_D_SH3.gro -o .../equi_HA_D_SH3/equi_HA_D_SH3_compact.pdb -s .../equi_HA_D_SH3/equi_HA_D_SH3.tpr -ur compact -pbc mol
>>> prot.display() #doctest: +ELLIPSIS
name         : 1y0m
sim_name     : equi_HA_D_SH3
coor_file    : .../equi_HA_D_SH3/equi_HA_D_SH3_compact.pdb
top_file     : .../top_D_SH3/SH3_D_neutral.top
tpr          : .../equi_HA_D_SH3/equi_HA_D_SH3.tpr
mdp          : .../equi_HA_D_SH3/equi_HA_D_SH3.mdp
xtc          : .../equi_HA_D_SH3/equi_HA_D_SH3.xtc
edr          : .../equi_HA_D_SH3/equi_HA_D_SH3.edr
log          : .../equi_HA_D_SH3/equi_HA_D_SH3.log
nt           : 0
ntmpi        : 0
sys_history  : 4
>>> #########################################
>>> ### Extract Potential Energy and Temp ###
>>> #########################################
>>> ener_pd = prot.get_ener(['Potential', 'Temp'])  #doctest: +ELLIPSIS
-Extract energy
gmx energy -f .../equi_HA_D_SH3/equi_HA_D_SH3.edr -o tmp.xvg
>>> ener_pd['Potential'].mean() #doctest: +ELLIPSIS
-2...

Note

An history of all command used could be saved.

Note

Files necessary for testing:../test/input/1y0m.pdb, ../test/input/5vav.pdb To do the unitary test, execute gmx5.py (-v for verbose mode)

add_ions(out_folder, name=None, ion_C=0.15, pname='NA', nname='CL', solv_name='SOL', check_file_out=True)

Add ion in a system to neutralise the sys_charge and to reach the ionic concentration ion_C.

Ion number are computed using the water number and the charge of the system:

  1. cation_num = int(ion_C x water_num)/55.5
  2. if cation_num+sys_charge >= 0 then anion_num = cation_num+sys_charge else cation_num = -sys_charge
Parameters:
  • out_folder (str) – path of the output file folder
  • name (str, optional, default=None) – generic name of the system
  • ion_C (float, optional, default=0.15) – ionic concentraton (Molar)
  • pname (str, optional, default="SOL") – cation name
  • pname – anion name
  • pname – solvant name
  • check_file_out (bool, optional, default=True) – flag to check or not if file has already been created. If the file is present then the command break.

Object requirement(s):

  • self.coor_file
  • self.top_file

Object field(s) changed:

  • self.coor_file
  • self.top_file
Example: 
>>> TEST_OUT = getfixture('tmpdir')
>>> prot = GmxSys(name='1y0m', coor_file=TEST_PATH+'/1y0m.pdb')
>>> prot.add_top(out_folder=TEST_OUT+'/add_ions/top_SH3/') #doctest: +ELLIPSIS
-Create topologie
gmx pdb2gmx -f .../input/1y0m.pdb -o 1y0m_pdb2gmx.pdb -p 1y0m_pdb2gmx.top -i 1y0m_posre.itp -water tip3p -ff charmm36-jul2017
Molecule topologie present in 1y0m_pdb2gmx.top , extract the topologie in a separate file: 1y0m_pdb2gmx.itp
Protein_chain_A
-ITP file: 1y0m_pdb2gmx.itp
-molecules defined in the itp file:
* Protein_chain_A
Rewrite topologie: 1y0m_pdb2gmx.top
>>> prot.create_box() #doctest: +ELLIPSIS
-Create pbc box
gmx editconf -f .../add_ions/top_SH3/1y0m_pdb2gmx.pdb -o .../add_ions/top_SH3/1y0m_pdb2gmx_box.pdb -bt dodecahedron -d 1.0
>>> prot.solvate_box(out_folder=TEST_OUT+'/add_ions/top_SH3_water/')
-Solvate the pbc box
Copy topologie file and dependancies
>>> prot.add_ions(out_folder=TEST_OUT+'/add_ions/top_SH3_water_ions/') #doctest: +ELLIPSIS
Copy topologie file and dependancies
-Create the tpr file  genion_1y0m_ion.tpr
gmx grompp -f .../gromacs/template/mini.mdp -c ../top_SH3_water/1y0m_water.pdb -r ../top_SH3_water/1y0m_water.pdb -p 1y0m_ion.top -po out_mini.mdp -o genion_1y0m_ion.tpr -maxwarn 1
-Add ions to the system with an ionic concentration of 0.15 M , sytem charge = 0.0 water num= 56...
Add ions : NA : 15   CL : 15
gmx genion -s genion_1y0m_ion.tpr -p 1y0m_ion.top -o 1y0m_ion.gro -np 15 -pname NA -nn 15 -nname CL
>>> prot.em(out_folder=TEST_OUT+'/add_ions/em_SH3_water_ions/', nsteps=100, constraints="none")
-Create the tpr file  1y0m.tpr
gmx grompp -f 1y0m.mdp -c ../top_SH3_water_ions/1y0m_ion.gro -r ../top_SH3_water_ions/1y0m_ion.gro -p ../top_SH3_water_ions/1y0m_ion.top -po out_1y0m.mdp -o 1y0m.tpr -maxwarn 1
-Launch the simulation 1y0m.tpr
gmx mdrun -s 1y0m.tpr -deffnm 1y0m -nt 0 -ntmpi 0 -nsteps -2 -nocopyright

Note

If name is not defined, the command will create a new .pdb and .top file name after the object name and adding “_ion”.

add_mdp(mdp_template, mdp_options, folder_out='', check_file_out=True)

Create the MD simulation input mdp file.

Read a template mdp file and replace define fields in mdp_options with the new value. In case the field name has a ‘-‘ , repalce it by : ‘_’.

Parameters:
  • mdp_template (str) – mdp file template
  • mdp_options (dict) – New parameters to use
  • folder_out (str, default="") – Path for output file
  • check_file_out (bool, optional, default=True) – flag to check or not if file has already been created. If the file is present then the command break.

Object requirement(s):

  • self.sim_name

Object field(s) changed:

  • self.mdp
add_ndx(ndx_cmd_input, ndx_name=None, folder_out='', check_file_out=True)

Create a ndx file using gmx make_ndx

Parameters:
  • ndx_cmd_input (str) – Input arguments for gmx make_ndx
  • ndx_name (str, default=None) – output name for the index file
  • folder_out (str, default="") – Path for output file
  • check_file_out (bool, optional, default=True) – flag to check or not if file has already been created. If the file is present then the command break.

Object requirement(s):

  • self.coor_file

Object field(s) changed:

  • self.ndx

Note

If name is not defined, will use the name of the object.

add_top(out_folder, name=None, ff='charmm36-jul2017', water='tip3p', check_file_out=True, pdb2gmx_option_dict={}, input_pdb2gmx='', posre_post='')

Launch the pdb2gmx command.

The objet variable self.coor_file has to be defined before launching this function. pdb2gmx will create a new coordinate file name+"_pdb2gmx.pdb", a topologie name+"_pdb2gmx.top" and several molecule itp and posre files. If name is not defined, it will use the object name.

Parameters:
  • out_folder (str) – path of the output file folder
  • name (str, optional, default=None) – generic name of the system
  • ff (str, optional, default="charmm36") – forcefield
  • water (str, optional, default="tip3p") – water model
  • check_file_out (bool, optional, default=True) – flag to check or not if file has already been created. If the file is present then the command break.
  • pdb2gmx_option_dict (dict, optional, default=None) – dictionnary of option for pdb2gmx, for example if you want to ignore input hydrogens use:{'ignh': None}. The ‘-‘ before the option is to avoid.
  • input_pdb2gmx (str, optional, default=None) – input for pdb2gmx request

Object requirement(s):

  • self.coor_file

Object field(s) changed:

  • self.coor_file
  • self.top_file
Example: 
>>> TEST_OUT = getfixture('tmpdir')
>>> prot = GmxSys(name='1y0m', coor_file=TEST_PATH+'/1y0m.pdb')
>>> #Basic usage :
>>> prot.add_top(out_folder=TEST_OUT+'/add_top/top_SH3') #doctest: +ELLIPSIS
-Create topologie
gmx pdb2gmx -f .../input/1y0m.pdb -o 1y0m_pdb2gmx.pdb -p 1y0m_pdb2gmx.top -i 1y0m_posre.itp -water tip3p -ff charmm36-jul2017
Molecule topologie present in 1y0m_pdb2gmx.top , extract the topologie in a separate file: 1y0m_pdb2gmx.itp
Protein_chain_A
-ITP file: 1y0m_pdb2gmx.itp
-molecules defined in the itp file:
* Protein_chain_A
Rewrite topologie: 1y0m_pdb2gmx.top
>>> #########################################
>>> # Use of different options for pdb2gmx: #
>>> #########################################
>>> # Ignore hydrogens: 'ignh': None
>>> # Define amino acid termini: 'ter': None
>>> # Needs to answer pdb2gmx request concerning termini
>>> # with: input_pdb2gmx ="1 \n 0"
>>> prot = GmxSys(name='1y0m', coor_file=TEST_PATH+'/1y0m.pdb')
>>> prot.add_top(out_folder=TEST_OUT+'/add_top/top_SH3_2/',
...     pdb2gmx_option_dict={'ignh': None, 'ter': None},
...     input_pdb2gmx="1 \n 0") #doctest: +ELLIPSIS
-Create topologie
gmx pdb2gmx -f .../input/1y0m.pdb -o 1y0m_pdb2gmx.pdb -p 1y0m_pdb2gmx.top -i 1y0m_posre.itp -water tip3p -ff charmm36-jul2017 -ignh -ter
Molecule topologie present in 1y0m_pdb2gmx.top , extract the topologie in a separate file: 1y0m_pdb2gmx.itp
Protein_chain_A
-ITP file: 1y0m_pdb2gmx.itp
-molecules defined in the itp file:
* Protein_chain_A
Rewrite topologie: 1y0m_pdb2gmx.top

Note

To avoid conflict with focefields, the environment variable $GMXLIB is change to GROMACS_MOD_DIRNAME+"/template/" where curently only charmm36 is present, if you want to use another forcefield, copy your forcefield folder in GROMACS_MOD_DIRNAME+"/template/", or change the current code.

add_tpr(name, r=None, po=None, folder_out='', check_file_out=True, **grompp_options)

Create a tpr file using gmx grompp.

Parameters:
  • name (str) – name of the simulation
  • r (str, default=None) – reference coordinate file for position restraints
  • po (str, default=None) – output file for the mdp file
  • folder_out (str, default="") – Path for output file
  • check_file_out (bool, optional, default=True) – flag to check or not if file has already been created. If the file is present then the command break.
  • **grompp_options – Optional arguments for gmx grompp

Object requirement(s):

  • self.mdp
  • self.coor_file
  • self.top_file

Object optional requirement(s):

  • self.ndx

Object field(s) changed:

  • self.tpr
static concat_coor(*coor_in_files, pdb_out)

Concat a list of coordinates file in one coordinate file:

Parameters:
  • coor_in_files (list of str) – list of pdb/gro files
  • pdb_out (str) – file to save the concat structure
Returns:

name of the new pdb file
Return type:

str

Note

This function does not use or affect the GmxSys object.

concat_edr(*edr_files_list, concat_edr_out, check_file_out=True)

Concat a list of energy file in one energy file:

Parameters:
  • xtc_in_files (list of str) – list of edr files
  • concat_edr_out (str) – file to save the concat energy

Object field(s) changed:

  • self.edr
concat_traj(*xtc_files_list, concat_traj_out, check_file_out=True)

Concat a list of trajectory file in one trajectory file:

Parameters:
  • xtc_in_files (list of str) – list of xtc files
  • concat_traj_out (str) – file to save the concat trajectory

Object field(s) changed:

  • self.xtc
convert_trj(name=None, ur='compact', pbc='mol', select='System', traj=True, specific_coor_out=None, check_file_out=True, **cmd_args)

Convert a trajectory or coordinate file using the commande gmx trjconv.

This is specially usefull when the protein is break across pbc. Using convert_trj() with default parameters will fix it.

Parameters:
  • name (str, optional, default=None) – generic name of the system
  • ur (str, default="compact") – unit-cell representation (“rect”, “tric”, “compact”)
  • pbc (str, default="mol") – PBC treatment (“none”, “mol”, “res”, “atom”, “nojump”, “cluster”, “whole”)
  • select (str, default="System") – group for output
  • specific_coor_out (str, optional, default=None) – specific output file
  • traj (bool, default=True) – Flag to convert trajectory or coordinates
  • check_file_out (bool, optional, default=True) – flag to check or not if file has already been created. If the file is present then the command break.
  • **cmd_args – Optional arguments for gmx trjconv

Object requirement(s):

  • self.tpr
  • self.coor_file or self.xtc

Object field(s) changed:

  • self.coor_file or self.xtc
Example: 
>>> TEST_OUT = getfixture('tmpdir')
>>> prot = GmxSys(name='1y0m', coor_file=TEST_PATH+'/1y0m.pdb')
>>> prot.add_top(out_folder=TEST_OUT+'/convert_trj/top_SH3/') #doctest: +ELLIPSIS
-Create topologie
gmx pdb2gmx -f .../input/1y0m.pdb -o 1y0m_pdb2gmx.pdb -p 1y0m_pdb2gmx.top -i 1y0m_posre.itp -water tip3p -ff charmm36-jul2017
Molecule topologie present in 1y0m_pdb2gmx.top , extract the topologie in a separate file: 1y0m_pdb2gmx.itp
Protein_chain_A
-ITP file: 1y0m_pdb2gmx.itp
-molecules defined in the itp file:
* Protein_chain_A
Rewrite topologie: 1y0m_pdb2gmx.top
>>> prot.create_box() #doctest: +ELLIPSIS
-Create pbc box
gmx editconf -f .../convert_trj/top_SH3/1y0m_pdb2gmx.pdb -o .../convert_trj/top_SH3/1y0m_pdb2gmx_box.pdb -bt dodecahedron -d 1.0
>>> prot.solvate_box(out_folder=TEST_OUT+'/convert_trj/top_SH3_water/')
-Solvate the pbc box
Copy topologie file and dependancies
>>> prot.em(out_folder=TEST_OUT+'/convert_trj/em_SH3_water/', nsteps=100, constraints="none")
-Create the tpr file  1y0m.tpr
gmx grompp -f 1y0m.mdp -c ../top_SH3_water/1y0m_water.pdb -r ../top_SH3_water/1y0m_water.pdb -p ../top_SH3_water/1y0m_water.top -po out_1y0m.mdp -o 1y0m.tpr -maxwarn 1
-Launch the simulation 1y0m.tpr
gmx mdrun -s 1y0m.tpr -deffnm 1y0m -nt 0 -ntmpi 0 -nsteps -2 -nocopyright
>>> prot.convert_trj(traj=False) #doctest: +ELLIPSIS
-Convert trj/coor
gmx trjconv -f .../convert_trj/em_SH3_water/1y0m.gro -o .../convert_trj/em_SH3_water/1y0m_compact.pdb -s .../convert_trj/em_SH3_water/1y0m.tpr -ur compact -pbc mol

Note

If name is not defined, the command will create a new pdb file name after the input one and adding “_compact.pdb” or “_compact.xtc”. If name is defined the pdb filed will be saved in the same directory as input file, the “_compact.pdb” or “_compact.xtc” will be added to name.

copy_box(nbox, name=None, check_file_out=True, **cmd_args)

Copy images of a given corrdinates in x, y, and z directions using gmx genconf.

nbox needs a list of 3 string for number x,y,z dimensions copy

This is specially usefull when the protein is break across pbc. Using convert_trj() with default parameters will fix it.

Parameters:
  • nbox (list of string) – list of 3 string for number of x, y, z dimensions copy
  • name (str, optional, default=None) – generic name of the system
  • check_file_out (bool, optional, default=True) – flag to check or not if file has already been created. If the file is present then the command break.
  • **cmd_args – Optional arguments for gmx genconf

Object requirement(s):

  • self.coor_file

Object field(s) changed:

  • self.coor_file
Example: 
>>> TEST_OUT = getfixture('tmpdir')
>>> prot = GmxSys(name='1y0m', coor_file=TEST_PATH+'/1y0m.pdb')
>>> prot.add_top(out_folder=TEST_OUT+'/copy_box/top_SH3/') #doctest: +ELLIPSIS
-Create topologie
gmx pdb2gmx -f .../input/1y0m.pdb -o 1y0m_pdb2gmx.pdb -p 1y0m_pdb2gmx.top -i 1y0m_posre.itp -water tip3p -ff charmm36-jul2017
Molecule topologie present in 1y0m_pdb2gmx.top , extract the topologie in a separate file: 1y0m_pdb2gmx.itp
Protein_chain_A
-ITP file: 1y0m_pdb2gmx.itp
-molecules defined in the itp file:
* Protein_chain_A
Rewrite topologie: 1y0m_pdb2gmx.top
>>> prot.create_box() #doctest: +ELLIPSIS
-Create pbc box
gmx editconf -f .../copy_box/top_SH3/1y0m_pdb2gmx.pdb -o .../copy_box/top_SH3/1y0m_pdb2gmx_box.pdb -bt dodecahedron -d 1.0
>>> prot.copy_box(nbox=[4,1,1])
-Copy pbc box using genconf

Note

If name is not defined, the command will create a new pdb file name after the input one and adding “_copy_box.pdb”. If name is defined the pdb filed will be saved in the same directory as input file, “_copy_box.pdb” will be added to name.

create_box(name=None, dist=1.0, box_type='dodecahedron', check_file_out=True)

Create pbc box using gmx editconf

Parameters:
  • name (str, optional, default=None) – generic name of the system
  • dist (float, default=1.0) – Distance between the solute and the box (nm)
  • box_type (str, default="dodecahedron") – Box type (“triclinic”, “cubic”, “dodecahedron”, “octahedron”)
  • check_file_out (bool, optional, default=True) – flag to check or not if file has already been created. If the file is present then the command break.

Object requirement(s):

  • self.coor_file

Object field(s) changed:

  • self.coor_file
Example: 
>>> TEST_OUT = getfixture('tmpdir')
>>> prot = GmxSys(name='1y0m', coor_file=TEST_PATH+'/1y0m.pdb')
>>> prot.add_top(out_folder=TEST_OUT+'/create_box/top_SH3/') #doctest: +ELLIPSIS
-Create topologie
gmx pdb2gmx -f .../input/1y0m.pdb -o 1y0m_pdb2gmx.pdb -p 1y0m_pdb2gmx.top -i 1y0m_posre.itp -water tip3p -ff charmm36-jul2017
Molecule topologie present in 1y0m_pdb2gmx.top , extract the topologie in a separate file: 1y0m_pdb2gmx.itp
Protein_chain_A
-ITP file: 1y0m_pdb2gmx.itp
-molecules defined in the itp file:
* Protein_chain_A
Rewrite topologie: 1y0m_pdb2gmx.top
>>> prot.create_box() #doctest: +ELLIPSIS
-Create pbc box
gmx editconf -f .../create_box/top_SH3/1y0m_pdb2gmx.pdb -o .../create_box/top_SH3/1y0m_pdb2gmx_box.pdb -bt dodecahedron -d 1.0

Note

If name is not defined, the command will create a new pdb file name after the input one and adding “_box.pdb”. If name is defined the pdb filed will be saved in the same directory as input file, the “_box.pdb” will be added to name.

create_peptide(sequence, out_folder, N_ter='None', C_ter='COOH', em_nsteps=1000, equi_nsteps=10000, posre_post='_pep')

Create a linear peptide structure and topologie:

  1. Create a peptide with pymol with one more residue G at the beginning of the peptide. This residue will then be change to an ACE. NH2 terminaison raise some issue with virtual sites and cannot be used.
  2. Create the topologie using add_top()
  3. Minimise the structure using em()
  4. Do a vacuum equilibration of the peptide using run_md_sim()
Parameters:
  • out_folder (str) – path of the output file folder
  • name (str, optional, default=None) – generic name of the system
  • ion_C (float, optional, default=0.15) – ionic concentraton (Molar)

Object requirement(s):

  • None

Object field(s) changed:

  • self.coor_file
  • self.top_file
Example: 
>>> TEST_OUT = getfixture('tmpdir')
>>> pep = GmxSys(name='SAM_pep')
>>> pep.create_peptide(sequence='SAM', out_folder=os.path.join(str(TEST_OUT), 'peptide'), em_nsteps=100, equi_nsteps=100) #doctest: +ELLIPSIS
-Make peptide: SAM
residue name:X
residue name:S
residue name:A
residue name:M
Succeed to save file .../peptide/SAM.pdb
-Create topologie
gmx pdb2gmx -f ../SAM.pdb -o SAM_pdb2gmx.pdb -p SAM_pdb2gmx.top -i SAM_posre.itp -water tip3p -ff charmm36-jul2017 -ignh -ter -vsite hydrogens
Molecule topologie present in SAM_pdb2gmx.top , extract the topologie in a separate file: SAM_pdb2gmx.itp
Protein_chain_P
-ITP file: SAM_pdb2gmx.itp
-molecules defined in the itp file:
* Protein_chain_P
Rewrite topologie: SAM_pdb2gmx.top
-Create pbc box
gmx editconf -f .../peptide/00_top/SAM_pdb2gmx.pdb -o .../peptide/00_top/SAM_pdb2gmx_box.pdb -bt dodecahedron -d 1.0
-Create the tpr file  SAM_pep.tpr
gmx grompp -f SAM_pep.mdp -c ../00_top/SAM_pdb2gmx_box.pdb -r ../00_top/SAM_pdb2gmx_box.pdb -p ../00_top/SAM_pdb2gmx.top -po out_SAM_pep.mdp -o SAM_pep.tpr -maxwarn 1
-Launch the simulation SAM_pep.tpr
gmx mdrun -s SAM_pep.tpr -deffnm SAM_pep -nt 0 -ntmpi 0 -nsteps -2 -nocopyright
-Create the tpr file  equi_vacuum_SAM.tpr
gmx grompp -f equi_vacuum_SAM.mdp -c ../01_mini/SAM_pep.gro -r ../01_mini/SAM_pep.gro -p ../00_top/SAM_pdb2gmx.top -po out_equi_vacuum_SAM.mdp -o equi_vacuum_SAM.tpr -maxwarn 1
-Launch the simulation equi_vacuum_SAM.tpr
gmx mdrun -s equi_vacuum_SAM.tpr -deffnm equi_vacuum_SAM -nt 0 -ntmpi 0 -nsteps -2 -nocopyright

Note

Pymol need to be installed to run the peptide creation

cyclic_peptide_top(out_folder, name=None, check_file_out=True)

Prepare a topologie for a cyclic peptide

  1. Create a peptide topologie with NH2 Cter and COO- Nter using``add_top()``.
  2. Delete useless termini atoms.
  3. Change atom types, names and charges.
  4. Add backbone bonds, angle and dihedral parameters.
  5. Finally compute the topologie with pdb2gmx add_top()
Parameters:
  • out_folder (str) – path of the output file folder
  • name (str, optional, default=None) – generic name of the system
  • check_file_out (bool, optional, default=True) – flag to check or not if file has already been created. If the file is present then the command break.

Object requirement(s):

  • self.coor_file

Object field(s) changed:

  • self.coor_file
  • self.top_file
Example: 
>>> TEST_OUT = getfixture('tmpdir')
>>> cyclic_pep = GmxSys(name='5vav', coor_file=TEST_PATH+'/5vav.pdb')
>>>
>>> #Basic usage :
>>> cyclic_pep.cyclic_peptide_top(out_folder=os.path.join(str(TEST_OUT),'cyclic/top')) #doctest: +ELLIPSIS
-Create topologie
gmx pdb2gmx -f ...input/5vav.pdb -o no_cyclic_5vav_pdb2gmx.pdb -p no_cyclic_5vav_pdb2gmx.top -i no_cyclic_5vav_posre.itp -water tip3p -ff charmm36-jul2017 -ignh -ter -vsite no
Molecule topologie present in no_cyclic_5vav_pdb2gmx.top , extract the topologie in a separate file: no_cyclic_5vav_pdb2gmx.itp
Protein_chain_A
-ITP file: no_cyclic_5vav_pdb2gmx.itp
-molecules defined in the itp file:
* Protein_chain_A
Rewrite topologie: no_cyclic_5vav_pdb2gmx.top
Protein_chain_A
Succeed to read file ...cyclic/top/no_cyclic_5vav_pdb2gmx.pdb ,  211 atoms found
Succeed to save file ...cyclic/top/5vav_pdb2gmx.pdb
>>> cyclic_pep.em(out_folder=TEST_OUT+'/cyclic/em/', nsteps=100, create_box_flag=True) #doctest: +ELLIPSIS
-Create pbc box
gmx editconf -f .../cyclic/top/5vav_pdb2gmx.pdb -o .../cyclic/top/5vav_pdb2gmx_box.pdb -bt dodecahedron -d 1.0
-Create the tpr file  5vav.tpr
gmx grompp -f 5vav.mdp -c ../top/5vav_pdb2gmx_box.pdb -r ../top/5vav_pdb2gmx_box.pdb -p ../top/5vav_pdb2gmx.top -po out_5vav.mdp -o 5vav.tpr -maxwarn 1
-Launch the simulation 5vav.tpr
gmx mdrun -s 5vav.tpr -deffnm 5vav -nt 0 -ntmpi 0 -nsteps -2 -nocopyright

Note

No options are allowed (forcefield, water model, termini capping) except for vsites.

Warning

Has not been tested with special residues like GLY or PRO !!

display()

Display defined attribute of the GmxSys object.

display_history()

Show all history

em(out_folder, name=None, nsteps=1000, posres='', create_box_flag=False, monitor=None, **mdp_options)

Minimize a system.

Parameters:
  • out_folder (str) – path of the output file folder
  • name (str, default=None) – name of the simulation to run
  • nsteps (int, default=1000) – number of minimisation steps
  • posres (str, default="") – option for the define variable in the mdp file, need to be define to have postion restraints
  • create_box_flag – flag to create or not a box to the input coor file.
  • mdp_options (dict) – Additional mdp parameters to use

Object requirement(s):

  • self.coor_file
  • self.top_file
  • self.nt
  • self.ntmpi
  • self.gpu_id

Object field(s) changed:

  • self.mdp
  • self.tpr
  • self.coor_file
  • self.xtc
em_2_steps(out_folder, name=None, no_constr_nsteps=1000, constr_nsteps=1000, posres='', create_box_flag=False, monitor=None, **mdp_options)

Minimize a system in two steps:

  1. minimisation without bond constraints
  2. minimisation using bond constraint for bonds involving hydrogen
Parameters:
  • out_folder (str) – path of the output file folder
  • name (str, default=None) – name of the simulation to run
  • no_constr_nsteps (int, default=1000) – number of minimisation steps in the first step
  • constr_nsteps (int, default=1000) – number of minimisation steps in the second step
  • posres (str, default="") – option for the define variable in the mdp file, need to be define to have postion restraints
  • create_box_flag (bool, optional, default=False) – flag to create or not a box to the input coor file.
  • mdp_options (dict) – Additional mdp parameters to use

Object requirement(s):

  • self.coor_file
  • self.top_file
  • self.nt
  • self.ntmpi
  • self.gpu_id

Object field(s) changed:

  • self.mdp
  • self.tpr
  • self.coor_file
  • self.xtc
em_equi_three_step_iter_error(out_folder, name=None, no_constr_nsteps=1000, constr_nsteps=1000, pdb_restr=None, nsteps_HA=100000, nsteps_CA=200000, nsteps_CA_LOW=400000, dt=0.005, dt_HA=0.002, maxwarn=0, iter_num=3, monitor=None, **mdp_options)

Minimize a system in 2 steps:

  1. minimisation without bond constraints
  2. minimisation using bond constraint for bonds involving hydrogen

Equilibrate a system in 3 steps:

  1. equilibration of nsteps_HA with position restraints on Heavy Atoms with dt = dt_HA
  2. equilibration of nsteps_CA with position restraints on Carbon Alpha with dt = dt
  3. equilibration of nsteps_CA_LOW with position restraints on Carbon Alpha with Low restraints with dt = dt

In case this process will crash (eg. LINCS WARNING …), the process will be rerun for iter_num time.

Parameters:
  • out_folder (str) – path of the output file folder
  • name (str, default=None) – name of the simulation to run
  • no_constr_nsteps (int, default=1000) – number of minimisation steps in the first step
  • constr_nsteps (int, default=1000) – number of minimisation steps in the second step
  • pdb_restr (str, default=None) – reference coordinate file for position restraints
  • nsteps_HA (int, default=100000) – number of equilibration steps with HA constraints
  • nsteps_CA (int, default=200000) – number of equilibration steps with CA constraints
  • nsteps_CA_LOW (int, default=400000) – number of equilibration steps with CA_LOW constraints
  • dt_HA (float, default=0.002) – integration time step for HA equilibration
  • dt (float, default=0.005) – integration time step for CA and CA_LOW equilibration
  • mdp_options (dict) – Additional mdp parameters to use

Object requirement(s):

  • self.coor_file
  • self.top_file
  • self.nt
  • self.ntmpi
  • self.gpu_id

Object field(s) changed:

  • self.mdp
  • self.tpr
  • self.coor_file
  • self.xtc
equi_three_step(out_folder, name=None, pdb_restr=None, nsteps_HA=100000, nsteps_CA=200000, nsteps_CA_LOW=400000, dt=0.005, dt_HA=0.002, maxwarn=0, monitor=None, **mdp_options)

Equilibrate a system in 3 steps:

  1. equilibration of nsteps_HA with position restraints on Heavy Atoms with dt = dt_HA
  2. equilibration of nsteps_CA with position restraints on Carbon Alpha with dt = dt
  3. equilibration of nsteps_CA_LOW with position restraints on Carbon Alpha with Low restraints with dt = dt
Parameters:
  • out_folder (str) – path of the output file folder
  • name (str, default=None) – name of the simulation to run
  • pdb_restr (str, default=None) – reference coordinate file for position restraints
  • nsteps_HA (int, default=100000) – number of equilibration steps with HA constraints
  • nsteps_CA (int, default=200000) – number of equilibration steps with CA constraints
  • nsteps_CA_LOW (int, default=400000) – number of equilibration steps with CA_LOW constraints
  • dt_HA (float, default=0.002) – integration time step for HA equilibration
  • dt (float, default=0.005) – integration time step for CA and CA_LOW equilibration
  • mdp_options (dict) – Additional mdp parameters to use

Object requirement(s):

  • self.coor_file
  • self.top_file
  • self.nt
  • self.ntmpi
  • self.gpu_id

Object field(s) changed:

  • self.mdp
  • self.tpr
  • self.coor_file
  • self.xtc
extend_equi_prod(tpr_file=None, nsteps=200000, dt=0.005, monitor=None)

Extend a simulation run.

Parameters:
  • tpr_file – path of the tpr file
  • nsteps (int, default=200000) – number of steps
  • dt (float, default=0.005) – integration time step use previously

Object requirement(s):

  • self.tpr
  • self.nt
  • self.ntmpi
  • self.gpu_id

Object field(s) changed:

  • self.tpr
  • self.sim_name
  • self.coor_file
  • self.xtc

Warning

The simulation can only run nsteps more steps. Should be improved to check how many steps have been run, and give a total number of step to compute.

get_all_output()

In a case of a simulation restart, outputs edr, log, gro and xtc files are called for example as self.sim_name+".partXXXX.edr" where XXXX is the iteration number of restart (eg. first restart: XXXX=0002).

This function return a dictionnary of all edr, log, coor_file and xtc list.

Object requirement(s):

  • self.sim_name
Returns:return dict containing edr, log, xtc and coor_file file list
Return type:dict
get_ener(selection_list, output_xvg='tmp.xvg', check_file_out=True, keep_ener_file=False)

Get energy of a system using gmx energy.

get_last_output()

In a case of a simulation restart, outputs edr, log, gro and xtc files are called for example as self.sim_name+".partXXXX.edr" where XXXX is the iteration number of restart (eg. first restart: XXXX=0002).

This function actualise the edr, log, coor_file and xtc variable.

Object requirement(s):

  • self.sim_name

Object field(s) changed:

  • self.edr
  • self.log
  • self.coor_file
  • self.xtc
get_simulation_time()

In a case of a simulation restart simulation, one would like to know how much simulation time has already been computed to reach a certain amount of time in the restart simulation. The command will check the cpt file using gmx check -f file.cpt.

Object requirement(s):

  • self.tpr

Object field(s) changed:

  • None
Returns:return simulation time (ns)
Return type:float
insert_mol_sys(mol_gromacs, mol_num, new_name, out_folder, check_file_out=True)

Insert a new molecule in a system:

Insert structure and topologie of mol_num copy of mol_gromacs molecule, in the system with 6 successive steps:

  1. Copy the molecule mol_num time
  2. Change the chain ID of mol_gromacs to “Y”, this step is necessary for vmd to recognize the inserted mol.
  3. Concat the two structure
  4. Insert the molecule in the solvant with a vmd script
  5. Update the topologie with the molecule and new water number
  6. If the charge is not null add ions to neutralize the system
Parameters:
  • mol_gromacs (GmxSys object) – molecule object to be inserted
  • mol_num (int) – molecule number to be inserted
  • new_name (str) – generic name of the system
  • out_folder (str) – path of the output file folder
  • check_file_out (bool, optional, default=True) – flag to check or not if file has already been created. If the file is present then the command break.

Object requirement(s):

  • self.coor_file
  • self.top_file

Object field(s) changed:

  • self.coor_file
  • self.top_file

Note

VMD don’t need anymore to be installed to run the peptide creation

prepare_top(out_folder, name=None, vsite='hydrogens', ignore_ZN=True)

Prepare the topologie of a protein:

  1. compute hisdine protonation with pdb2pqr
  2. Change Histidine resname according to the protonation
  3. Correct cystein resname
  4. Correct chain ID’s
  5. Zinc Finger: Add Zinc in the pdb and change residue type of CYS and HIS coordinating the Zinc
  6. Finally compute the topologie with pdb2gmx add_top()
Parameters:
  • out_folder (str) – path of the output file folder
  • name (str, optional, default=None) – generic name of the system
  • vsite (str, optional, default="hydrogens") – option for topologie’s bonds constraints (“none”, “hydrogens”, “all”)
  • ignore_ZN (bool, optional, default=False) – option for not adding parameters to ZINC finger

Object requirement(s):

  • self.coor_file

Object field(s) changed:

  • self.coor_file
  • self.top_file
Example: 
>>> TEST_OUT = getfixture('tmpdir')
>>> # Create the topologie of a protein and do a minimisation:
>>> prot = GmxSys(name='1y0m', coor_file=TEST_PATH+'/1y0m.pdb')
>>> prot.prepare_top(out_folder=TEST_OUT+'/prepare_top/top_SH3/') #doctest: +ELLIPSIS
Succeed to read file .../input/1y0m.pdb ,  648 atoms found
Succeed to save file tmp_pdb2pqr.pdb
pdb2pqr.py --ff CHARMM --ffout CHARMM --chain tmp_pdb2pqr.pdb 00_1y0m.pqr
Succeed to read file 00_1y0m.pqr ,  996 atoms found
Chain: A  Residue: 0 to 60
Succeed to save file 01_1y0m_good_his.pdb
-Create topologie
gmx pdb2gmx -f 01_1y0m_good_his.pdb -o 1y0m_pdb2gmx.pdb -p 1y0m_pdb2gmx.top -i 1y0m_posre.itp -water tip3p -ff charmm36-jul2017 -ignh -vsite hydrogens
Molecule topologie present in 1y0m_pdb2gmx.top , extract the topologie in a separate file: 1y0m_pdb2gmx.itp
Protein_chain_A
-ITP file: 1y0m_pdb2gmx.itp
-molecules defined in the itp file:
* Protein_chain_A
Rewrite topologie: 1y0m_pdb2gmx.top

Note

No options are allowed (forcefield, water model, termini capping) except for vsites.

Note

Starting file need to be a pdb, this should be changed.

production(out_folder, name=None, nsteps=400000, dt=0.005, maxwarn=0, monitor=None, **mdp_options)

Run a production run.

Parameters:
  • out_folder (str) – path of the output file folder
  • name (str, default=None) – name of the simulation to run
  • nsteps (int, default=400000) – number of minimisation steps
  • dt (float, default=0.005) – number of minimisation steps
  • mdp_options (dict) – Additional mdp parameters to use

Object requirement(s):

  • self.coor_file
  • self.top_file
  • self.nt
  • self.ntmpi
  • self.gpu_id

Object field(s) changed:

  • self.mdp
  • self.tpr
  • self.coor_file
  • self.xtc
run_md_sim(out_folder, name, mdp_template, mdp_options, pdb_restr=None, maxwarn=0, monitor=None)

Run a simulation using 3 steps:

  1. Create a mdp file
  2. Create a tpr file using gmx grompp
  3. Launch the simulation using gmx mdrun
Parameters:
  • out_folder (str) – path of the output file folder
  • name (str) – name of the simulation to run
  • mdp_template (str) – mdp file template
  • mdp_options (dict) – New parameters to use
  • pdb_restr (str, default=None) – reference coordinate file for position restraints
  • maxwarn (int, default=0) – Maximum number of warnings when using gmx grompp
  • monitor – option to monitor a simulation, if not none monitor should contains two values: function the function to be ran while simulation is running and input parameters for the function

Object requirement(s):

  • self.coor_file
  • self.top_file
  • self.nt
  • self.ntmpi
  • self.gpu_id

Object field(s) changed:

  • self.mdp
  • self.tpr
  • self.coor_file
  • self.xtc
run_simulation(check_file_out=True, cpi=None, nsteps=-2, rerun=False, monitor=None)

Launch the simulation using gmx mdrun

Parameters:
  • check_file_out (bool, optional, default=True) – flag to check or not if file has already been created. If the file is present then the command break.
  • cpi (str, default=None) – checkpoint file, if defined, it will restart a simulation and run nsteps
  • nsteps (int, default=-2) – Number of steps to run, (-2 : will use mdp parameter)
  • rerun (dict, default=None) – option to rerun a simulation (eg. recompute energy)
  • monitor – option to monitor a simulation, if not none monitor should contains two values: function the function to be ran while simulation is running and input parameters for the function

Object requirement(s):

  • self.tpr
  • self.sim_name
  • self.nt
  • self.ntmpi
  • self.gpu_id

Additional requirement(s) for rerun:

  • self.xtc

Object field(s) changed:

  • self.coor_file
  • self.xtc
  • self.edr
  • self.log

Note

The function must be launched in the path where the tpr is present.

Note

If cpi file is defined the simulation will restart with the -noappend option, if cpi is not defined, but the .cpt file exist, it will restart with “append”.

save_state()

Save last state

solvate_add_ions(out_folder, name=None, ion_C=0.15, create_box_flag=True, box_dist=1.1)

Solvate a system with three succesive steps:

  1. Create box using create_box()
  2. Add water using solvate_box()
  3. Add ions using add_ions()
Parameters:
  • out_folder (str) – path of the output file folder
  • name (str, optional, default=None) – generic name of the system
  • ion_C (float, optional, default=0.15) – ionic concentraton (Molar)

Object requirement(s):

  • self.coor_file
  • self.top_file

Object field(s) changed:

  • self.coor_file solvate_box
  • self.top_file
Example: 
>>> TEST_OUT = getfixture('tmpdir')
>>> prot = GmxSys(name='1y0m', coor_file=TEST_PATH+'/1y0m.pdb')
>>> prot.add_top(out_folder=TEST_OUT+'/solvate_add_ions/top_SH3/') #doctest: +ELLIPSIS
-Create topologie
gmx pdb2gmx -f .../input/1y0m.pdb -o 1y0m_pdb2gmx.pdb -p 1y0m_pdb2gmx.top -i 1y0m_posre.itp -water tip3p -ff charmm36-jul2017
Molecule topologie present in 1y0m_pdb2gmx.top , extract the topologie in a separate file: 1y0m_pdb2gmx.itp
Protein_chain_A
-ITP file: 1y0m_pdb2gmx.itp
-molecules defined in the itp file:
* Protein_chain_A
Rewrite topologie: 1y0m_pdb2gmx.top
>>> prot.solvate_add_ions(out_folder=TEST_OUT+'/solvate_add_ions/top_SH3_water_ions/') #doctest: +ELLIPSIS
-Create pbc box
gmx editconf -f .../solvate_add_ions/top_SH3/1y0m_pdb2gmx.pdb -o .../solvate_add_ions/top_SH3/1y0m_pdb2gmx_box.pdb -bt dodecahedron -d 1.1
-Solvate the pbc box
Copy topologie file and dependancies
Copy topologie file and dependancies
-Create the tpr file  genion_1y0m_water_ion.tpr
gmx grompp -f .../gromacs/template/mini.mdp -c 1y0m_water.pdb -r 1y0m_water.pdb -p 1y0m_water_ion.top -po out_mini.mdp -o genion_1y0m_water_ion.tpr -maxwarn 1
-Add ions to the system with an ionic concentration of 0.15 M , sytem charge = 0.0 water num= 62...
Add ions : NA : 16   CL : 16
gmx genion -s genion_1y0m_water_ion.tpr -p 1y0m_water_ion.top -o 1y0m_water_ion.gro -np 16 -pname NA -nn 16 -nname CL
>>> prot.em(out_folder=TEST_OUT+'/solvate_add_ions/em_SH3_water_ions/', nsteps=100, constraints = "none")
-Create the tpr file  1y0m.tpr
gmx grompp -f 1y0m.mdp -c ../top_SH3_water_ions/1y0m_water_ion.gro -r ../top_SH3_water_ions/1y0m_water_ion.gro -p ../top_SH3_water_ions/1y0m_water_ion.top -po out_1y0m.mdp -o 1y0m.tpr -maxwarn 1
-Launch the simulation 1y0m.tpr
gmx mdrun -s 1y0m.tpr -deffnm 1y0m -nt 0 -ntmpi 0 -nsteps -2 -nocopyright

Note

If name is not defined, it will use the object name.

solvate_box(out_folder, name=None, radius=0.21, cs='share/gromacs/top/spc216.gro', check_file_out=True)

Solvate the pbc box with water or another mol defined with cs using the gmx solvate command.

Parameters:
  • out_folder (str) – path of the output file folder
  • name (float, optional, default=0.21) – generic name of the system
  • radius – default van der Waals distance (nm)
  • cs (str, optional, default=``WATER_GRO``) – solvant coordinate file
  • check_file_out (bool, optional, default=True) – flag to check or not if file has already been created. If the file is present then the command break.

Object requirement(s):

  • self.coor_file
  • self.top_file

Object field(s) changed:

  • self.coor_file
  • self.top_file
Example: 
>>> TEST_OUT = getfixture('tmpdir')
>>> prot = GmxSys(name='1y0m', coor_file=TEST_PATH+'/1y0m.pdb')
>>> prot.add_top(out_folder=TEST_OUT+'/solv_box/top_SH3/') #doctest: +ELLIPSIS
-Create topologie
gmx pdb2gmx -f .../input/1y0m.pdb -o 1y0m_pdb2gmx.pdb -p 1y0m_pdb2gmx.top -i 1y0m_posre.itp -water tip3p -ff charmm36-jul2017
Molecule topologie present in 1y0m_pdb2gmx.top , extract the topologie in a separate file: 1y0m_pdb2gmx.itp
Protein_chain_A
-ITP file: 1y0m_pdb2gmx.itp
-molecules defined in the itp file:
* Protein_chain_A
Rewrite topologie: 1y0m_pdb2gmx.top
>>> prot.create_box() #doctest: +ELLIPSIS
-Create pbc box
gmx editconf -f .../solv_box/top_SH3/1y0m_pdb2gmx.pdb -o .../solv_box/top_SH3/1y0m_pdb2gmx_box.pdb -bt dodecahedron -d 1.0
>>> prot.solvate_box(out_folder=TEST_OUT+'/solv_box/top_SH3_water/')
-Solvate the pbc box
Copy topologie file and dependancies

Note

If name is not defined, the command will create a new .pdb and .top file name after the object name and adding “_water”.