[gmx-users] comm-mode option for applied external electric field.

2019-02-01 Thread ARNAB MUKHERJEE
Hi,

I am trying to understand the comm-mode option in .mdp. I am simulating a
DNA protein system with electric field. Now, if I understand correctly the
comm-mode option, if I set it to linear, and comm-grps to System, it will
remove the velocity of COM velocity of the system. But here I want to
compute the velocity of COM of protein (which is positively charged and
moves in the direction of the applied field). Here I am using Langevin
dynamics. So my question here is comm-mode option, as it removes the
velocity of COM, does it also rescale the coordinates? If I run the same
simulation by using the option comm-mode = Linear, and for other case
comm-mode = None, in these 2 cases I should get different trajectories,
right?

Also, since for my case the protein and ions are driven by the field, I set
comm-grps = DNA_W, which is the group containing all DNA and water atoms
(here the DNA is position restrained). Is this the right thing to do here?
Or should I use comm-grps = none? I am confused here.

I would highly appreciate any help!

Thank you in advance.

Regards,

Arnab
-- 
Gromacs Users mailing list

* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting!

* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

* For (un)subscribe requests visit
https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a 
mail to gmx-users-requ...@gromacs.org.


[gmx-users] Problem in velocity plot while doing MD simulation of DNA protein with E field

2019-01-18 Thread ARNAB MUKHERJEE
Hi,

I am trying to study the motion of protein along a DNA with application of
electric field. I have position restrained the DNA, and the DNA axis is
along Z, and I am applying an electric field also along +Z. Now since the
protein is strongly positively charged it moves along the DNA in the
direction of the field (along +Z).

Now when I extract the z coordinate of the central C_alpha atom of the
protein, as expected I see for each cycle (since my box is periodic in all
directions) the z coordinate increases with time. But the problem is when I
extract and plot the z component of velocity of the central C_alpha atom of
protein, it fluctuates about 0. I don't understand how can the z component
of the velocity be even negative, when z coordinate shows positive slope
all the time.

In order to extract the position and velocities, I used the following
commands :

gmx_mpi traj -f traj_comp.xtc -s md_run.tpr -n index1.ndx -ox coor-CB.xv

gmx_mpi traj -f traj.trr -s md_run.tpr -n index1.ndx -ov vel-CB.xvg

I hope that the commands are correct.

Below I have pasted my .mdp file

title   = NVT equilibration with position restraint on all solute
(topology modified)
; Run parameters
integrator  = sd; leap-frog integrator
nsteps  = 250   ; 1 * 50 = 500 ps
;nsteps  = 5000
dt  = 0.02  ; 1 fs
; Output control
nstxout = 1000  ; save coordinates every 10 ps
nstvout = 1000  ; save velocities every 10 ps
nstcalcenergy   = 50
nstenergy   = 1000  ; save energies every 1 ps
nstxtcout   = 2500
;nstxout-compressed  = 5000   ; save compressed coordinates every 1.0 ps
 ; nstxout-compressed replaces nstxtcout
;compressed-x-grps  = System  ; replaces xtc-grps
nstlog  = 1000  ; update log file every 1 ps
; Bond parameters
continuation= no   ; first dynamics run
constraint_algorithm = lincs ; holonomic constraints
constraints = none  ; all bonds (even heavy atom-H bonds)
constrained
;lincs_iter = 2 ; accuracy of LINCS
lincs_order = 4 ; also related to accuracy
epsilon_r   = 15
; Neighborsearching
cutoff-scheme   = Verlet
ns_type = grid  ; search neighboring grid cels
nstlist = 10; 20 fs
rvdw_switch = 1.0
rlist   = 1.2   ; short-range neighborlist cutoff (in nm)
rcoulomb= 1.2   ; short-range electrostatic cutoff (in nm)
rvdw= 1.2   ; short-range van der Waals cutoff (in nm)
vdwtype = Cut-off   ; Twin range cut-offs rvdw >= rlist
;vdw-modifier= Force-switch
;Electrostatics
coulombtype = PME   ; Particle Mesh Ewald for long-range
electrostatics
pme_order   = 4 ; cubic interpolation
fourierspacing  = 0.12  ; grid spacing for FFT

ld-seed = -1

; Temperature coupling is on
tcoupl  = v-rescale
tc_grps = System
tau_t   = 2.0
ref_t   = 300

;energygrps = DNA Protein W ION

;freezegrps = Frozen-DNA-atoms
;freezedim = Y Y Y

E-x = 1 0 1 ;
E-y = 1 0 1 ;
E-z = 1 0.28 1 ;

; Pressure coupling is off
pcoupl  = no; no pressure coupling in NVT
; Periodic boundary conditions
pbc = xyz   ; 3-D PBC
; Dispersion correctiion
DispCorr= no; account for cut-off vdW scheme
; Velocity generation
gen_vel = yes   ; assign velocities from Maxwell
distribution
gen_temp= 300   ; temperature for Maxwell distribution
gen_seed= -1; generate a random seed
; COM motion removal
; These options remove motion of the protein/bilayer relative to the
solvent/ions
nstcomm = 50
comm-mode   = Linear
comm-grps   = DNA_W
;
refcoord_scaling = all

Earlier I had "coom-grps = System". So I was thinking that because it was
removing the velocity of COM of the system, hence the velocity of protein
that I was getting was only due to the thermal fluctuation, since it is
removing the velocity due to the electric field. Hence now I used
"coom-grps = DNA_W" where "DNA_W" is the group containing all the DNA and
water atoms, but still when I plot the z component of velocity of the
central C_alpha atom of protein, it again shows fluctuation about 0, I
don't understand why.

Can anyone please help me understand the problem here?

Thank you in advance,

Regards,

Arnab Mukherjee
-- 
Gromacs Users mailing list

* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting!

* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

* For (un)subscribe requests visit
https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a 
mail to gmx-users-requ...@gromacs.org.


Re: [gmx-users] Position restraints while performing pulling run for Umbrella sampling

2018-12-21 Thread ARNAB MUKHERJEE
Hello everyone,

I hope I have been able to write it clearly the problem I am facing. It
would be very helpful for me if someone can please guide me how to solve
this problem.

Thank you in advance,

Regards,

Arnab

On Thu, Dec 20, 2018 at 8:20 PM ARNAB MUKHERJEE 
wrote:

