QM/MM protocols

Table of Contents

Maestro

Download the tar.gz from the Schrödinger homepage, unpack where you want it and go to /path/to/schrodinger (could be /home/erikh/Programs/Maestro_2017-1_Linux-x86_64_Academic).

  1. run install ./INSTALL - note what dir you install it to...
  2. Add export SCHRODINGER=path/to/installdir/ to .bash_profile (install dir could be opt/schrodinger2017-1)
  3. NOTE: On fedora, I experience problems with some libraries. This was solved by export LD_PRELOAD=/usr/lib64/libstdc++.so.6:/lib64/libgcc_s.so.1

ComQum

A complete is also given in Ryde Comqum

Setup of ComQum Calculations

A ComQum QM/MM calculation requires a protein setup first. If you have not done this, go to the protein setup section.

Creating a spherical system from the protein setup

During the Simulated annealing step of the protein setup, we employ an octahedral water cap. This needs to be spherical in QM/MM calculations (since we do not have periodic boundary conditions in QM/MM). In a first step we therefore create a spherical system using the script below. Note: If the old procedure has been followed, the simulated annealing is done with a spherical water cap and therefore this step can be skipped.

  • run the program cpptraj with cpptraj prmtop ptraj.in. Use the ptraj.in file below ptraj.in

    trajin mdrest3
trajout mdrest3-wrap restart
image origin center familiar com :98 # central residue
go

    Note: The central residue can be found with changepdb as described in part 9. of the protein setup

  • Run changeparm

    changeparm <<EOF
prmtop
p
pdb3.pdb
m
mdrest3-wrap
q
EOF

  • The system can now be cut, e.g., with changepdb

    changepdb <<EOF
pdb3.pdb
m
p
cup
<press enter>
0 0 0
45     # Here different radii can be tested (ensure the system is still spherical)
n
pdb45-cut.pdb
n
EOF

  • rename the file pdb45-cut.pdb to comqum.pdb

Define QM system

  • Prepare syst1 file with QM system

    The syst1 file contains the QM system. It is structures as

    Title
xx--yy pdb atom numbers for atoms in the QM region.

x junction atoms
y

    Example syst1 file

    This is a title: LPMO system
1-27 # First residue in QM region (numbers refer to pdb file) 
1259 # First atom of next residue in QM region (junction atom) 
1261-1271 # here comes the residue 
2576 # more residues (junction) 
2578-2592 # more... 
3358-3364 # This is a metal and a few waters 

27 # Junctions - Caps are always CH3 groups junctions should be between C-C bonds 
1259 
2576

  • run pdbtocomqum

    pdbtocomqum

    1. type pdbtocomqum in a terminal
    2. select logfile
    3. n (no short contacts)
    4. enter title
    5. specify QM system (best done by file "syst1", see below for the syntax of this file)
    6. enter cut off radius for system 2 (normally 6 is fine). System 2 is the MM system that is optimized
    7. N (do not remove amino acids from system 2)
    8. N (do not include more amino acids in system 2)
    9. enter cut-off radius for system 3 (1000 is fine). System 3 is the MM system that is fixed
    10. use new format
    11. enter junction atoms. This can also be through the syst1 file. Check that the junctions are known (i.e. not -1)
    12. do not remove charge (n)
    13. do not redistribute charges (n)

  • Make turbomole and amber files.

    Run comqumtoturbo and comqumtoamber

  • Make partop3 prmcrd3 with tleap. Copy the leap.in file (and all .in and .dat files + the leapprc file) from the equilibrium run and modify it:

    x=loadpdb pdb.in3
saveamberparm x prmtop3 prmcrd3

    Delete the line with solventcap (e.g. like solvatecap x TIP3PBOX {0.0 0.0 0.0} 40). Then run

    tleap -s -f leap.in

    Typical problems with tleap:

    1. Leap connects everything that is not separated by "TER". Therefore, make sure to separate groups that should not be connected with TER in the input file

leap.in sample file

  source leaprc.ff14SB
  source leaprc.GLYCAM_06j-1
  AddPath /lunarc/nobackup/users/erikh/lpmo/pmo-complex/Leap
  loadAmberPrep wam.in
  loadAmberPrep hic.in
  loadAmberPrep cu2.in
  loadAmberPrep o2m.in
  loadAmberParams wam.dat
  loadAmberParams hic.dat
  loadAmberParams cu2.dat
  loadAmberParams o2.dat
  x=loadpdb pdb.in3
  bond x.41.SG x.167.SG
  saveamberparm x prmtop3 prmcrd3
  quit

Typical problems with tleap:

  1. Leap connects everything that is not separated by "TER". Therefore, make sure to separate groups that should not be connected with TER in the input file

  • Modify charges of non-standard residues

    If there was non-standard residues in the equilibrium calculations, they should be modified, also in the comQum calculation. Use the script changeparm

    changeparm <<EOF
prmtop3
m
comqum.pdb
w

