Understanding & Creating BTCL Scripts, Instructions

Section 1 - Subsets and Atom Specifications

Section 1 of this tutorial makes use of a BTCL input file called seg1.inp, which illustrates the use of subsets and atom specs.

IMPORTANT: To run this tutorial, you need to copy all the files that are in $BIOSYM/tutorial/discover/BTCL into a directory in which you have write-permission.

The script in seg1.inp obtains arrays of atoms (which are actually references into the System database) directly and by using subsets. For most applications, the direct approach is fine. Subsets capture information in actual database constructs and are more useful in dealing with relationships, e.g., in recording pairs of atoms corresponding to a set of distances. For a subset of this type, the define statement would appear less redundant: subset define distance instead of subset define subset, where the second subset keyword just indicates that the subset being defined is generic.

Run the script and note the effects of the various commands. Things to try:

What happens if you change echo [object coord] (the last line of the file) to echo $coord?
Can you devise some more involved atom specs and run the resulting input file successfully? For ideas, you can look at the documentation of the select command.

Section 2

Section 2 of this tutorial makes use of seg2.inp, which illustrates some BTCL language constructs while accessing an archive file.

The script begins with the definition of a procedure called compareArchive. An alternative to including the full procedure definition would be to define it in a separate file--perhaps named compareArchive.tcl. The script could then incorporate the procedure by including the command source compareArchive.tcl.

The compareArchive procedure takes two direct inputs: a subset of atoms (backbone atoms to be used in an rms comparison) and an archive filename. The procedure also declares the input variable cmprIndx to be global, which means that it can be accessed at the top level and within the procedure.

The compareArchive procedure begins by saving all current system coordinates for subsequent restoration. It also stores coordinates for the input subset of atoms separately and computes the energy of the current system.

The procedure next uses the if and info exists constructs to determine whether or not to proceed. It also appends an else construct to the if. This is not required, but it is often useful, as is the similar elseif construct. If some conditions are satisfied, the readFile command is used to bring in a new set of coordinates from the archive file. Notice the difference between the readFile and the begin and reset commands, which delete everything and start from scratch. The archive coordinates are used to perform another energy calculation, and both archived and current coordinates for the subset of atoms are used for an rms comparison. The energy difference and the rms comparison value are printed out, and the original coordinates are restored.

The main script defines the subset of heavy atoms and initially calls the compareArchive procedure without setting the cmprIndx variable. It then tries again with $cmprIndx = 0, using a catch to avoid aborting the run if there is a problem with the use of index 0. The catch construct intercepts system errors and suspends automatic error processing. It is useful for testing risky commands and for replacing low-level error messages with higher-level ones. Finally, the script uses a for loop to try the values 1 and 2 for cmprIndx.

Run the script and verify that it works as it ought to. Things to try:

What happens if you try to read a non-existent frame?

What happens if you change:

	catch {compareArchive sbset energy_test}

to:

	catch {compareArchive sbset energfy_test}

Can you use a while loop instead of a for loop? Documentation on conditional loops in BTCL is available.

Section 3

Section 3 of this tutorial makes use of seg3.inp, which illustrates some basic uses of the molGeom and restraint commands, and also some TCL list commands (lappend, llength).

The script obtains a distance (using molGeom) and the energy, then alters and remeasures the distance (using molGeom), and recalculates the energy. It then minimizes the conformation and reobtains the distance and energy. Next, it applies a distance restraint, minimizes, and reobtains the distance and energy. Finally, it outputs the sequence of energies and distances obtained.

Note the difference in the restraint command between the Discover 94.0 and Discover 95.0 programs. The restraint create command is simpler in the Discover 95.0 program because more defaults are provided. The defaults can be altered using the restraint function and restraint target commands.

The Discover 95.0 program also supports restraint files. To read a restraint file, use the command:

	readRestraintFile filename

If the optional filename is not provided, the file-reader will look for $PROJECT.rstrnt.

Run the script and note the effects of the various commands. Things to try:

Use the molGeom command to modify some torsions instead of distances.
Try working with vectors instead of scalars.
Write a script to read in the file cram7helix.rstrnt and carry out some processing. The files cram7helix.car and cram7helix.mdf are available in $BIOSYM/tutorial/discover/BTCL--you should copy them to a directory in which you have write-permission.

Section 4

Section 4 of this tutorial is the creation of your own BTCL script. You should try writing a script that will:

Generate new conformations from the butane input conformation. (You should have copied butane.car and butane.mdf from $BIOSYM/tutorial/discover/BTCL into a directory in which you have write-permission and in which your script is located.)
Calculate the energy and heavy-atom rms deviation between a new conformation and each previously stored conformation.
Store the new conformation if its energy or heavy-atom rms deviation differs by more than a defined tolerance from each previously stored conformation.

Sections one, two, and three of this lesson provide examples of the BTCL constructs you will need to construct your script.

Start your script with the definition of a procedure to be called compareArchive. This procedure should take four inputs:

subset - the subset consisting of the backbone atoms
filename - archive filename
enerTol - energy difference tolerance
rmsTol - rms deviation tolerance

The procedure should record the energies of the stored conformations in a list stored in a global variable called energies.

Start the procedure by saving all current system coordinates, so that they can be restored later. Also store coordinates for the input subset of atoms separately and compute the energy of the current system.

Use the if and info exists constructs to determine whether there is anything in the list contained in the energies global variable.

If it is empty, store the conformation.
If it is not empty, loop over the stored energy values looking for an energy that differs from the current energy by less than enerTol.
- For each such energy, use the readFile command to bring in the corresponding set of coordinates from the archive file.
- Use the archived and current coordinates of the subset of atoms for an rms comparison.
  - If the rms comparison result does not exceed rmsTol, restore the coordinates and return.
  - If the rms comparison result exceeds rmsTol, continue looping through the list in the energies variable.
If the loop completes, restore the coordinates, record the conformation using:
```
  writeFile archive file = $filename frame = 0
```
and add the current energy to the energies variable, using lappend.

At this point the procedure is complete.

In the main script (i.e., after the begin command), define the subset of heavy atoms and then loop over a number of conformations. You may want to generate conformations using a command such as:

  molGeom  set \
    torsion  relative  $incr \
    "*:1:C1"  "*:1:C2"  "*:1:C3"  "*:1:C4"

where $incr is the number of degrees by which you want to increment the torsion. For each conformation, call the compareArchive procedure that you defined at the start of your command input file.

Run your script and examine the resulting .arc file. Things to try:

Add a minimization before each call to the compareArchive procedure. How does this affect the resulting .arc file? (Use a new .arc filename for this run so that you start afresh.)
Add a restraint before looping over the conformations. How does this affect the resulting .arc file? (Use a new .arc filename for this run so that you start afresh.)

If you have trouble writing a script that does what you want, a sample file is supplied. It is the file named eng.inp and is one of those that you copied from $BIOSYM/tutorial/discover/BTCL.

Section 5

Section 5 of this tutorial makes use of poly_benz.inp, which illustrates additional capabilities of the molGeom command.

The script obtains the principal axis of a connected pair of benzene rings, finds a least-squares line through the atoms of the structure, and rotates the structure about that line.

Run the script and use the Insight program to view the initial and final conformations. Experiment with the effects of other keywords for the molGeom command.

Main access page

Advanced-Use (BTCL) access

BTCL - Tutorial Access.

Lesson on understanding & creating BTCL scripts, Introduction