VMD-L Mailing List

From: Vermaas, Josh (vermaasj_at_msu.edu)
Date: Sun Jul 23 2023 - 16:28:50 CDT

If you are just looking to get your ligands in the middle of your box, and you already have an XTC file handy (which has periodic box information in it), this is what I’d do in PBCtools-space.

package require pbctools

set pdbs {"A.pdb" "B.pdb"}

set xtcs {"A.xtc" "B.xtc"}



foreach pdb $pdbs xtc $xtcs {

#Making a new molecule from a PDB has the side effect of making a

#new frame from this PDB, which may or may not have PBC information in it.

mol new $pdb

#So, I’m going to delete the coordinate information

animate delete

mol addfile $xtc step 80 waitfor all

#Make whole the molecules that GROMACS breaks when it writes the XTC with join

pbc join residue -first 0 -last 0

#Unwrap the trajectory to eliminate jumps across a trajectory

pbc unwrap

#Re-wrap so that some selection you care about ends up in the middle.

#I recommend wrapping around one of the two ligands.

pbc wrap -all -compound residue -centersel "seltextgoeshere" -center com -all

animate write dcd "${pdb}_wrapped.dcd" waitfor all [molinfo top]

}



-Josh


From: <owner-vmd-l_at_ks.uiuc.edu> on behalf of "Ackad, Edward" <eackad_at_siue.edu>
Date: Sunday, July 23, 2023 at 11:19 AM
To: "vmd-l_at_ks.uiuc.edu" <vmd-l_at_ks.uiuc.edu>
Subject: Re: vmd-l: keeping two ligands in a box

Thanks Josh!
Can you suggest a way to avoid this or get the results I want a different way?
Thansk,
Eddie


_________________________________________________________
Edward Ackad, Ph.D<https://urldefense.com/v3/__http:/www.siue.edu/*7Eeackad__;JQ!!DZ3fjg!9w7mwSC6ggUWy_7PfV8Wr5yxVLVru571WCG5hDD0YbXteqC7dIRJLOOTpL6ptwNFS6qN1EZ5pI5d9Oihlw$>
Associate Professor of Physics
Computational Nanophotonics/Biophysics
Southern Illinois University Edwardsville
(618) 650-2390

________________________________
From: Vermaas, Josh <vermaasj_at_msu.edu>
Sent: Sunday, July 23, 2023 3:45 AM
To: Ackad, Edward <eackad_at_siue.edu>; vmd-l_at_ks.uiuc.edu <vmd-l_at_ks.uiuc.edu>
Subject: Re: vmd-l: keeping two ligands in a box


During the trajectory after you unwrap it, are there any timesteps when the position is >9999 or < -999? The error message you are getting is from one of the fixed-width columns in your pdb being overrun with a value that is too big or small. The position is probably the problem, as that is the only piece of a multiframe pdb that changes with time, and unwrapped positions can grow really fast. With an unwrapping command before the pdb writing command in the first loop, you’ll never make it to the second half of your script.



-Josh



From: <owner-vmd-l_at_ks.uiuc.edu> on behalf of "Ackad, Edward" <eackad_at_siue.edu>
Date: Saturday, July 22, 2023 at 11:37 AM
To: "vmd-l_at_ks.uiuc.edu" <vmd-l_at_ks.uiuc.edu>
Subject: vmd-l: keeping two ligands in a box



Hi all,

I've been having trouble getting a script to work for the following case. I have two compounds in different solvents. I want to see their interaction over time so I need to wrap them in the box. Here is the script I've been trying to get working but it complains about the pdb write or it ends up breaking all the bonds and completely distorting the compounds. I try writing the pdb twice but it seems if there are too many frames then the first pdb write gives:

Info) Opened coordinate file A.pdbA.pdb for writing.

PDB WRITE ERROR: Position, occupancy, or b-factor (beta) for atom 14

                 cannot be written in PDB format.

                 File will be truncated.

ERROR) write_timestep returned nonzero





Script:



set pdbs {"A.pdb" "B.pdb"}

set xtcs {"A.xtc" "B.xtc"}



# join the molecules and get rid of the solvent

for {set mymol 0} {$mymol < 2 } { incr mymol } {

    set nameis [lindex $pdbs $mymol]

    set namextc [lindex $xtcs $mymol]

    mol new "$nameis" type pdb

    mol addfile "$namextc" type xtc first 0 last -1 step 80 filebonds 1 autobonds 1 waitfor all

    puts "Joining: $nameis"

    pbc join residue -first 0 -last 0

    pbc unwrap

    set sel [atomselect top "residue 0 1"]

    $sel writepdb "${nameis}A.pdb"

    puts "Time to animate $nameis 1"

    animate write dcd "${nameis}.dcd" waitfor all sel $sel [molinfo top]

    puts "Time to savepdb $nameis 2"

    set sel2 [atomselect top "residue 0 1" frame 0]

    $sel2 writepdb "${nameis}.pdb"

    puts "done! $nameis"

    mol delete top

}



puts "Starting Wrap"

set outfile [open "status.txt" "w"]

for {set mymol 0} {$mymol < 2 } { incr mymol } {

    set nameis [lindex $pdbs $mymol]

    mol new "${nameis}A.pdb"

    mol addfile "$nameis.dcd" type dcd first 0 last -1 step 1 filebonds 1 autobonds 1 waitfor all

    puts "done loading $nameis dcd"

    set all [atomselect top "all"]

    set cellDims [measure minmax $all]



    set cellLengths [vecsub [lindex $cellDims 1] [lindex $cellDims 0]]

    pbc wrap -all -compound residue -centersel "all" -center com -all

    puts "doing animate"

    puts $outfile "doing animate"

    # Save the wrapped trajectory

    animate write dcd "${nameis}_wrapped.dcd" waitfor all [molinfo top]

    mol delete top

}

close $outfile



Thank you in advance!

_________________________________________________________
Edward Ackad, Ph.D<https://urldefense.com/v3/__http:/www.siue.edu/*7Eeackad__;JQ!!DZ3fjg!9LK1oTjo3CyP1NwARoA0Sn1HKwousyqZ6tEyUnIeBIuZvkgZs8U1OtQS0F6R12MmDowyXAlMjGtb1kNaNA$>
Associate Professor of Physics
Computational Nanophotonics/Biophysics
Southern Illinois University Edwardsville
(618) 650-2390