q
EOF

  • Add junctions to comqum.dat

  • Run cqprep

    cqprep

    1. run cqprep and follow the instructions (a cap is inserted - use the data from the equilibrium run).
    2. prmtop1 is generated from prmtop3. What is also done: In prmtop1 is a (rudementary) force field for the QM system, so that the MM energy of the QM system can be calculated (given in sander.out1).
    3. Charges from prmtop1 are removed and pdbout1 is written out (what is in this file?)
    4. Charges from system 1 are removed from prmtop2 and pdbout3 is written (what is in this file?)
    5. Test run sander.in1 (write down energies)
    6. Test run sander.in2 (write down energies)
    7. Define is started (select DFT method, bases etc. manually)

  • Clean up gzip *; cqgunzip; mkdir Gz; mv *.gz Gz; cp * Gz; gzip Gz/*&

  • Protein free runs..

    1. Replace $protein fixed with $protein free in comqum.dat
    2. Insert $esp_fit kollman in the control file
    3. Start comqum normally
    4. In case of problems run:
      dscf >logd (or ridft >logd)
ridft - proper mulliken (comment out the pointcharge line in the control file)
fixcharge58_turboin >> fixc
fixcharge_amberin >> fixc
fixcharge>>fixc
    5. The are sometimes problems with the charges (one problem can occur if the QM system is slightly different from the system used to obtained the RESP charges for non-standard residues). This can be solved by
      • For all residues where there are differences between the current QM system and the system used to obtained charges, replace the charges of these residues with standard AMBER charges (see e.g. another residue of the same amino acid in the pdb file) in the system pdb file
      • Insert the new charges in the prmtop3 file with changeparm (see above)
      • run 
        fixcharge_amberin
fixcharge_turboin
fixcharge
fixcharge_amberout

Note: Comqum expects that all residues are neutral (as done in the AMBER force field). Some force fields woks with non-neutral residues (e.g. the GLYCAM force field). If some residues are not neutral, it can be necessary to add a correction via comqum.dat. This should only be done with care! normally, a error in the fixcharge is rather due to a problem in the setup. An example of where a charge-correction is necessary is when a force field, e.g., for a substrate is defined as several residues with different charge. In this case, a charge correction can be necessary if the QM region only includes some of the residues in the substrate (since ComQum expects neutral residues).

$residue_charge_corr
n # numer of residue to be corrected
m c # m= residue id, c = correction

Example

$residue_charge_corr
1
240 0.1940

Run then

fixcharge_amberin
fixcharge_turboin
fixcharge
fixcharge_amberout
  • MM analysis

In some (rare) cases it is possible to get a large differences by energies in the MM parts: This is often (but not always) caused by vdW terms. Changeparm can be used to analyse the contributions from individual terms. This is done as follows:

Collect the comqum calculations in two directories called dir1 and dir2

changeparm
path_to_dir1/prmtop2
qmd # This is the analysis command
path_to_dir1/prmtop1
m
path_to_dir1/prmcrd3
m
path_to_dir1/prmcrd1
m
path_to_dir2/prmcrd3
m
path_to_dir2/prmcrd1

ComQum-X

Here goes comqum-X description...

Big-QM

The set-up big-QM calculations proceeds via. the changepdb program. The big QM calculations construct a system that includes:

  1. All residues within 4.5 Å to 6 Å
  2. All buried charges
  3. All junctions are moved 3 residues away from the "usual" small QM system

big-qm with changepdb

  1. Start from an pdb ("pdb3.pdb") file that contains the QM/MM optimised system. Copy the files from the comqum optimization along with the syst1 file that was used to set-up the system to a new folder "big-qm".
  2. Rename syst1 to s1.
  3. Run changepdb on command "bc" on pdb3.pdb

Polarizable Embedding

Calculations of potentials for polarizable embedding calculations in DALTON or DIRAC are done with the polarizable embedding assistant script (PEAS). The script and a general installation guide is given below. The script requires

  • Python (above version 2.7)
  • Installations of either MOLCAS or DALTON (required to calculate localized multipoles and polarizabilities through the LoProp protocol). Note that for DALTON, additional python scripts are required

A general installation guide can be found at on the gitlab page: peas gitlab

Installation on aurora

  • Get peas from gitlab 
    git clone git@gitlab.com:pe-software/peas.git peas-aug2018
  • Follow the install instructionspython 
    setup.py build
  • Go to build dir cd build and install at an appropriate place e.g.
    python setup.py install --prefix=/home/erikh/programs/peas
  • Set environment variables (see examples below)

g-fortran/ompi 

module purge
module load GCC/5.4.0-2.26 OpenMPI/1.10.3 # Alternative with gfortran/ompi compiler
module load Python/2.7.12 
export PEASHOME=/home/erikh/programs/peas
export LD_LIBRARY_PATH=$LD_LIBRARY_PATH:$PEASHOME/lib:$LD_LIBRARY_PATH
export INCLUDE=$INCLUDE:$PEASHOME/lib
export LD_RUN_PATH=$LD_RUN_PATH:$PEASHOME/lib
export MANPATH=$MANPATH:$PEASHOME/share/man
export PKG_CONFIG_PATH=$PEASHOME/lib/pkgconfig:$PKG_CONFIG_PATH
export PYTHONPATH=$PEASHOME/lib/python2.7/site-packages    
PATH=$PATH:$PEASHOME/bin
export PATH

Compatible with Molcas intel 

module purge
module load foss/2016b
module load CMake/3.5.2
module load HDF5/1.8.18
module load GSL/2.1
module load Python/2.7.12
module load iomkl/2017.01    
export PEASHOME=/home/erikh/programs/peas
export LD_LIBRARY_PATH=$LD_LIBRARY_PATH:$PEASHOME/lib:$LD_LIBRARY_PATH
export INCLUDE=$INCLUDE:$PEASHOME/lib
export LD_RUN_PATH=$LD_RUN_PATH:$PEASHOME/lib
export MANPATH=$MANPATH:$PEASHOME/share/man
export PKG_CONFIG_PATH=$PEASHOME/lib/pkgconfig:$PKG_CONFIG_PATH
export PYTHONPATH=$PEASHOME/lib/python2.7/site-packages    
PATH=$PATH:$PEASHOME/bin
export PATH
  • Put the variables into a script peas.sh and place it /somewhere/peas.sh (usage: source /somewhere/peas.sh)

Running peas (examples)

  • Set environment variables
    source /home/erikh/programs/.peas.sh

results matching ""

    No results matching ""