> Hi,
>
> I am intending to perform pulling run for umbrella sampling, on my coarse
> grained DNA protein system where I have a bundle containing 20 DNAs and
> proteins. I want to pull 1 DNA out of this bundle to perform Umbrella
> Sampling. I was going through this tutorial
> http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin/gmx-tutorials/umbrella/05_pull.html.
> So since I want to pull 1 DNA out of the bundle, I would like to position
> restrain the atoms of the other 19 DNAs. But as mentioned here
> http://www.gromacs.org/Documentation/How-tos/Position_Restraints, I need
> to add the position restraints under the molecule type of DNA, in my .top
> file, but then shouldn't it restrain the atoms of all the 20 DNAs? How can
> I restrain the atoms of other 19 DNAs, except the DNA that I am pulling?
>
> I would highly appreciate any help.
>
> Thank you very much,
>
> Regards,
>
> Arnab
>
-- 
Gromacs Users mailing list

* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting!

* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

* For (un)subscribe requests visit
https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a 
mail to gmx-users-requ...@gromacs.org.


[gmx-users] Position restraints while performing pulling run for Umbrella sampling

2018-12-20 Thread ARNAB MUKHERJEE
Hi,

I am intending to perform pulling run for umbrella sampling, on my coarse
grained DNA protein system where I have a bundle containing 20 DNAs and
proteins. I want to pull 1 DNA out of this bundle to perform Umbrella
Sampling. I was going through this tutorial
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin/gmx-tutorials/umbrella/05_pull.html.
So since I want to pull 1 DNA out of the bundle, I would like to position
restrain the atoms of the other 19 DNAs. But as mentioned here
http://www.gromacs.org/Documentation/How-tos/Position_Restraints, I need to
add the position restraints under the molecule type of DNA, in my .top
file, but then shouldn't it restrain the atoms of all the 20 DNAs? How can
I restrain the atoms of other 19 DNAs, except the DNA that I am pulling?

I would highly appreciate any help.

Thank you very much,

Regards,

Arnab
-- 
Gromacs Users mailing list

* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting!

* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

* For (un)subscribe requests visit
https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a 
mail to gmx-users-requ...@gromacs.org.


[gmx-users] Position restraints with Martini force field

2018-11-22 Thread ARNAB MUKHERJEE
Hi everyone,

I am simulating a DNA protein system with Martini forcefield. I want to use
position restraints on all DNA atoms. So I used gmx genrestr command to
generate .itp file for position restraints of all the DNA atoms. So without
including position restraints, my .top file looks like this :
 #include "martini_v2.1-dna.itp"
#include "martini_v2.0_ions.itp"

#define RUBBER_BANDS

#include "Nucleic_A+Nucleic_B.itp"
#include "Protein_A.itp"

[ system ]
; name
Martini system from DNA-CG.gro

[ molecules ]
; namenumber
Nucleic_A+Nucleic_B  1
Protein_A1
W   207523
NA   238
CL   21

Here all the .itp files are from Martini force field. It works fine without
problem.

Now when I include the position restraints, my .top file looks like the
following :

#include "martini_v2.1-dna.itp"
#include "martini_v2.0_ions.itp"

#define RUBBER_BANDS

#include "Nucleic_A+Nucleic_B.itp"
#include "posre_DNA.itp"
#include "Protein_A.itp"

[ system ]
; name
Martini system from DNA-CG.gro

[ molecules ]
; namenumber
Nucleic_A+Nucleic_B  1
Protein_A1
W   207523
NA   238
CL   21

So the only difference is that I included here the position restraints .itp
file generated using gmx genrestr (posre_DNA.itp) just below the DNA .itp
file (Nucleic_A+Nucleic_B.itp).

My doubt is that is this the right way to use position restrains? Is the
location in the .top file where I included  posre_DNA.itp file, correct?

I would highly appreciate any help.

Thank you very much,

Regards,

Arnab
-- 
Gromacs Users mailing list

* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting!

* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

* For (un)subscribe requests visit
https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a 
mail to gmx-users-requ...@gromacs.org.


Re: [gmx-users] Problem during equilibration using semi-isotropic pressure coupling

2018-10-11 Thread ARNAB MUKHERJEE
Hi Justin,

Thank you very much for your reply. I would try also position restrain
instead of freezing. But it is still strange to me why I got average
pressure of 1 bar after running only for 5 ns for the completely frozen
finite DNA when I did isotropic pressure coupling. But for the infinite
frozen DNA with semi isotropic pressure coupling it shows problem.

Thank you very much,

Regards,

Arnab
-- 
Gromacs Users mailing list

* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting!

* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

* For (un)subscribe requests visit
https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a 
mail to gmx-users-requ...@gromacs.org.


[gmx-users] Problem during equilibration using semi-isotropic pressure coupling

2018-10-09 Thread ARNAB MUKHERJEE
Hi everyone,

I have a system of infinite CG (Martini) DNA (length of the DNA same as
length of the periodic box) protein system. I am trying to equilibrate this
system to 1 bar pressure. Since my DNA is same as the size of the box, I
don't want the box to fluctuate along Z direction. Hence I am using semi
isotropic pressure coupling, with the box fluctuating only along XY.
I have built 2 systems. One where the DNA is completely frozen, and the 2nd
one only DNA end points are frozen. For both the systems I ran
equilibration NPT run for over 20 ns, and still the average pressure
doesn't show 1 bar, which is the target pressure. For the only end points
frozen DNA protein system, the average pressure shows 1.5 bar, and for the
completely frozen DNA system, the average pressure shows around 5.5 bar,
after running for over 20 ns for both the systems.
Here is how my .mdp file looks :

 title   = NVT equilibration with position restraint on all solute
(topology modified)
; Run parameters
integrator  = sd; leap-frog integrator
tinit   = 0
init-step   = 1600
nsteps  = 1200
dt  = 0.001 ; 1 fs
; Output control
nstxout = 0 ; save coordinates every 10 ps
nstvout = 0 ; save velocities every 10 ps
nstcalcenergy   = 50
nstenergy   = 1000  ; save energies every 1 ps
nstxtcout   = 2500
;nstxout-compressed  = 5000   ; save compressed coordinates every 1.0 ps
 ; nstxout-compressed replaces nstxtcout
