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:
- Load the structure into SYBYL
- Save the structure as a PDB file
- Zap the current molecule from SYBYL
- Load in the molecule in the PDB file in its place
- Check the atom and bond types of any heteroatoms to ensure that they are correct; if not then edit them as required
- Save the structure as a mol2 file You are now ready to use this new mol2 file to set up your CHARMM input.
- Load Insight-II via the command:
- insightii
- 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}
- 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}
- 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}
- 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}
- 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
- 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
- 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.
- Check in a UNIX shell to ensure that each of the following files are present in your directory:
- jobname.rtf
- jobname.psf
- jobname.crd
- If all looks to be okay, save your molecule in standard Insight-II format:
- File >> Save_Folder
- {Execute}
- Quit InsightII:
- Session >> Quit
- {Execute}
- 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 GROUP ATOM 1CB CT -0.10 ATOM 1SG ST -0.00 GROUP ATOM 2SG ST -0.00 ATOM 2CB CT -0.10 DELETE ATOM 1HG DELETE ATOM 2HG 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 - ssh to the cluster via the following syntax:
- ssh kumgmc01.msg.ku.edu
- Obtain a template charmm directory for your calculation:
- cp -r /home/ghl/charmm/template .
- Change the directory name from "template" to your own job name:
- mv template jobname
- 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
- 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
- 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
- Run the CHARMM job via:
- chms jobname
- 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
- 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. - 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