Preparing and Running CHARMM Torsional Analyses for Small Peptides

Pre-requisites: this tutorial presumes that you're starting with a suitably protonated structure file in mol2 format. If you instead have your structure in PDB format, you can start at step 4. If you generated your structure in SYBYL, there will frequently be information in the mol2 file that will confuse Insight-II. This can be eliminated by converting the structure to PDB and back again as follows:

  1. Load the structure into SYBYL
  2. Save the structure as a PDB file
  3. Zap the current molecule from SYBYL
  4. Load in the molecule in the PDB file in its place
  5. Check the atom and bond types of any heteroatoms to ensure that they are correct; if not then edit them as required
  6. Save the structure as a mol2 file You are now ready to use this new mol2 file to set up your CHARMM input.
  7. Load Insight-II via the command:
    • insightii
  8. Load your structure file into Insight-II
    • Molecule >> Get
    • In the "Get Molecule" window, click on the "SYBYL_Mol2" button
    • If you are not in the directory containing your file, click on the appropriate directory (and subdirectories if necessary) in the "Parameters" window until your file is listed as an option
    • Click on the desired file.
    • In the "Get Molecule" window, click on {Execute}
  9. Select the force field:
    • FF >> Select
    • In the "Parameters" window, click the "Charmm" button
    • In the "Parameters" window, click "CHARMM.frc"
    • In the "Select Forcefield" button, click the "Clear_Charges" button
    • In the "Select Forcefield" button, click {Execute}
  10. Set force field parameters:
    • FF >> Potentials
    • In the "Potentials Forcefield" button, click on each of the three "Fix" buttons
    • In the "Potentials Forcefield" button, click {Execute}
  11. Check your atom types via:
    • Molecule >> Label
    • In the "Parameters" window, click "Atom"
    • In the "Label Molecule" window, the "Label option" button should be set to "Single_property", the "Label level" to "Atom" and the "Label action" to "On"
    • In the "Label Molecule" window, click on the "Label Property" box
    • From the "Parameters" window, select "Potential"
    • In the "Label Molecule" window, click {Execute}
  12. Verify the following:
    • Amide N's should be "NP", except for proline ("NX")
    • Amide H's should be "H"
    • Carboxyl C's should be "C"
    • Carboxyl O's should be "O"
    • Alkyl (and alpha) C's should be "CT"
    • Alkyl H's should be "HA"
    • For other atom types, check with the instructor
  13. If an atom type is invalid:
    • MSI >> Builder
    • Atom >> Potential
    • Click on the atom type you wish to change an atom to
    • Click on the specific atom with the mis-assigned type
    • Do so for all mis-assigned atoms
    • {Cancel} to quit
  14. Save your molecule in a format that CHARMM can read:
    • Molecule >> Put
    • In the "Put Molecule" window, click on the "CHARMm" button
    • Click in the "Mol File Name" text area, and enter an appropriate name
    • Make sure that each of the "Transformed", "Create_PSF File" and "Create_RTF File" switches are turned on.
  15. Check in a UNIX shell to ensure that each of the following files are present in your directory:
    • jobname.rtf
    • jobname.psf
    • jobname.crd
    where "jobname" is the name you gave to the system in step 7
  16. If all looks to be okay, save your molecule in standard Insight-II format:
    • File >> Save_Folder
    • {Execute}
  17. Quit InsightII:
    • Session >> Quit
    • {Execute}
  18. Edit your "jobname.rtf" file and delete all of the following lines found near the end of the file:
    PRES DISU        -0.20
    ! This patch is used to provide a disulphide bridge between two CYS residues
     ATOM 1CB  CT    -0.10
     ATOM 1SG  ST    -0.00
     ATOM 2SG  ST    -0.00
     ATOM 2CB  CT    -0.10
    BOND 1SG 2SG
    ANGLE  1CB 1SG 2SG      1SG 2SG 2CB
    DIHE   1CA 1CB 1SG 2SG      1HB1 1CB 1SG 2SG      1HB2 1CB 1SG 2SG
    DIHE   1CB 1SG 2SG 2CB
    DIHE   1SG 2SG 2CB 2CA      1SG 2SG 2CB 2HB1     1SG 2SG 2CB 2HB2
    IC     1CA 1CB 1SG 2SG    0.0  0.0  180.0  0.0  0.0
    IC     1CB 1SG 2SG 2CB    0.0  0.0  180.0  0.0  0.0
    IC     1SG 2SG 2CB 2CA    0.0  0.0  180.0  0.0  0.0
    but remember to leave the final "END" line in place
  19. ssh to the cluster via the following syntax:
    • ssh
  20. Obtain a template charmm directory for your calculation:
    • cp -r /home/ghl/charmm/template .
  21. Change the directory name from "template" to your own job name:
    • mv template jobname
  22. Enter the new directory, make a note of the full directory path, and register this directory into the "datadir.def" file:
    • cd jobname
    • pwd
    • (record the output from the previous step)
    • vi datadir.def
    • change "/home/ghl/charmm/template" to "/home/steph/jobname" on both lines 4 and 5
    • write the new datadir.def file and quit
  23. Copy the files you generated from Insight-II into the "data" directory
    • cp /usr/people3/steph/jobname/jobname.crd ./data
    • cp /usr/people3/steph/jobname/jobname.psf ./data
    • cp /usr/people3/steph/jobname/jobname.rtf ./data
    Note: obviously you should exchange "/usr/people3/steph/jobname" for whatever is the full path name of the directory where you created these files.
  24. There are two template input files for you to use: anneal.inp and analysis.inp, where the first is a standard setup for simulated annealing and the second one is suitable for carrying out torsional analysis. Suppose you are going to do simulated annealing:
    • cp anneal.inp jobname.inp
    • vi jobname.inp
    • Change "/home/ghl/charmm/template/datadir.def" to reflect the proper directory
    • Change all other instances of "test" to "jobname" via the following vi line command: :1,$s/template/jobname/g
    • Save the file and quit
  25. Run the CHARMM job via:
    • chms jobname
    And answer the question regarding which processor to use in an appropriate fashion (i.e., select a low-load processor).
  26. Running the analysis job is the same as in steps 23,24, except that in order to extract specific torsions in a format that you can use for analysis, you will need to run a post-processing command. First, look at your output file to find the specific number of the to be analyzed via:
    • vi jobname.out
  27. The output file should be large. In order to find the torsional data, execute the following vi command: / DIHEDRAL VALUES FOR STEP 
    Note that the "/" is important -- this is what tells vi to search for a character string. You should then see a list of torsions for the first data step. The second column in this list is the torsion number. Record this number for your torsion of interest, and exit from the vi session.
  28. Run the following post-processing command: exdih
    It will ask you first for the name of the output file to process (give the full file name) and then for the torsion number determined in step 27. It will output the relevant torsional data in a file that looks like "jobname.out.xx.list", where "jobname.out" is your output file and "xx" is the torsion number.

Gerry Lushington   8/23/2004

David Johnson
Access to more than 40 computational chemistry software programs and databases
High-performance computational tools accelerate drug discovery and minimize costs
Analyze complex, multidimensional data sets to rapidly generate biological insights
KU Today