;compressed-x-grps  = System  ; replaces xtc-grps
nstlog  = 1000  ; update log file every 1 ps
; Bond parameters
continuation= no   ; first dynamics run
constraint_algorithm = lincs ; holonomic constraints
constraints = none  ; all bonds (even heavy atom-H bonds)
constrained
;lincs_iter = 2 ; accuracy of LINCS
lincs_order = 4 ; also related to accuracy
epsilon_r   = 15
; Neighborsearching
cutoff-scheme   = Verlet
ns_type = grid  ; search neighboring grid cels
nstlist = 10; 20 fs
rvdw_switch = 1.0
rlist   = 1.2   ; short-range neighborlist cutoff (in nm)
rcoulomb= 1.2   ; short-range electrostatic cutoff (in nm)
rvdw= 1.2   ; short-range van der Waals cutoff (in nm)
vdwtype = Cut-off   ; Twin range cut-offs rvdw >= rlist
;vdw-modifier= Force-switch
;Electrostatics
coulombtype = PME   ; Particle Mesh Ewald for long-range electrostatics
pme_order   = 4 ; cubic interpolation
fourierspacing  = 0.12  ; grid spacing for FFT
; Temperature coupling is on
tcoupl  = v-rescale
tc_grps = System
tau_t   = 2.0
ld-seed  =  -1

;energygrps = DNA W_ION_Protein
;energygrp-excl = DNA DNA
freezegrps = DNA
freezedim = Y N Y

; Pressure coupling is off
;pcoupl = no; no pressure coupling in NVT
Pcoupl = parrinello-rahman
Pcoupltype  = semiisotropic
;tau_p   = 5.0
tau_p   = 2.5
compressibility = 3e-4 0
ref_p   = 1.0 1.0
; Periodic boundary conditions
pbc = xyz   ; 3-D PBC
; Dispersion correctiion
DispCorr= no; account for cut-off vdW scheme
; Velocity generation
gen_vel = no; assign velocities from Maxwell distribution
gen_temp= 300   ; temperature for Maxwell distribution
gen_seed= -1; generate a random seed
; COM motion removal
; These options remove motion of the protein/bilayer relative to the
solvent/ions
nstcomm = 50
comm-mode   = Linear
comm-grps   = System
;
;refcoord_scaling = com
refcoord_scaling = all

As you see I used Parrinello Rahman Pressure coupling. Initially I had used
tau_p = 5.0 ps, and then I reduced it to 2.5 ps but there is still not much
change. I also did equilibration run for Finite DNA system with isotropic
pressure coupling and using tau_p = 5.0 ps, I got an average pressure of 1
bar just after 5ns NPT run. While calculating the average pressure for both
these systems, the RMSD shows around 10 bar.

For the infinite DNA (completely frozen DNA system) I have tried using
Berensden thermostat with tau_p = 1ps, and after 4 ns run, I have average
pressure around 7 bar. I am continuing this run to check if it shows better
results.

I do not understand what is problem occuring while equlibrating for the
infinite DNA systems using semi-isotropic pressure coupling? Why it cannot
reach the target pressure?

I would highly appreciate any help.

Thank you very much,

Regards,

Arnab
-- 
Gromacs Users mailing list

* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting!

* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

* For (un)subscribe requests visit
https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a 
mail to gmx-users-requ...@gromacs.org.


Re: [gmx-users] Problem while equilibrating system with frozen groups

2018-09-11 Thread ARNAB MUKHERJEE
Hi Peter,

Thank you very much for your suggestion. I will try using atom restraints
and check if the pressure coupling works properly.

Thanks and regards,

Arnab.
On Sep 11, 2018 10:33 AM, "Peter Kroon"  wrote:

> Hey Arnab,
>
>
> Why do you freeze atoms instead of putting position restraints? Why not
> let the loose entirely? Freezing creates large artefacts.
>
> Second, 4 ns equilibration is *way* too short.
>
> Lastly, see http://www.gromacs.org/Documentation/Terminology/Pressure
>
>
> Peter
>
>
> On 10-09-18 19:01, ARNAB MUKHERJEE wrote:
> > Hi,
> >
> > I am simulating DNA (infinite) protein Martini system with electric
> field.
> > I have 2 systems, for 1 complete DNA is frozen and for the 2nd one, only
> 4
> > atoms of the DNA is frozen. When I perform NPT equilibration for the
> > complete DNA frozen system, it shows me this error :
> >
> >  Fatal error:
> > 4442 particles communicated to PME rank 24 are more than 2/3 times the
> > cut-off out of the domain decomposition cell of their charge group in
> > dimension y.
> > This usually means that your system is not well equilibrated.
> >
> > But it runs fine for the other system with 4 DNA atoms frozen. For the
> > completely frozen DNA, if I freeze only 2 dimensions, letting 1 dimension
> > free it runs without complaining. I do not understand why it shows
> problem
> > for the completely frozen DNA with 3 dimensions frozen? Since the DNA is
> > infinite I use semi-isotropic pressure coupling keeping the
> compressibility
> > along Z 0. So the pressure coupling options in my .mdp file looks like
> this
> > :
> >
> > freezegrps = DNA
> > freezedim = Y N Y
> >
> > ; Pressure coupling is off
> > ;pcoupl = no; no pressure coupling in NVT
> > Pcoupl = parrinello-rahman
> > Pcoupltype  = semiisotropic
> > tau_p   = 5.0
> > compressibility = 3e-4 0
> > ref_p   = 1.0 1.0
> > refcoord_scaling = all
> >
> > I tried it for finite DNA too and it shows same problem for completely
> > frozen DNA.
> >
> > I ran the completely Frozen DNA system keeping 1 dimension free, and
> > freezing the other 2 dimensions. I ran the equilibration run for 500 ps.
> I
> > plotted the pressure and the average pressure is 9.15 bar, while the
> target
> > pressure that I had put in the .mdp file is 1 bar. So there is a large
> > difference. I was searching regarding this problem, and I found that I
> need
> > to use "energygrp-excl" for frozen groups. So I ran the completely frozen
> > DNA system with  energygrp-excl = DNA, using cutoff-scheme = group (since
> > energygrp-excl is not supported in Verlet in my gromacs version 5.0.6),
> and
> > the average pressure I get is again 9.28 bar. For the system with 4 DNA
> > atoms frozen, the average pressure shows 5.17 bar, which is still
> different
> > from the target pressure.
> >
> > I am now running these equilibration runs for longer (4 ns), to see if
> that
> > helps to get the average pressure closer to 1 bar. Normally when all
> atoms
> > are free what I found was 0.5 ns was enough to get average pressure 1
> bar,
> > and equilibrate the system. But I don't understand what is the problem I
> am
> > facing here with frozen atoms?
> >
> > I turn off the electric field during the equilibration runs.
> >
> > I would highly appreciate any help!
> >
> > Thank you very much,
> >
> > Regards,
> >
> > Arnab
>
>
>
> --
> Gromacs Users mailing list
>
> * Please search the archive at http://www.gromacs.org/
> Support/Mailing_Lists/GMX-Users_List before posting!
>
> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>
> * For (un)subscribe requests visit
> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
> send a mail to gmx-users-requ...@gromacs.org.
>
-- 
Gromacs Users mailing list

* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting!

* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

* For (un)subscribe requests visit
https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a 
mail to gmx-users-requ...@gromacs.org.


[gmx-users] Problem while equilibrating system with frozen groups

2018-09-10 Thread ARNAB MUKHERJEE
Hi,

I am simulating DNA (infinite) protein Martini system with electric field.
I have 2 systems, for 1 complete DNA is frozen and for the 2nd one, only 4
atoms of the DNA is frozen. When I perform NPT equilibration for the
complete DNA frozen system, it shows me this error :

 Fatal error:
4442 particles communicated to PME rank 24 are more than 2/3 times the
cut-off out of the domain decomposition cell of their charge group in
dimension y.
This usually means that your system is not well equilibrated.

But it runs fine for the other system with 4 DNA atoms frozen. For the
completely frozen DNA, if I freeze only 2 dimensions, letting 1 dimension
free it runs without complaining. I do not understand why it shows problem
for the completely frozen DNA with 3 dimensions frozen? Since the DNA is
infinite I use semi-isotropic pressure coupling keeping the compressibility
along Z 0. So the pressure coupling options in my .mdp file looks like this
:

freezegrps = DNA
freezedim = Y N Y

; Pressure coupling is off
;pcoupl = no; no pressure coupling in NVT
Pcoupl = parrinello-rahman
Pcoupltype  = semiisotropic
tau_p   = 5.0
compressibility = 3e-4 0
ref_p   = 1.0 1.0
refcoord_scaling = all

I tried it for finite DNA too and it shows same problem for completely
frozen DNA.

I ran the completely Frozen DNA system keeping 1 dimension free, and
freezing the other 2 dimensions. I ran the equilibration run for 500 ps. I
plotted the pressure and the average pressure is 9.15 bar, while the target
pressure that I had put in the .mdp file is 1 bar. So there is a large
difference. I was searching regarding this problem, and I found that I need
to use "energygrp-excl" for frozen groups. So I ran the completely frozen
DNA system with  energygrp-excl = DNA, using cutoff-scheme = group (since
energygrp-excl is not supported in Verlet in my gromacs version 5.0.6), and
the average pressure I get is again 9.28 bar. For the system with 4 DNA
atoms frozen, the average pressure shows 5.17 bar, which is still different
from the target pressure.

I am now running these equilibration runs for longer (4 ns), to see if that
helps to get the average pressure closer to 1 bar. Normally when all atoms
are free what I found was 0.5 ns was enough to get average pressure 1 bar,
and equilibrate the system. But I don't understand what is the problem I am
facing here with frozen atoms?

I turn off the electric field during the equilibration runs.

I would highly appreciate any help!

Thank you very much,

Regards,

Arnab
-- 
Gromacs Users mailing list

* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting!

* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

* For (un)subscribe requests visit
https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a 
mail to gmx-users-requ...@gromacs.org.


Re: [gmx-users] Cannot compute velocity of COM of a group of atoms

2018-08-08 Thread ARNAB MUKHERJEE
Hi,

Thank you very much for your suggestion.

By "nsteps" I meant the # of step interval between which I would record the
velocity. So my simulation run has 10 million steps (200 ns) and I chose
nstvout = 1000.

But I was wondering is it possible to do a mdrun -rerun (changing the .mdp
file), and compute the velocity for the old trajectory (.trr) file, rather
than running a new simulation again?

Thanks and Regards,

Arnab

On Wed, Aug 8, 2018 at 4:15 PM, Justin Lemkul  wrote:

>
>
> On 8/7/18 1:47 PM, ARNAB MUKHERJEE wrote:
>
>> Dear all,
>>
>> I have simulated a system of DNA and Protein, and I want to calculate the
>> velocity of the center of mass of protein as a function of time. So I used
>> the following command :
>>
>> gmx traj -f traj_comp.trr -s md_run-E-Field.tpr -n index.ndx -ov
>> test-vel.xvg -com
>>
>> I have pasted below the file that it generates :
>>
>> # GROMACS:  gmx traj, VERSION 5.0.6
>> # Executable:   /applic/applications/GROMACS/5.0.6/bin/gmx
>> # Library dir:  /applic/applications/GROMACS/5.0.6/share/gromacs/top
>> # Command line:
>> #   gmx traj -f traj_comp.trr -s md_run-E-Field.tpr -n index.ndx -ov
>> test-vel.xvg -com
>> # gmx is part of G R O M A C S:
>> #
>> # Georgetown Riga Oslo Madrid Amsterdam Chisinau Stockholm
>> #
>> @title "Center of mass velocity"
>> @xaxis  label "Time (ps)"
>> @yaxis  label "Velocity (nm/ps)"
>> @TYPE xy
>> @ view 0.15, 0.15, 0.75, 0.85
>> @ legend on
>> @ legend box on
>> @ legend loctype view
>> @ legend 0.78, 0.8
>> @ legend length 2
>> @ s0 legend "Protein X"
>> @ s1 legend "Protein Y"
>> @ s2 legend "Protein Z"
>> "test-vel.xvg" 24L, 731C  24,1
>> Bot
>>
>> So it doesn't have the velocities. Is it due to the fact that in my .mdp
>> file, I had set nstvout  = 0 ? I have pasted below the initial part of my
>> .mdp file
>>
>> title   = NVT equilibration with position restraint on all solute
>> (topology modified)
>> ; Run parameters
>> integrator  = md; leap-frog integrator
>> nsteps  = 500   ; 1 * 50 = 500 ps
>> ;nsteps  = 5000
>> dt  = 0.02  ; 1 fs
>> ; Output control
>> nstxout = 0 ; save coordinates every 10 ps
>> nstvout = 0 ; save velocities every 10 ps
>> nstcalcenergy   = 50
>> nstenergy   = 1000  ; save energies every 1 ps
>> nstxtcout   = 2500
>> ;nstxout-compressed  = 5000   ; save compressed coordinates every 1.0 ps
>>   ; nstxout-compressed replaces nstxtcout
>> ;compressed-x-grps  = System  ; replaces xtc-grps
>> nstlog  = 1000  ; update log file every 1 ps
>> ; Bond parameters
>> continuation= no   ; first dynamics run
>> constraint_algorithm = lincs ; holonomic constraints
>> constraints = none  ; all bonds (even heavy atom-H bonds)
>> constrained
>> ;lincs_iter = 2 ; accuracy of LINCS
>> lincs_order = 4 ; also related to accuracy
>>
>> Is there a way I can still compute the velocities, or do I need to run the
>> simulation again, putting nstvout = nsteps ?
>>
>
> Yes, you need a non-zero value for nstvout, but you shouldn't set it to
> the same value as nsteps, because then you'll only get one frame. You need
> some sensible output frequency that allows you to collect relevant
> statistics of the quantity of interest.
>
> -Justin
>
> --
> ==
>
> Justin A. Lemkul, Ph.D.
> Assistant Professor
> Virginia Tech Department of Biochemistry
>
> 303 Engel Hall
> 340 West Campus Dr.
> Blacksburg, VA 24061
>
> jalem...@vt.edu | (540) 231-3129
> http://www.thelemkullab.com
>
> ==
>
> --
> Gromacs Users mailing list
>
> * Please search the archive at http://www.gromacs.org/Support
> /Mailing_Lists/GMX-Users_List before posting!
>
> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>
> * For (un)subscribe requests visit
> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
> send a mail to gmx-users-requ...@gromacs.org.
>
-- 
Gromacs Users mailing list

* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting!

* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

* For (un)subscribe requests visit
https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a 
mail to gmx-users-requ...@gromacs.org.


[gmx-users] Cannot compute velocity of COM of a group of atoms

2018-08-07 Thread ARNAB MUKHERJEE
Dear all,

I have simulated a system of DNA and Protein, and I want to calculate the
velocity of the center of mass of protein as a function of time. So I used
the following command :

gmx traj -f traj_comp.trr -s md_run-E-Field.tpr -n index.ndx -ov
test-vel.xvg -com

I have pasted below the file that it generates :

# GROMACS:  gmx traj, VERSION 5.0.6
# Executable:   /applic/applications/GROMACS/5.0.6/bin/gmx
# Library dir:  /applic/applications/GROMACS/5.0.6/share/gromacs/top
# Command line:
#   gmx traj -f traj_comp.trr -s md_run-E-Field.tpr -n index.ndx -ov
test-vel.xvg -com
# gmx is part of G R O M A C S:
#
# Georgetown Riga Oslo Madrid Amsterdam Chisinau Stockholm
#
@title "Center of mass velocity"
@xaxis  label "Time (ps)"
@yaxis  label "Velocity (nm/ps)"
@TYPE xy
@ view 0.15, 0.15, 0.75, 0.85
@ legend on
@ legend box on
@ legend loctype view
@ legend 0.78, 0.8
@ legend length 2
@ s0 legend "Protein X"
@ s1 legend "Protein Y"
@ s2 legend "Protein Z"
"test-vel.xvg" 24L, 731C  24,1
Bot

So it doesn't have the velocities. Is it due to the fact that in my .mdp
file, I had set nstvout  = 0 ? I have pasted below the initial part of my
.mdp file

title   = NVT equilibration with position restraint on all solute
(topology modified)
; Run parameters
integrator  = md; leap-frog integrator
nsteps  = 500   ; 1 * 50 = 500 ps
;nsteps  = 5000
dt  = 0.02  ; 1 fs
; Output control
nstxout = 0 ; save coordinates every 10 ps
nstvout = 0 ; save velocities every 10 ps
nstcalcenergy   = 50
nstenergy   = 1000  ; save energies every 1 ps
nstxtcout   = 2500
;nstxout-compressed  = 5000   ; save compressed coordinates every 1.0 ps
 ; nstxout-compressed replaces nstxtcout
;compressed-x-grps  = System  ; replaces xtc-grps
nstlog  = 1000  ; update log file every 1 ps
; Bond parameters
continuation= no   ; first dynamics run
constraint_algorithm = lincs ; holonomic constraints
constraints = none  ; all bonds (even heavy atom-H bonds)
constrained
;lincs_iter = 2 ; accuracy of LINCS
lincs_order = 4 ; also related to accuracy

Is there a way I can still compute the velocities, or do I need to run the
simulation again, putting nstvout = nsteps ?

I would highly appreciate any help.

Thank you very much,

Regards,

Arnab
-- 
Gromacs Users mailing list

* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting!

* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

* For (un)subscribe requests visit
https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a 
mail to gmx-users-requ...@gromacs.org.


[gmx-users] Error : symtab get_symtab_handle 367 not found

2018-06-19 Thread ARNAB MUKHERJEE
Hi,

I am simulating a Martini coarse grained DNA protamine system. I have 1
infinite DNA + 1 protamine. Since while using the martini python script for
coarse graining it removes the 2 phosphate atoms of the terminal base
pairs, so in order to have the right infinite DNA system, while building
it, I had build 2 extra base pairs, which I removed them manually and
renumbered the topology .itp file atom numbers to match with the .gro file.
So my periodic cuboidal box size in Z direction (the DNA is aligned along
Z) is # of BPs*3.38 Angstrom, so that the pbc along Z makes the DNA
infinite.

The problem is when I simulate this system, during the short NPT
equilibration, it shows this error :

  Fatal error:
symtab get_symtab_handle 367 not found

The strange thing is if I submit a short run, in interactive mode for eg,
it runs fine, as long as it stays connected. But when I submit the complete
run in the cluster, it shows this error, before even starting Step 0. I
tried to google about this error, but didn't find much info. I have been
running this simulation in version 5.0.6. I checked also using newer
version of gromacs 5.1.4, it shows running status, but in the .log file it
shows this :

 Started mdrun on rank 0 Thu Jun 14 01:29:37 2018
   Step   Time Lambda
  00.00.0


Not all bonded interactions have been properly assigned to the domain
decomposition cells

Are the 2 different errors that I get in the different versions, connected?

To test if there is a problem in the force field (.itp file) since I had to
modify it manually to renumber the atoms, I ran a finite 50 BP DNA with
modified .itp file for the 50 BP DNA, and it runs fine.

