Re: Problems with wrapAll

From: Vermaas, Joshua (Joshua.Vermaas_at_nrel.gov)
Date: Fri Jul 22 2016 - 11:55:30 CDT

Hi Sonu,

If this molecule is currently being simulated in vacuum, you could also try starting from a cubic box (lets say 40x40x40), and using pressure control to shrink the box down. I'm still not sure why you think the density means each molecule has that dimension though. I would rather imagine that there is a big mesh of these polymers, and that as a whole the mesh has a certain density, rather than individual chains being squeezed into a compact space.

-Josh

On 07/22/2016 10:37 AM, MEHRAN MB wrote:
Hi Sonu,

Basically, here are what Josh suggested:
* There are also female users on this forum so it would better to use terms like "Dear Sir or Madam".
* you need to increase the cutoffs length which seems it was a typo and you have already fixed that issue in your second email .CONF files.
* you cannot (as long as I know) have simulation box smaller than your molecule size. The simulation box, at least, should be equal to [ molecule_size + cuttoff length ]. This is the minimum size of the simulation box that you have to start with. Then, if you want to shrink your molecule, one way to do this is to do steering molecular dynamic using Collective Variable (colvar ) in NAMD. You can define radius of gyration as a variable and use harmonic restraint with moving centre to force your molecule to have target smaller radius of gyration. I have not done this before but it should be straight forward if you read Collective variable section in NAMD user guide. http://www.ks.uiuc.edu/Research/namd/2.10/ug/

here is what you should probably do:

add following two lines in your .CONF file

colvars on
colvarsConfig $name.in<http://name.in>

which $name is name of collective variable configuration file.
then in that file you should have something like

colvar {
   name Rg # arbitrary name I chose Rg
   gyration {
      atoms {
          "select your atom here. look at http://www.ks.uiuc.edu/Research/namd/2.10/ug/node56.html#SECTION000133100000000000000
      }
   }
}

# look at http://www.ks.uiuc.edu/Research/namd/2.10/ug/node58.html#SECTION000135300000000000000
harmonic {
   colvars Rg
   centers Ri # Ri your initial radius of gyration
   forceConstant 1.0
   targetCenters Rf # target radius of gyration
   targetNumSteps Ns # Number of steps for steering
}

hope it helps,

kind regards,
Mehran

On Fri, Jul 22, 2016 at 7:58 AM, Sonu Kumar <sonu2kd_at_gmail.com<mailto:sonu2kd_at_gmail.com>> wrote:

That was a mistake. Btw I don't quite get how to do what you suggest. Actually the molecule is larger than the pbc and its polymer.
I have written both the conf files. Can you suggest what I must add. I have further elaborated it down.

I have a polymer molecule of 1806 atoms with a density of 0.9.
The molecule has been minimized for 5 ns and has dimensions 38*32*26.

As per density I want to force it in a box of 25*25*25. However, upon wrapping it using pbctools and generating the pdb which is further used in the md simulation it is showing error probably because of bad contacts.

Can you please suggest me how to perform md simulation on this molecule by forcing it in a box smaller than this. It shall be helpful if you can just tell the part of script that I must add in my configuration file. I have pasted the two configuration files I am using. Please tell what I should add in either of these two to achieve what I want.

FIRST CONF FILE :

#My configuration file for minimization of the initial pdb :(please tell what I should add):

structure Poly.psf
coordinates Poly.pdb
outputName first

set temperature 298
firsttimestep 0

paraTypeCharmm on
parameters par_all22_prot.prm

temperature $temperature
exclude scaled1-4
1-4scaling 1.0
cutoff 12.0
switching on
switchdist 10.0
pairlistdist 14.0

# Integrator Parameters
timestep 2.0 ;# 2fs/step
rigidBonds all ;# needed for 2fs steps
nonbondedFreq 1
fullElectFrequency 2
stepspercycle 10

langevin on ;# do langevin dynamics
langevinDamping 1 ;# damping coefficient (gamma) of 5/ps
langevinTemp $temperature
langevinHydrogen no ;# don't couple langevin bath to hydrogens

restartfreq 500 ;# 500steps = every 1ps
dcdfreq 500
xstFreq 500
outputEnergies 100
outputPressure 100

if {1} {
minimize 1000
reinitvels $temperature
}

run 2500000

SECOND CONF FILE :(FOR PBC)

structure final.psf
coordinates minimized.pdb
outputName sample

set temperature 298

firsttimestep 0

margin 0

paraTypeCharmm on
parameters par_all22_prot.prm

temperature 0
reassignFreq 1000
reassignIncr 10
reassignHold 300

if {1} {
cellBasisVector1 24.954999923706055 0 0
cellBasisVector2 0 24.954999923706055 0
cellBasisVector3 0 0 24.954999923706055
cellOrigin 12.9974365234375 13.76833438873291 11.9735107421875
}
wrapWater on
wrapAll on

# Force-Field Parameters
exclude scaled1-4
1-4scaling 1.0
cutoff 15.0
switching on
switchdist 13.0
pairlistdist 17.5

# Integrator Parameters
timestep 1.0 ;# 2fs/step
rigidBonds none ;# needed for 2fs steps
#rigidTolerance 0.0005
nonbondedFreq 1
fullElectFrequency 2
stepspercycle 10

#PME (for full-system periodic electrostatics)
if {1} {
PME yes
PMEGridSpacing 1.0

#manual grid definition
PMEGridSizeX 32
PMEGridSizeY 32
PMEGridSizeZ 32
}

langevin off ;# do langevin dynamics
langevinDamping 1 ;# damping coefficient (gamma) of 5/ps
langevinTemp $temperature
langevinHydrogen no ;# don't couple langevin bath to hydrogens

restartfreq 500 ;# 500steps = every 1ps
dcdfreq 500
xstFreq 500
outputEnergies 100
outputPressure 100

if {1} {
minimize 1000
reinitvels $temperature
}

run 100000

This archive was generated by hypermail 2.1.6 : Sun Dec 31 2017 - 23:20:35 CST