From: Axel Kohlmeyer (akohlmey_at_gmail.com)
Date: Tue Dec 13 2011 - 14:16:12 CST

rajan,

how about using this script (just whipped together in a few minutes)?

this can be run in textmode, but you have to install
a recent alpha test version of VMD or update the
topotools plugin with the latest version from here:

http://sites.google.com/site/akohlmey/software/topotools#TOC-Installation

you cannot have a selection that crosses multiple molecules,
but can easily combine multiple molecules, even if the atoms
overlap and *then* detect the overlapping atoms using "exwithin"
and make it remove entire molecules with "same fragment as".

HTH,
    axel.

=============================

#!/usr/local/bin/vmd -e

package require topotools 1.2

# minimum required distance between atoms
set mindist 2.5

# load both molecules
set solute [mol new TerthioPSFGEN50_D.psf waitfor all]
mol addfile MDRun50_D.pdb mol $solute waitfor all

set solvent [mol new 3cBoxCH2Cl2.psf waitfor all]
mol addfile 3cEqOPLSCH2Cl2.pdb mol $solvent waitfor all

set combined [TopoTools::mergemols "$solute $solvent"]

# transfer box dimensions, if available.
molinfo $combined set {a b c alpha beta gamma} [molinfo $solvent get
{a b c alpha beta gamma}

set sel [atomselect $solute all]
# define a macro that make it easy to determine the solute
atomselect macro solute "index < [$sel num]"
$sel delete

set sel [atomselect $combined "not (same fragment as (exwithin
$mindist of solute))"]
$sel writepsf merged.psf
$sel writepdb merged.pdb

# uncomment if run in text mode.
#after 1000 {quit}

====================================

On Tue, Dec 13, 2011 at 1:46 PM, Rajan Vatassery <rajan_at_umn.edu> wrote:
> Hi,
>        I have been trying to solvate a molecule with CH2Cl2 (the OPLS 3-site
> model), but I haven't been able to get the solvate (version 1.3) plugin
> to work. My problem is a bit easier than what the solvate plugin is
> built for so I am constrained to this plugin.
>        Specifically, I have a molecule that I have minimized in vacuum. I know
> the dimensions of this molecule (from minmax) so I created a solvent box
> that would encompass the molecule with the solvent padding that I
> wanted. Ideally I would like to merge the two pdb files together, and
> have a simple algorithm to remove close contacts then writepdb and
> writepsf. However, I cannot figure out a way to do this. Here's what I
> have tried:
> *****
> Solvate plugin command:
>
> solvate TerthioPSFGEN50_D.psf MDRun50_D.pdb -o solvateD -s CU -minmax
> {{-27.30 -27.30 -27.30} {27.3 27.32 27.32}} -x 0 -y 0 -z 0 +x 0 +y 0 +z
> 0 -b 2 -spdb 3cEqOPLSCH2Cl2.pdb -spsf 3cBoxCH2Cl2.psf -stop
> 3cEqOPLSCH2Cl2.top -ws 85.976 -ks {{type CH2}}
>
> My files are:
> Solute PSF: TerthioPSFGEN50_D.psf
> Solute PDB: MDRun50_D.pdb
> Solvent PDB: 3cEqOPLSCH2Cl2.pdb
> Solvent PSF: 3cBoxCH2Cl2.psf
> Solvent topology file: 3cEqOPLSCH2Cl2.top
> Side length of the cubic solvent box is: 85.976 A
> Each solvent molecule has a segid of QQQ and only one atom of type CH2
>
> This plugin does not seem to work with my solvent because it either does
> not encompass the solute with solvent as expected, and if I play around
> with the -ws or any of the other dimensions to make the
> solvent encompass the solute, I get an error telling me about a bug.
> This error comes up (I believe) exclusively when the solvent box is
> replicated in a periodic fashion. I have tried just about every
> combination of the numerical values in the command line, and nothing
> seems to work.
>
> *****
> Outside of the solvate.tcl script, I've tried a simple scripting
> attempt:
>
> It seems that the operative line in the solvate.tcl script is:
> set sel [atomselect top "segid QQQ and $keysel and same residue as \
>                   (x < $xmin or x > $xmax or \
>                    y < $ymin or y > $ymax or \
>                    z < $zmin or z > $zmax or \
>                    within $rwat of (not segid QQQ))"]
>
> so I tried to duplicate the "within" portion of this. What I realized is
> that I think it is impossible to get two different IDs to interact with
> each other.
>        I first loaded VMD with the command:
> vmd MDRun50_D.pdb TerthioPSFGEN50_D.psf
> and then I loaded the files 3cEqOPLSCH2Cl2.pdb and 3cBoxCH2Cl2.psf with
> the "new molecule" GUI tool. This gave me the solute (MDRun50_D.pdb with
> ID 0) and the solvent (3cEqOPLSCH2Cl2.pdb with ID 1) in the same VMD
> Main window.
>        Is there any way that I can write the command:
> set sel [atomselect 1 "within 2 of molid 0"]
> in some other fashion? In the solvate.tcl script, line:
> within $rwat of (not segid QQQ)
> seems to be doing exactly what I need, but because the two PDBs are in
> different IDs, I cannot use the same syntax. Am I stuck merging the
> files by hand and then doing this elimination myself? I'll post all the
> relevant files on
> http://www.tc.umn.edu/~rajan/VMD-l.html
> in a few minutes.
>
> Thanks for any help I can get.
>
> Rajan
>

-- 
Dr. Axel Kohlmeyer
akohlmey_at_gmail.com  http://goo.gl/1wk0
College of Science and Technology
Temple University, Philadelphia PA, USA.