I am not able to understand what is the problem. I am pasting my input .mdp
parameters that I used for the run. I used semiisotropic pressure coupling
as I want to keep the Z dimension of the box constant. I have also frozen 4
atoms in the 2 ends of the DNA since I want to keep the DNA aligned along
Z, since I later want to apply an E field along Z direction.

 title   = NVT equilibration with position restraint on all solute
(topology modified)
; Run parameters
integrator  = md; leap-frog integrator
;nsteps = 3000  ; 1 * 50 = 500 ps
nsteps  = 50
dt  = 0.001 ; 1 fs
; Output control
nstxout = 0 ; save coordinates every 10 ps
nstvout = 0 ; save velocities every 10 ps
nstcalcenergy   = 50
nstenergy   = 1000  ; save energies every 1 ps
nstxtcout   = 2500
;nstxout-compressed  = 5000   ; save compressed coordinates every 1.0 ps
 ; nstxout-compressed replaces nstxtcout
;compressed-x-grps  = System  ; replaces xtc-grps
nstlog  = 1000  ; update log file every 1 ps
; Bond parameters
continuation= no   ; first dynamics run
constraint_algorithm = lincs ; holonomic constraints
constraints = none  ; all bonds (even heavy atom-H bonds)
constrained
;lincs_iter = 2 ; accuracy of LINCS
lincs_order = 4 ; also related to accuracy
epsilon_r   = 15
; Neighborsearching
cutoff-scheme   = Verlet
ns_type = grid  ; search neighboring grid cels
nstlist = 10; 20 fs
rvdw_switch = 1.0
rlist   = 1.2   ; short-range neighborlist cutoff (in nm)
rcoulomb= 1.2   ; short-range electrostatic cutoff (in nm)
rvdw= 1.2   ; short-range van der Waals cutoff (in nm)
vdwtype = Cut-off   ; Twin range cut-offs rvdw >= rlist
;vdw-modifier= Force-switch
;Electrostatics
coulombtype = PME   ; Particle Mesh Ewald for long-range electrostatics
pme_order   = 4 ; cubic interpolation
fourierspacing  = 0.12  ; grid spacing for FFT
; Temperature coupling is on
tcoupl  = v-rescale
tc_grps = System
tau_t   = 1.0
ref_t   = 300

;energygrps = DNA W_ION_Protein
;energygrp-excl = DNA DNA
freezegrps = DNA-Frozen-Atoms
freezedim = Y Y Y

; Pressure coupling is off
;pcoupl = no; no pressure coupling in NVT
Pcoupl = parrinello-rahman
Pcoupltype  = semiisotropic
tau_p   = 5.0
compressibility = 3e-4 0
ref_p   = 1.0 1.0
; Periodic boundary conditions
pbc = xyz   ; 3-D PBC
; Dispersion correctiion
DispCorr= no; account for cut-off vdW scheme
; Velocity generation
gen_vel = yes   ; assign velocities from Maxwell distribution
gen_temp= 300   ; temperature for Maxwell distribution
gen_seed= -1; generate a random seed
; COM motion removal
; These options remove motion of the protein/bilayer relative to the
solvent/ions
nstcomm = 50
comm-mode   = Linear
comm-grps   = System
;
refcoord_scaling = com
;refcoord_scaling = all

I would highly appreciate any help.

Thank you in advance,

Regards,

Arnab Mukherjee
-- 
Gromacs Users mailing list

* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting!

* Can't p

[gmx-users] Error : symtab get_symtab_handle 367 not found

2018-06-18 Thread ARNAB MUKHERJEE
Hi,

I am simulating a Martini coarse grained DNA protamine system. I have 1
infinite DNA + 1 protamine. Since while using the martini python script for
coarse graining it removes the 2 phosphate atoms of the terminal base
pairs, so in order to have the right infinite DNA system, while building
it, I had build 2 extra base pairs, which I removed them manually and
renumbered the topology .itp file atom numbers to match with the .gro file.
So my periodic cuboidal box size in Z direction (the DNA is aligned along
Z) is # of BPs*3.38 Angstrom, so that the pbc along Z makes the DNA
infinite.

The problem is when I simulate this system, during the short NPT
equilibration, it shows this error :

  Fatal error:
symtab get_symtab_handle 367 not found

The strange thing is if I submit a short run, in interactive mode for eg,
it runs fine, as long as it stays connected. But when I submit the complete
run in the cluster, it shows this error, before even starting Step 0. I
tried to google about this error, but didn't find much info. I have been
running this simulation in version 5.0.6. I checked also using newer
version of gromacs 5.1.4, it shows running status, but in the .log file it
shows this :

 Started mdrun on rank 0 Thu Jun 14 01:29:37 2018
   Step   Time Lambda
  00.00.0


Not all bonded interactions have been properly assigned to the domain
decomposition cells

Are the 2 different errors that I get in the different versions, connected?

To test if there is a problem in the force field (.itp file) since I had to
modify it manually to renumber the atoms, I ran a finite 50 BP DNA with
modified .itp file for the 50 BP DNA, and it runs fine.

I am not able to understand what is the problem. I am pasting my input .mdp
parameters that I used for the run. I used semiisotropic pressure coupling
as I want to keep the Z dimension of the box constant. I have also frozen 4
atoms in the 2 ends of the DNA since I want to keep the DNA aligned along
Z, since I later want to apply an E field along Z direction.

 title   = NVT equilibration with position restraint on all solute
(topology modified)
; Run parameters
integrator  = md; leap-frog integrator
;nsteps = 3000  ; 1 * 50 = 500 ps
nsteps  = 50
dt  = 0.001 ; 1 fs
; Output control
nstxout = 0 ; save coordinates every 10 ps
nstvout = 0 ; save velocities every 10 ps
nstcalcenergy   = 50
nstenergy   = 1000  ; save energies every 1 ps
nstxtcout   = 2500
;nstxout-compressed  = 5000   ; save compressed coordinates every 1.0 ps
 ; nstxout-compressed replaces nstxtcout
