* read ligand structures from MOL2 database and
* place in ASV Integrase binding site, followed
* by energy minimization and analysis
*

stream mmff_setup.str

!read CHARMM topology file that contains MMFF information
OPEN UNIT 10 CARD READ NAME top_all22_prot_mmff.inp
READ RTF CARD UNIT 10
CLOSe unit 10

!output data file
open unit 50 write form name test_cplx.dat

!MOL2 database of compounds in DOCKed orientations
open unit 60 read  form name mmff_db_2.mol2

!list of MOL2 structure numbers (corresponding
!to mmff_db_2.mol2) to extract from database
open unit 61 read  form name mmff_db_2.list

!mmff_db_2.mol2 contains all compounds of collowing numbers
!9169
!78502
!103840
!117022
!127625
!179171
!227189
!356465
!382139
!642649

set i 1

label gen_loop

! obtain compound number from mmff_db_2.list
get cmpd unit 61

!read compound from database and generate
read db @cmpd unit 60
generate

!read Integrase protein structure
open unit 62 read  form name asvin.mrk
read merck unit 62 append

!generate entire system
generate

update inbfrq -1 ihbfrq 0 cutnb 14. ctofnb 12. ctonnb 10. -
atom switch vatom vswitch cdie e14fac 0.75

!constrain all protein residues greater than 8 A
!from the protein
define active sele .byres. (ires 1 .around. 8.0) end

cons harm force 2.0 mass sele .not. active end

mini SD nstep 10 nprint 2

!determine interaction energy between protein and drug
!and output to database
inter sele ires 1 show end sele .not. ires 1 end

!obtain number of possible hbonds between drug and Integrase
define hetero sele .not. (type H* .or. type C* .or. type P* -
 .or. type S*) end

coor dist cut 4.0 sele (ires 1 .and. hetero) end -
sele ((.not. ires 1) .and. hetero) end

!write interaction energies and number of possible hbonds
!to data file
write title unit 50
* @cmpd ?ener ?elec ?vdw ?npair
*

!delete all atoms (and internal connectivity in
!charmm) prior to analyzing next compound
dele atom sele all end

calc i = @i + 1

if i le 5 goto gen_loop

stop