;compressed-x-grps  = System  ; replaces xtc-grps
nstlog  = 1000  ; update log file every 1 ps
; Bond parameters
continuation= no   ; first dynamics run
constraint_algorithm = lincs ; holonomic constraints
constraints = none  ; all bonds (even heavy atom-H bonds)
constrained
;lincs_iter = 2 ; accuracy of LINCS
lincs_order = 4 ; also related to accuracy
epsilon_r   = 15
; Neighborsearching
cutoff-scheme   = Verlet
ns_type = grid  ; search neighboring grid cels
nstlist = 10; 20 fs
rvdw_switch = 1.0
rlist   = 1.2   ; short-range neighborlist cutoff (in nm)
rcoulomb= 1.2   ; short-range electrostatic cutoff (in nm)
rvdw= 1.2   ; short-range van der Waals cutoff (in nm)
vdwtype = Cut-off   ; Twin range cut-offs rvdw >= rlist
;vdw-modifier= Force-switch
;Electrostatics
coulombtype = PME   ; Particle Mesh Ewald for long-range electrostatics
pme_order   = 4 ; cubic interpolation
fourierspacing  = 0.12  ; grid spacing for FFT
; Temperature coupling is on
tcoupl  = v-rescale
tc_grps = System
tau_t   = 1.0
ref_t   = 300

;energygrps = DNA W_ION_Protein
;energygrp-excl = DNA DNA
freezegrps = DNA-Frozen-Atoms
freezedim = Y Y Y

; Pressure coupling is off
;pcoupl = no; no pressure coupling in NVT
Pcoupl = parrinello-rahman
Pcoupltype  = semiisotropic
tau_p   = 5.0
compressibility = 3e-4 0
ref_p   = 1.0 1.0
; Periodic boundary conditions
pbc = xyz   ; 3-D PBC
; Dispersion correctiion
DispCorr= no; account for cut-off vdW scheme
; Velocity generation
gen_vel = yes   ; assign velocities from Maxwell distribution
gen_temp= 300   ; temperature for Maxwell distribution
gen_seed= -1; generate a random seed
; COM motion removal
; These options remove motion of the protein/bilayer relative to the
solvent/ions
nstcomm = 50
comm-mode   = Linear
comm-grps   = System
;
refcoord_scaling = com
;refcoord_scaling = all

I would highly appreciate any help.

Thank you in advance,

Regards,

Arnab Mukherjee
-- 
Gromacs Users mailing list

* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting!

* Can't p

[gmx-users] Error : symtab get_symtab_handle 367 not found

2018-06-14 Thread ARNAB MUKHERJEE
Hi,

I am simulating a Martini coarse grained DNA protamine system. I have 1
infinite DNA + 1 protamine. Since while using the martini python script for
coarse graining it removes the 2 phosphate atoms of the terminal base
pairs, so in order to have the right infinite DNA system, while building
it, I had build 2 extra base pairs, which I removed them manually and
renumbered the topology .itp file atom numbers to match with the .gro file.
So my periodic cuboidal box size in Z direction (the DNA is aligned along
Z) is # of BPs*3.38 Angstrom, so that the pbc along Z makes the DNA
infinite.

The problem is when I simulate this system, during the short NPT
equilibration, it shows this error :

  Fatal error:
symtab get_symtab_handle 367 not found

The strange thing is if I submit a short run, in interactive mode for eg,
it runs fine, as long as it stays connected. But when I submit the complete
run in the cluster, it shows this error, before even starting Step 0. I
tried to google about this error, but didn't find much info. I have been
running this simulation in version 5.0.6. I checked also using newer
version of gromacs 5.1.4, it shows running status, but in the .log file it
shows this :

 Started mdrun on rank 0 Thu Jun 14 01:29:37 2018
   Step   Time Lambda
  00.00.0


Not all bonded interactions have been properly assigned to the domain
decomposition cells

Are the 2 different errors that I get in the different versions, connected?

To test if there is a problem in the force field (.itp file) since I had to
modify it manually to renumber the atoms, I ran a finite 50 BP DNA with
modified .itp file for the 50 BP DNA, and it runs fine.

I am not able to understand what is the problem. I am pasting my input .mdp
parameters that I used for the run. I used semiisotropic pressure coupling
as I want to keep the Z dimension of the box constant. I have also frozen 4
atoms in the 2 ends of the DNA since I want to keep the DNA aligned along
Z, since I later want to apply an E field along Z direction.

 title   = NVT equilibration with position restraint on all solute
(topology modified)
; Run parameters
integrator  = md; leap-frog integrator
;nsteps = 3000  ; 1 * 50 = 500 ps
nsteps  = 50
dt  = 0.001 ; 1 fs
; Output control
nstxout = 0 ; save coordinates every 10 ps
nstvout = 0 ; save velocities every 10 ps
nstcalcenergy   = 50
nstenergy   = 1000  ; save energies every 1 ps
nstxtcout   = 2500
;nstxout-compressed  = 5000   ; save compressed coordinates every 1.0 ps
 ; nstxout-compressed replaces nstxtcout
;compressed-x-grps  = System  ; replaces xtc-grps
nstlog  = 1000  ; update log file every 1 ps
; Bond parameters
continuation= no   ; first dynamics run
constraint_algorithm = lincs ; holonomic constraints
constraints = none  ; all bonds (even heavy atom-H bonds)
constrained
;lincs_iter = 2 ; accuracy of LINCS
lincs_order = 4 ; also related to accuracy
epsilon_r   = 15
; Neighborsearching
cutoff-scheme   = Verlet
ns_type = grid  ; search neighboring grid cels
nstlist = 10; 20 fs
rvdw_switch = 1.0
rlist   = 1.2   ; short-range neighborlist cutoff (in nm)
rcoulomb= 1.2   ; short-range electrostatic cutoff (in nm)
rvdw= 1.2   ; short-range van der Waals cutoff (in nm)
vdwtype = Cut-off   ; Twin range cut-offs rvdw >= rlist
;vdw-modifier= Force-switch
;Electrostatics
coulombtype = PME   ; Particle Mesh Ewald for long-range electrostatics
pme_order   = 4 ; cubic interpolation
fourierspacing  = 0.12  ; grid spacing for FFT
; Temperature coupling is on
tcoupl  = v-rescale
tc_grps = System
tau_t   = 1.0
ref_t   = 300

;energygrps = DNA W_ION_Protein
;energygrp-excl = DNA DNA
freezegrps = DNA-Frozen-Atoms
freezedim = Y Y Y

; Pressure coupling is off
;pcoupl = no; no pressure coupling in NVT
Pcoupl = parrinello-rahman
Pcoupltype  = semiisotropic
tau_p   = 5.0
compressibility = 3e-4 0
ref_p   = 1.0 1.0
; Periodic boundary conditions
pbc = xyz   ; 3-D PBC
; Dispersion correctiion
DispCorr= no; account for cut-off vdW scheme
; Velocity generation
gen_vel = yes   ; assign velocities from Maxwell distribution
gen_temp= 300   ; temperature for Maxwell distribution
gen_seed= -1; generate a random seed
; COM motion removal
; These options remove motion of the protein/bilayer relative to the
solvent/ions
nstcomm = 50
comm-mode   = Linear
comm-grps   = System
;
refcoord_scaling = com
;refcoord_scaling = all

I would highly appreciate any help.

Thank you in advance,

Regards,

Arnab Mukherjee
-- 
Gromacs Users mailing list

* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting!

* Can't p

[gmx-users] Problem to compute interaction energy between different groups

2017-11-24 Thread ARNAB MUKHERJEE
Hello,

I have simulated a system of DNA and Protamine (Coarse grained, Martini). I
want to look at the interaction energy between different groups like DNA -
Protamine, DNA - ions, etc. I have defined energy groups in the mdp file,
and I pass the index file with the groups while building the tpr file.
After the run is complete, when I run gmx energy command, no matter what
interaction energy I choose, it always gives me 0 (for all time steps). But
the total Potential energy works. I am not able to understand why?

This is how my .mdp file looks like :

title   = NVT equilibration with position restraint on all solute
(topology modified)
; Run parameters
integrator  = md; leap-frog integrator
nsteps  = 6000  ; 1 * 50 = 500 ps
;nsteps  = 5000
dt  = 0.01  ; 1 fs
; Output control
nstxout = 0 ; save coordinates every 10 ps
nstvout = 0 ; save velocities every 10 ps
nstcalcenergy   = 50
nstenergy   = 1000  ; save energies every 1 ps
nstxtcout   = 2500
;nstxout-compressed  = 5000   ; save compressed coordinates every 1.0 ps
 ; nstxout-compressed replaces nstxtcout
;compressed-x-grps  = System  ; replaces xtc-grps
nstlog  = 1000  ; update log file every 1 ps
; Bond parameters
continuation= no   ; first dynamics run
constraint_algorithm = lincs ; holonomic constraints
constraints = none  ; all bonds (even heavy atom-H bonds)
constrained
;lincs_iter = 2 ; accuracy of LINCS
lincs_order = 4 ; also related to accuracy
epsilon_r   = 15
; Neighborsearching
cutoff-scheme   = Verlet
ns_type = grid  ; search neighboring grid cels
nstlist = 10; 20 fs
rvdw_switch = 1.0
rlist   = 1.2   ; short-range neighborlist cutoff (in nm)
rcoulomb= 1.2   ; short-range electrostatic cutoff (in nm)
rvdw= 1.2   ; short-range van der Waals cutoff (in nm)
vdwtype = Cut-off   ; Twin range cut-offs rvdw >= rlist
;vdw-modifier= Force-switch
;Electrostatics
coulombtype = PME   ; Particle Mesh Ewald for long-range
electrostatics
pme_order   = 4 ; cubic interpolation
fourierspacing  = 0.12  ; grid spacing for FFT
; Temperature coupling is on
tcoupl  = v-rescale
tc_grps = System
tau_t   = 1.0
ref_t   = 300

energygrps = DNA Protein W ION

freezegrps = DNA
freezedim = Y Y Y

; Pressure coupling is off
pcoupl  = no; no pressure coupling in NVT
; Periodic boundary conditions
pbc = xyz   ; 3-D PBC
; Dispersion correctiion
DispCorr= no; account for cut-off vdW scheme
; Velocity generation
gen_vel = yes   ; assign velocities from Maxwell
distribution
gen_temp= 300   ; temperature for Maxwell distribution
gen_seed= -1; generate a random seed
; COM motion removal
; These options remove motion of the protein/bilayer relative to the
solvent/ions
nstcomm = 50
comm-mode   = Linear
comm-grps   = System
;
refcoord_scaling = com


I would greatly appreciate any help.

Thanks in advance,

Thanks and Regards,

Arnab Mukherjee
-- 
Gromacs Users mailing list

* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting!

* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

* For (un)subscribe requests visit
https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a 
mail to gmx-users-requ...@gromacs.org.


[gmx-users] Freeze a group of atoms in gromacs

2017-02-07 Thread ARNAB MUKHERJEE
Hello,

I am simulating a course grained DNA with counterions. I want to know if
there is a way in gromacs to freeze all the degrees of freedom of DNA, so
as to completely restrict its motion. So, the DNA in my system will be
fixed while the ions will be moving. How can I achieve this in gromacs?

Thanks a lot,

Regards,

Arnab Mukherjee
-- 
Gromacs Users mailing list

* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting!

* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

* For (un)subscribe requests visit
https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a 
mail to gmx-users-requ...@gromacs.org.


[gmx-users] Having WCA bonded interaction with FENE potential

2017-01-24 Thread ARNAB MUKHERJEE
Hello,

I am simulating a charged bead-spring modelled polyelectrolyte with FENE
bond stretching potential and Deby Huckel potential for electrostatics. I
need to combine FENE bonded interaction combined with WCA potential,
otherwise the beads overlap with each other since there is no replusive
term. Before this I ran the same system with no electrostatics, so to take
into account WCA for bonded interaction I did No. of bonded exclusions
nexcl = 0, and the results of the simulation seem to be ok.

Now for the same system with electrostatics, if I do nexcl = 0, then there
will also be electrostatic repulsion between the bonded beads, which I do
not want. So, how can I have WCA interaction between the bonded atoms,
without having electrostatic interaction between them?

I know in worst case, I can build a table for the bond stretching
potential, where I can add the WCA term with FENE potential. But is there
an easier way to do this? As FENE potential is already there as one of the
bonded interactions in gromacs.

Thanks a lot,

Regards,

Arnab Mukherjee
-- 
Gromacs Users mailing list

* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting!

* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

* For (un)subscribe requests visit
https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a 
mail to gmx-users-requ...@gromacs.org.


[gmx-users] Cosine bending potential

2016-12-03 Thread ARNAB MUKHERJEE
Hi,

What is the function no. (funct) for the functional form of the bending
potential U_bend = k_theta*(cos theta - cos theta_0)^2 ?

Thanks and Regards,

Arnab
-- 
Gromacs Users mailing list

* Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting!

* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

* For (un)subscribe requests visit
https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a 
mail to gmx-users-requ...@gromacs.org.