[gmx-users] units for atomic coordinates in traj.xtc? nm?

2011-03-25 Thread Jennifer Williams

A very quick question.

Can someone confirm that atomic coordinates are stored as in the  
traj.xtc file in units of nm (as opposed to Angstroms)?


I am trying to figure out whether or not in the g_msd program reads in  
nm or Angstroms. The standard output of the xy graph is nm2 and ps. I  
can't see that the units are being divided by 10 in the g_msd code so  
I am guessing the traj.xtc feeds in the units in nm.  Just want to be  
sure...


Thanks




--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

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


Re: [gmx-users] Instantaneous Square Displacement

2011-03-09 Thread Jennifer Williams



Thanks for the advice. I has also found that last source on google and  
have been thinking how I could apply this.


I assume that if I plotted [rt - r0]^2 against time then (for a single  
molecule) I would get peaks on the graph when the molecules are mobile  
and dips (where [rt-r0]^2 is close to zero) for the periods where  
molecules are stuck.


This would seem to make sense only for a single molecule (as Mark  
suggested) as averaged over all molecules, the peaks and troughs would  
average out and I wouldn't really be illustrating anything.


I have made a rough sketch of what the plot should look like and for a  
few moves of the molecules (using only 15 positions) it seems to make  
sense. However I don't really find it intuitive given that iSD isn't  
widely used and for more than a few moves it gets crowded and  
difficult to interpret.


So as I see it, the graph would sample only one molecule over a short  
time period.


Any more thoughts on this?

Thanks

Jenny

Quoting "Justin A. Lemkul" :




Mark Abraham wrote:

On 8/03/2011 3:01 AM, Jennifer Williams wrote:


Hi,

I am writing a paper where I describe that gas molecules move  
inside a pore and then stick for long periods of time in  
occlusions in the pore wall.


A reviewer has mentioned that I could illustrate this effect by  
using "instantaneous square-displacement".


I have already produced MSD vs time plots and used them to obtain  
the self diffusion coefficient. Can someone shed some light on how  
I can obtain the instantaneous square displacement in gromacs?


I have no idea what "ISD" means, and Google doesn't know either :)  
Perhaps they want to see the diffusion of a single molecule?




Searching for "instantaneous square displacement" turns up very  
little (3 results), but the last seems to be what you need, as long  
as this person is correct:


http://smartech.gatech.edu/bitstream/handle/1853/13994/bai_xianming_200612_phd.pdf?sequence=1

Section 2.3.3.

-Justin

--


Justin A. Lemkul
Ph.D. Candidate
ICTAS Doctoral Scholar
MILES-IGERT Trainee
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at  
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the www  
interface or send it to gmx-users-requ...@gromacs.org.

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






Dr. Jennifer Williams
Institute for Materials and Processes
School of Engineering
University of Edinburgh
Sanderson Building
The King's Buildings
Mayfield Road
Edinburgh, EH9 3JL, United Kingdom
Phone: ++44 (0)131 650 4 861


--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

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


Re: [gmx-users] using one thermostat for entire system ??

2011-03-09 Thread Jennifer Williams



Thanks Mark for the advice.

I have just rerun a test simulation with each of my gas species  
coupled separately to a thermostat and have got similar values for the  
quantities of interest (diffusion coefficients). However I am not sure  
that this will satisfy the reviewer without a bit more justification  
of why I chose to use a single thermostat for my system.


I remembered reading some gromacs information on the application of  
thermostats (on which I based my decision to apply one thermostat for  
the entire system). I have just managed to locate it on the FAQ page:


http://www.gromacs.org/Documentation/Terminology/Thermostats

What Not To Do
Some hints on practices that generally not a good idea to use:
?	Do not use separate thermostats for every component of your system.  
Some molecular dynamics thermostats only work well in the  
thermodynamic limit.  A group must be of sufficient size to justify  
its own thermostat.  If you use one thermostat for, say, a small  
molecule, another for protein, and another for water, you are likely  
introducing errors and artifacts that are hard to predict. In  
particular, do not couple ions in aqueous solvent in a separate group  
from that solvent.


My system consists of
55 molecules of CO2
38 molecules of N2
3 molecules of O2
One large molecule of MCM-41 (mostly frozen but with mobile surface  
groups) consisting of 4284 atoms.


My take on this was that the gas molecules were in small supply  
compared to the other part of the system so I should avoid using  
separate thermostats. Was I mistaken in making this assumption?


Is there any feeling for what "sufficient size" as referred to above  
is? Does anyone know of any papers I could reference in my explanation?


Any comments welcome

Thanks

Jenny






Quoting Mark Abraham :


On 8/03/2011 2:50 AM, Jennifer Williams wrote:

Hi,

I simulated the diffusion of small gases (CO2, N2) in a framework  
structure which was mostly frozen with some mobile surface groups.  
I applied a temperature thermostat to the entire system (i.e I  
didn't couple the gas molecules and framework separately). I have  
now been asked the following by a reviewer:


"Please comment on any artefacts that might arise as a result of  
non-equipartition of energies.  For example, what is the calculated  
temperature of each of the gas species and mobile species?"


I have tried to explore the temperature of each species separately  
but g_energy will only give me the energy of the system as a whole.  
Are there any other tools within gromacs to look at the temperature  
of each species in turn from the md runs I already have?


Does anyone have a suggestion as to how I can address the  
reviewer's comment and proove that using one thermostat for the  
whole system is OK (that is what I am hoping!).  Is there anything  
in particular that I should be looking closely at?


g_energy can only report the global temperature as saved in the .edr  
file, computed during the simulation. The only way of decomposing  
the temperature is to save velocities to the .trr file with nstvout.  
(Some colleagues had to do this recently.) Then, judicious use of  
trjconv to take subsets and matching topology hacking (tpbconv  
doesn't quite work) mean you can use mdrun -rerun on the subset .trr  
and matching .tpr to get a group-wise temperature. That can  
demonstrate whether temperature differentials (ie. heat flows) exist.


If you didn't save velocities, you're out of luck.

Mark
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at  
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the www  
interface or send it to gmx-users-requ...@gromacs.org.

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






Dr. Jennifer Williams
Institute for Materials and Processes
School of Engineering
University of Edinburgh
Sanderson Building
The King's Buildings
Mayfield Road
Edinburgh, EH9 3JL, United Kingdom
Phone: ++44 (0)131 650 4 861


--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

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


[gmx-users] using one thermostat for entire system ??

2011-03-07 Thread Jennifer Williams

Hi,

I simulated the diffusion of small gases (CO2, N2) in a framework  
structure which was mostly frozen with some mobile surface groups. I  
applied a temperature thermostat to the entire system (i.e I didn't  
couple the gas molecules and framework separately). I have now been  
asked the following by a reviewer:


"Please comment on any artefacts that might arise as a result of  
non-equipartition of energies.  For example, what is the calculated  
temperature of each of the gas species and mobile species?"


I have tried to explore the temperature of each species separately but  
g_energy will only give me the energy of the system as a whole. Are  
there any other tools within gromacs to look at the temperature of  
each species in turn from the md runs I already have?


Does anyone have a suggestion as to how I can address the reviewer's  
comment and proove that using one thermostat for the whole system is  
OK (that is what I am hoping!).  Is there anything in particular that  
I should be looking closely at?


Thanks

Jenny


--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

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


[gmx-users] Instantaneous Square Displacement

2011-03-07 Thread Jennifer Williams


Hi,

I am writing a paper where I describe that gas molecules move inside a  
pore and then stick for long periods of time in occlusions in the pore  
wall.


A reviewer has mentioned that I could illustrate this effect by using  
"instantaneous square-displacement".


I have already produced MSD vs time plots and used them to obtain the  
self diffusion coefficient. Can someone shed some light on how I can  
obtain the instantaneous square displacement in gromacs?


Thanks

Jenny

--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

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


[gmx-users] compiling error in tools

2011-01-25 Thread Jennifer Williams


Hello

I know this is possibly an issue for the IT support at my Uni but I  
was wondering if someone could shed some light on what may have gone  
wrong so at least I can point them in the right direction.


As I am getting very small diffusion coefficients when using g_msd, I  
want to increase the number of decimal points the diffusion  
coefficient is displayed to.

In the file gmx_msd.c within src/tools I changed

 fprintf(out,"# D[%10s] = %.4f (+/- %.4f) (1e-5 cm^2/s)\n",

to

 fprintf(out,"# D[%10s] = %.8f (+/- %.8f) (1e-5 cm^2/s)\n",

then, as I have always done, I just typed in make g_msd to recompile  
this tool.

I get the following error message:

src/tools> make g_msd
mpicc -DHAVE_CONFIG_H -I. -I../../src  -I../../include  
-DGMXLIBDIR=\"/home/jwillia4/GRO/share/top\"  
-I/home/jwillia4/GRO/include  -O3 -fomit-frame-pointer  
-finline-functions -Wall -Wno-unused -funroll-all-loops -MT g_msd.o  
-MD -MP -MF .deps/g_msd.Tpo -c -o g_msd.o g_msd.c

mv -f .deps/g_msd.Tpo .deps/g_msd.Po
/bin/sh ../../libtool --tag=CC   --mode=link mpicc  -O3  
-fomit-frame-pointer -finline-functions -Wall -Wno-unused  
-funroll-all-loops  -L/home/jwillia4/GRO/lib  -lgslcblas  -o g_msd  
g_msd.o libgmxana_mpi.la ../mdlib/libmd_mpi.la ../gmxlib/libgmx_mpi.la  
-lgsl -lnsl -lfftw3f -lm
mpicc -O3 -fomit-frame-pointer -finline-functions -Wall -Wno-unused  
-funroll-all-loops -o g_msd g_msd.o  -L/home/jwillia4/GRO/lib  
./.libs/libgmxana_mpi.a  
/home/jwillia4/GRO/gromacs-4.0.7/src/mdlib/.libs/libmd_mpi.a  
../mdlib/.libs/libmd_mpi.a  
/home/jwillia4/GRO/gromacs-4.0.7/src/gmxlib/.libs/libgmx_mpi.a  
../gmxlib/.libs/libgmx_mpi.a /home/jwillia4/GRO/lib/libgslcblas.so  
/home/jwillia4/GRO/lib/libgsl.so -lnsl  
/home/jwillia4/GRO/lib/libfftw3f.a -lm   -Wl,--rpath  
-Wl,/home/jwillia4/GRO/lib -Wl,--rpath -Wl,/home/jwillia4/GRO/lib
/home/jwillia4/GRO/gromacs-4.0.7/src/gmxlib/.libs/libgmx_mpi.a(main.o): In  
function `init_par':

main.c:(.text+0xc0): undefined reference to `ompi_mpi_comm_world'
main.c:(.text+0xc8): undefined reference to `ompi_mpi_comm_world'
/home/jwillia4/GRO/gromacs-4.0.7/src/gmxlib/.libs/libgmx_mpi.a(main.o): In  
function `init_multisystem':

main.c:(.text+0xb99): undefined reference to `ompi_mpi_comm_world'
main.c:(.text+0xbd9): undefined reference to `ompi_mpi_comm_world'
main.c:(.text+0xbee): undefined reference to `ompi_mpi_comm_world'
/home/jwillia4/GRO/gromacs-4.0.7/src/gmxlib/.libs/libgmx_mpi.a(network.o):network.c:(.text+0x15): more undefined references to `ompi_mpi_comm_world'  
follow
/home/jwillia4/GRO/gromacs-4.0.7/src/gmxlib/.libs/libgmx_mpi.a(network.o): In  
function `gmx_sumi_sim':

network.c:(.text+0x96): undefined reference to `ompi_mpi_op_sum'
network.c:(.text+0x9b): undefined reference to `ompi_mpi_int'
/home/jwillia4/GRO/gromacs-4.0.7/src/gmxlib/.libs/libgmx_mpi.a(network.o): In  
function `gmx_sumf_sim':

network.c:(.text+0x4e6): undefined reference to `ompi_mpi_op_sum'
network.c:(.text+0x4eb): undefined reference to `ompi_mpi_float'
/home/jwillia4/GRO/gromacs-4.0.7/src/gmxlib/.libs/libgmx_mpi.a(network.o): In  
function `gmx_sumd_sim':

network.c:(.text+0x9a6): undefined reference to `ompi_mpi_op_sum'
network.c:(.text+0x9ab): undefined reference to `ompi_mpi_double'
/home/jwillia4/GRO/gromacs-4.0.7/src/gmxlib/.libs/libgmx_mpi.a(network.o): In  
function `gmx_sumi':

network.c:(.text+0xe8f): undefined reference to `ompi_mpi_op_sum'
network.c:(.text+0xe94): undefined reference to `ompi_mpi_int'
network.c:(.text+0xec0): undefined reference to `ompi_mpi_int'
network.c:(.text+0xed7): undefined reference to `ompi_mpi_op_sum'
network.c:(.text+0xedc): undefined reference to `ompi_mpi_int'
network.c:(.text+0x12fe): undefined reference to `ompi_mpi_op_sum'
network.c:(.text+0x1303): undefined reference to `ompi_mpi_int'
/home/jwillia4/GRO/gromacs-4.0.7/src/gmxlib/.libs/libgmx_mpi.a(network.o): In  
function `gmx_sumf':

network.c:(.text+0x134e): undefined reference to `ompi_mpi_float'
network.c:(.text+0x1354): undefined reference to `ompi_mpi_op_sum'
network.c:(.text+0x1380): undefined reference to `ompi_mpi_float'
network.c:(.text+0x1397): undefined reference to `ompi_mpi_op_sum'
network.c:(.text+0x139c): undefined reference to `ompi_mpi_float'
network.c:(.text+0x183e): undefined reference to `ompi_mpi_op_sum'
network.c:(.text+0x1843): undefined reference to `ompi_mpi_float'
/home/jwillia4/GRO/gromacs-4.0.7/src/gmxlib/.libs/libgmx_mpi.a(network.o): In  
function `gmx_sumd':

network.c:(.text+0x1890): undefined reference to `ompi_mpi_op_sum'
network.c:(.text+0x1895): undefined reference to `ompi_mpi_double'
network.c:(.text+0x18c2): undefined reference to `ompi_mpi_double'
network.c:(.text+0x18d7): undefined reference to `ompi_mpi_op_sum'
network.c:(.text+0x18dc): undefined reference to `ompi_mpi_double'
network.c:(.text+0x1d86): undefined reference to `ompi_mpi_op_sum'
network.c:(.text+0x1d8b): undefined reference to `ompi_mpi_double'
/home/jwillia4/GRO/gromacs-4.0.7/src/gmxlib/.libs/libgmx_m

[gmx-users] energy group exlusions for frozen structures-use of maxwarn -1?

2010-10-28 Thread Jennifer Williams

Hello.

I am using a porous structure (MOF). I am interested in how guest  
molecules inside the porous structure interact with each other and the  
pores of the structure. As my structure represents a rigid crystalline  
framework it is usually kept rigid/frozen and only non bonded  
parameters (LJ and electrostatics) are taken into account for this  
structure.


(I think) I read somewhere on the forum the energy group exclusions  
should be applied for frozen atoms (also if I don?t do this I get some  
v. large energies).


So I use the following in my .mdp file

; Selection of energy groups
energygrps   = MOF

; ENERGY GROUP EXCLUSIONS
; Pairs of energy groups for which all non-bonded interactions are excluded
energygrp_excl   = MOF MOF

However on running grompp I get the following warning:

WARNING 1 [file MCM_TDI.top, line 32874]:
  Can not exclude the lattice Coulomb energy between energy groups

Why is this? My frozen structure has partial charges but is overall  
charge neutral.
To get around this I have to use maxwarn -1. I am a bit wary of using  
maxwarn. Can someone advise whether it is sensible to use maxwarn in  
this case?


Thanks

Jenny



--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

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


Re: [gmx-users] genconf to insert small molecules triclinic unit cell

2010-10-22 Thread Jennifer Williams


Hi Justin,

The .pdb I sent with my last e-mail was a result of trying to add 500  
molecules. The max of 168 were added and then genconf couldn't add  
anymore despite the large vacant space. Trying to add more molecules  
to this .pdb won't work. This is what makes me think the space was  
somehow inaccessible.


If I carry out an mdrun on this structure, the NO molecules distribute  
themselves into this empty space. If I take the structure file  
following MD I can add an extra 24 molecules with genbox.


The worry is not for these small NO molecules where I can use MD to  
create more space. I am also trying to insert larger organic molecules  
and if parts of the unit cell are seen as inaccessible I won't be able  
to get the molecules in there in the first place.



Jenny




Quoting "Justin A. Lemkul" :




Jennifer Williams wrote:


Hi Justin,

Thanks for the response (again!).

I would be happy if genbox inserted molecules randomly (within the   
box defined by my .pdf file) but the output doesn't look random but  
 subject to some strange (symmetry?)constraints which don't allow   
random insertions into some areas of my unit cell.


The molecules will crowd into one space leaving another portion of   
the cell completely empty. The inserted molecules assume the shape   
of a box whilst the rest of my structure is a triclinic cell.


One pic is probably worth a thousand words. I've attached a .pdb of  
 the final structure. I'd appreciate if you could tell me if this  
is  normal output for genbox,




I see nothing wrong with it.  Looks like genbox started adding
molecules on one "side" of your box, and continued adding until it
reached 100, then stopped. Looks like you simply have more than enough
space to add more of your NO molecules, that's all.  The molecules are
all within the unit cell, so I see no issue.

If you want a more homogeneous distribution of NO, you may have to come
up with a different method.

-Justin


Thanks very much

Jenny




Quoting "Justin A. Lemkul" :




Jennifer Williams wrote:

Hi,

I have a strange problem when I try to insert small molecules   
into  my cell using genconf. Here the cell is a crystalline   
metal-organic  framework with lots of space for the molecules   
inside. The edges of  the metal organic framework defined the pbcs.


When I view the final structure, the inserted molecules appear to  
  have a different symmetry to that of the MOF/cell. They   
notieably  avoid some empty regions inside the box and are   
inserted outside  regions of my triclinic cell. The inserted   
molecules don?t occupy a  triclinic shape at all.


I have used this feature before on a similar triclinic cell and   
it  worked perfectly. I can?t tell what I am doing wrong this time.


I outline the steps I take below.

1. I am trying to insert a number of molecules of NO into a
triclinic crystal cell. The pdb file of the crystal cell has the   
 following cell parameters.


CRYST1   25.885   25.8856.806  90.00  90.00 120.00 P 1   1

This is the pdb file of the small molecule I am trying to insert:

CRYST1   25.885   25.885   27.2232  90.00  90.00 120.00 P 1   1
HETATM1  NNO NOO 1   1.150   0.000   0.000  1.00 0.00  
N
HETATM2  ONO NOO 1   0.000   0.000   0.000  1.00 0.00  
O

CONECT12
END

I use the following commands:

To centre the cell:
editconf -c -f  MOF.pdb -o MOF_centered.pdb

To make a 1x1x4 box ( I need to increase the length in z) so that  
 I  can use sensible cut-offs later on


genconf -nbox 1 1 4 -f MOF_centered.pdb -o MOF_4c.pdb


To add the NO molecules:

genbox ?cp MOF_4c.pdb -ci NOO.pdb -nmol 100 -o inserted_NO.pdb

Any ideas appreciated as I have been going around in circles.




When using genbox -ci, there is no guarantee of any type of symmetry.
It inserts molecules randomly, wherever it finds space.  I don't think
you're doing anything wrong, I think you're just encountering a
limitation of genbox.

-Justin


I am using gromacs 4.0.7

Thanks

Jenny




--


Justin A. Lemkul
Ph.D. Candidate
ICTAS Doctoral Scholar
MILES-IGERT Trainee
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the www
interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists








--


Justin A. Lemkul
Ph.D. Candidate
ICTAS Doctoral Scholar
MILES-IGERT Trainee
D

[gmx-users] genconf to insert small molecules triclinic unit cell

2010-10-22 Thread Jennifer Williams

Hi,

I have a strange problem when I try to insert small molecules into my  
cell using genconf. Here the cell is a crystalline metal-organic  
framework with lots of space for the molecules inside. The edges of  
the metal organic framework defined the pbcs.


When I view the final structure, the inserted molecules appear to have  
a different symmetry to that of the MOF/cell. They notieably avoid  
some empty regions inside the box and are inserted outside regions of  
my triclinic cell. The inserted molecules don?t occupy a triclinic  
shape at all.


I have used this feature before on a similar triclinic cell and it  
worked perfectly. I can?t tell what I am doing wrong this time.


I outline the steps I take below.

1. I am trying to insert a number of molecules of NO into a triclinic  
crystal cell. The pdb file of the crystal cell has the following cell  
parameters.


CRYST1   25.885   25.8856.806  90.00  90.00 120.00 P 1   1

This is the pdb file of the small molecule I am trying to insert:

CRYST1   25.885   25.885   27.2232  90.00  90.00 120.00 P 1   1
HETATM1  NNO NOO 1   1.150   0.000   0.000  1.00 0.00   
   N
HETATM2  ONO NOO 1   0.000   0.000   0.000  1.00 0.00   
   O

CONECT12
END

I use the following commands:

  To centre the cell:
editconf -c -f  MOF.pdb -o MOF_centered.pdb

To make a 1x1x4 box ( I need to increase the length in z) so that I  
can use sensible cut-offs later on


genconf -nbox 1 1 4 -f MOF_centered.pdb -o MOF_4c.pdb


To add the NO molecules:

genbox ?cp MOF_4c.pdb -ci NOO.pdb -nmol 100 -o inserted_NO.pdb

Any ideas appreciated as I have been going around in circles.


I am using gromacs 4.0.7

Thanks

Jenny


--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] Polar Hydrogen missing using PRODRG

2010-10-14 Thread Jennifer Williams


Hi,

So, I did as suggested. Ran PRODRG the first time using the JME Editor  
and took the pdb file from this. Then I added the ADDHYD OAE  
underneath (actually it doesn't matter which O (OAD or OAE) is  
protonated as long as it is only one of them).


The first bit of output looked encouraging...

PRODRG> PDB mode detected.
PRODRG> Molecule complexity index: 2.00.
PRODRG> Added hydrogen on atom OAE .
PRODRG>  18 hydrogen(s) added.
PRODRG> Using charge groups.
PRODRG> Net charge on molecule:   0.000
PRODRG>  20 partial charges  0 ambiguous
PRODRG> Cannot assign type to atom ' OAE'.
ERRDRG> Error in GROMOS atom names/types.
PRODRG> Program terminated unsuccessfully, sorry!

So it understood that I wanted a hydrogen and the net charge is at  
last 0(usually -1 without that hydrogen) but then it crashes! I don't  
understand why it can't assign a type to OAE. IS it because a hydrogen  
is now attached so the hybridisation has changed from sp2 to sp3?


Do I need to change the identity of OAE before putting a hydrogen on?  
If so to what?


Thanks

Jenny


REMARK
REMARK
REMARK  This file was generated by PRODRG version 071121.0636
REMARK  PRODRG written/copyrighted by Daan van Aalten
REMARK  and Alexander Schuettelkopf
REMARK
REMARK  Questions/comments to d...@davapc1.bioch.dundee.ac.uk
REMARK
REMARK  When using this software in a publication, cite:
REMARK  A. W. Schuettelkopf and D. M. F. van Aalten (2004).
REMARK  PRODRG - a tool for high-throughput crystallography
REMARK  of protein-ligand complexes.
REMARK  Acta Crystallogr. D60, 1355--1363.
REMARK
REMARK
HETATM1  CAA DRG 1   2.360   1.610   1.840  1.00 20.00  
C
HETATM2  CAN DRG 1   3.580   2.340   1.260  1.00 20.00  
C
HETATM3  CAB DRG 1   3.330   2.650  -0.210  1.00 20.00  
C
HETATM4  CAJ DRG 1   4.820   1.450   1.460  1.00 20.00  
C
HETATM5  CAL DRG 1   6.160   2.130   1.140  1.00 20.00  
C
HETATM6  CAG DRG 1   6.870   1.740  -0.010  1.00 20.00  
C
HETATM7  CAI DRG 1   8.100   2.360  -0.300  1.00 20.00  
C
HETATM8  CAF DRG 1   6.640   3.140   2.010  1.00 20.00  
C
HETATM9  CAH DRG 1   7.870   3.760   1.720  1.00 20.00  
C
HETATM   10  CAM DRG 1   8.600   3.370   0.570  1.00 20.00  
C
HETATM   11  CAO DRG 1   9.930   4.110   0.310  1.00 20.00  
C
HETATM   12  CAC DRG 1  10.120   4.580  -1.140  1.00 20.00  
C
HETATM   13  CAK DRG 1  11.150   3.250   0.700  1.00 20.00  
C
HETATM   14  OAE DRG 1  12.100   3.870   1.240  1.00 20.00  
O
HETATM   15  OAD DRG 1  11.130   2.020   0.480  1.00 20.00  
O

CONECT12
CONECT2134
CONECT32
CONECT425
CONECT5468
CONECT657
CONECT76   10
CONECT859
CONECT98   10
CONECT   1079   11
CONECT   11   10   12   13
CONECT   12   11
CONECT   13   11   14   15
CONECT   14   13
CONECT   15   13
END
ADDHYD OAE








Quoting "Justin A. Lemkul" :



You need to run PRODRG twice.  The first time, do not use ADDHYD.  Make
note of the atom names that PRODRG assigns.  Then run PRODRG a second
time, using ADDHYD with the atom name that PRODRG uses for the O atom
you want protonated.

-Justin

Jennifer Williams wrote:



Hi Justin,

Thanks for your answer. The ADDHYD in the FAQ sounds like exactly   
what I need however I can't get this command to work. I draw in my   
molecule using the JME editor and get this .pdb file in the PRODRG   
window.


CC(C)Cc1ccc(C(C)C(=O)O)cc1
JME 2002.05 Thu Oct 14 10:48:38 BST 2010

15 15V2000
   0.0.0. C
   1.21242.10000. C
   7.27464.20000. C
   8.48700.70000. O
   9.69952.80000. O
   3.63732.10000. C
   4.84970.0. C
   4.84972.80000. C
   6.06220.70000. C
   2.42490.0. C
   8.48702.10000. C
   3.63730.70000. C
   6.06222.10000. C
   1.21240.70000. C
   7.27462.80000. C
 1 14  1  0
 2 14  1  0
 3 15  1  0
 4 11  2  0
 5 11  1  0
 6  8  1  0
 6 12  2  0
 7  9  2  0
 7 12  1  0
 8 13  2  0
 9 13  1  0
10 12  1  0
10 14  1  0
11 15  1  0
13 15  1  0
M  END

If I then try adding ADDHYD O to the window and clicking "run   
PRODRG" I get the following:


ERRDRG> Atom 'O' referenced in instruction was not found.
PRODRG> Program terminated unsuccessfully, sorry!

I've then tried renaming O to OM as there are 2 Os in my structure   
and I only want one protonated. However, no other name is   
recognized so I get this message:


ERRDRG> Currently o

Re: [gmx-users] Polar Hydrogen missing using PRODRG

2010-10-14 Thread Jennifer Williams



Hi Justin,

Thanks for your answer. The ADDHYD in the FAQ sounds like exactly what  
I need however I can't get this command to work. I draw in my molecule  
using the JME editor and get this .pdb file in the PRODRG window.


CC(C)Cc1ccc(C(C)C(=O)O)cc1
JME 2002.05 Thu Oct 14 10:48:38 BST 2010

 15 15V2000
0.0.0. C
1.21242.10000. C
7.27464.20000. C
8.48700.70000. O
9.69952.80000. O
3.63732.10000. C
4.84970.0. C
4.84972.80000. C
6.06220.70000. C
2.42490.0. C
8.48702.10000. C
3.63730.70000. C
6.06222.10000. C
1.21240.70000. C
7.27462.80000. C
  1 14  1  0
  2 14  1  0
  3 15  1  0
  4 11  2  0
  5 11  1  0
  6  8  1  0
  6 12  2  0
  7  9  2  0
  7 12  1  0
  8 13  2  0
  9 13  1  0
 10 12  1  0
 10 14  1  0
 11 15  1  0
 13 15  1  0
M  END

If I then try adding ADDHYD O to the window and clicking "run PRODRG"  
I get the following:


ERRDRG> Atom 'O' referenced in instruction was not found.
PRODRG> Program terminated unsuccessfully, sorry!

I've then tried renaming O to OM as there are 2 Os in my structure and  
I only want one protonated. However, no other name is recognized so I  
get this message:


ERRDRG> Currently only N, C, O, S, P, Cl, Br, F, I are supported.
PRODRG> Program terminated unsuccessfully, sorry!

Next I tried using the PDB generated by PRODRG and repasting it with  
the ADDHYD command into the PRODRG window.


Now the pdb looks like this:

REMARK
REMARK
REMARK  This file was generated by PRODRG version 071121.0636
REMARK  PRODRG written/copyrighted by Daan van Aalten
REMARK  and Alexander Schuettelkopf
REMARK
REMARK  Questions/comments to d...@davapc1.bioch.dundee.ac.uk
REMARK
REMARK  When using this software in a publication, cite:
REMARK  A. W. Schuettelkopf and D. M. F. van Aalten (2004).
REMARK  PRODRG - a tool for high-throughput crystallography
REMARK  of protein-ligand complexes.
REMARK  Acta Crystallogr. D60, 1355--1363.
REMARK
REMARK
HETATM1  CAA DRG 1   2.360   1.610   1.840  1.00 20.00  
C
HETATM2  CAN DRG 1   3.580   2.340   1.260  1.00 20.00  
C
HETATM3  CAB DRG 1   3.330   2.650  -0.210  1.00 20.00  
C
HETATM4  CAJ DRG 1   4.820   1.450   1.460  1.00 20.00  
C
HETATM5  CAL DRG 1   6.160   2.130   1.140  1.00 20.00  
C
HETATM6  CAG DRG 1   6.870   1.740  -0.010  1.00 20.00  
C
HETATM7  CAI DRG 1   8.100   2.360  -0.300  1.00 20.00  
C
HETATM8  CAF DRG 1   6.640   3.140   2.010  1.00 20.00  
C
HETATM9  CAH DRG 1   7.870   3.760   1.720  1.00 20.00  
C
HETATM   10  CAM DRG 1   8.600   3.370   0.570  1.00 20.00  
C
HETATM   11  CAO DRG 1   9.930   4.110   0.310  1.00 20.00  
C
HETATM   12  CAC DRG 1  10.120   4.580  -1.140  1.00 20.00  
C
HETATM   13  CAK DRG 1  11.150   3.250   0.700  1.00 20.00  
C
HETATM   14  OAE DRG 1  12.100   3.870   1.240  1.00 20.00  
O
HETATM   15  OAD DRG 1  11.130   2.020   0.480  1.00 20.00  
O

CONECT12
CONECT2134
CONECT32
CONECT425
CONECT5468
CONECT657
CONECT76   10
CONECT859
CONECT98   10
CONECT   1079   11
CONECT   11   10   12   13
CONECT   12   11
CONECT   13   11   14   15
CONECT   14   13
CONECT   15   13
END

I've tried changing the HETATM to OM and using ADDHYD OM and then  
played around with various other things like renaming the atom symbol  
at the end of the line and even the OAD, DRG but I always get the  
following:


ERRDRG> Atom 'OM' referenced in instruction was not found.
PRODRG> Program terminated unsuccessfully, sorry!

I am using the PRODRG beta server available online at

http://davapc1.bioch.dundee.ac.uk/cgi-bin/prodrg_beta

Will I fare any better if I use the source code or is there something  
obvious I am doing wrong? Does the ADDHYD have to be added in any  
particular location of the pdb file? I have tried a few random  
positions including at the very top, very bottom and just beneath the  
O to be protonated. Nothing worked.


Any ideas much appreciated,

Jenny



Quoting "Justin A. Lemkul" :




Jennifer Williams wrote:



Hi,

I am trying to get a .itp file for a simple molecule (Ibuprofen).   
This contains a COOH group.


THe problem is that PRODRG removes the polar hydrogen of the COOH   
from the .pdb file and generates a .itp file without it.


I am not concerned about the other aromatic hydrogens but I really   
need

[gmx-users] Polar Hydrogen missing using PRODRG

2010-10-13 Thread Jennifer Williams



Hi,

I am trying to get a .itp file for a simple molecule (Ibuprofen). This  
contains a COOH group.


THe problem is that PRODRG removes the polar hydrogen of the COOH from  
the .pdb file and generates a .itp file without it.


I am not concerned about the other aromatic hydrogens but I really  
need the polar hydrogen modelled explicitly.


Is there a fix for this within PRODRG?

I know I can use pdb2gmx to add hydrogens but the residue names of the  
pdb file are not recognized by the existing forcefields so I usually  
bypass using pdb2gmx.



Thanks

Jenny



Dr. Jennifer Williams
Institute for Materials and Processes
School of Engineering
University of Edinburgh
Sanderson Building
The King's Buildings
Mayfield Road
Edinburgh, EH9 3JL, United Kingdom
Phone: ++44 (0)131 650 4 861


--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

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


[gmx-users] Convert LJ 9-6 to LJ 12-6?

2010-09-15 Thread Jennifer Williams


Hi

I want to check some of my gromacs results using another computer  
program which is set up to read in LJ 12-6 parameters only.


I have a set of 9-6 LJ parameters. Can gromacs interconvert LJ-9-6 to  
12-6? If not can anyone point me in the right direction to a paper,  
bit of code etc.


Changing the source code of the other program to read in LJ 9-6 is  
something I'd rather avoid if there is another solution.



Thanks





--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

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


[gmx-users] pbc in one direction only for analysis?

2010-08-31 Thread Jennifer Williams


Dear gromacs users

I have surface groups anchored on a cylindrical pore wall (similar to  
a carbon nanotube). The pore runs along the z direction. I am trying  
to determine to what extent my surface groups are clustered together  
and was thinking of using g_mindist and/or g_rdf for this analysis.


I usually work with periodic boundary conditions in all 3 directions  
for MD runs.


The problem I have is that I want to use periodic boundaries only in  
the z direction and not the x and y for the anaysis. I.e I dont want  
the surface groups in my pore to see other surface groups in  
neighbouring cells in the x or y direction.


From what I can see, the options with g_rdf is either to have pbc in  
all 3 directions or not at all. Does anyone know a way around this?


I am analyzing one pdb file and not a trajectory.

I am using gromacs 4.0.7

Thanks



--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

Can't post? Read http://www.gromacs.org/mailing_lists/users.php


[gmx-users] Can I make make_ndx interactive?

2010-08-03 Thread Jennifer Williams


Hi,


Is there some way of making the make_ndx file interactive so that I  
can include it in a script?


I would usually type in  "a C_O and O_C" to select my atoms

I've tried altering the example on the webpage for making commands  
interactive:


make_ndx -flags 

[gmx-users] Re: Can I make make_ndx interactive?

2010-08-03 Thread Jennifer Williams


 I realised "flags" is not part of the command.

!! Sorry for the silly question. It's the end of a long day and I'm tired!



Quoting Jennifer Williams :



Hi,


Is there some way of making the make_ndx file interactive so that I can
include it in a script?

I would usually type in  "a C_O and O_C" to select my atoms

I've tried altering the example on the webpage for making commands
interactive:

make_ndx -flags <


--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

Can't post? Read http://www.gromacs.org/mailing_lists/users.php


[gmx-users] no non bonded parameters from PRODRG

2010-07-30 Thread Jennifer Williams



Hi,

I am trying to create .top files for some complex organic molecules  
using the PRODRG server and the gromos 96 forcefield. PRODRG gives me  
bonded parameters (stretching, angle bending, torsions) etc but no  
lennard jones parameters (sigma, epsilon) or charges.


Where can I get the non bonded terms from? Does it make sense to  
extract these values from the OPLS forcefield or elsewhere when I am  
using gromas 96 forcefield to get the bonded parameters. I am wary of  
cobbelling forcefield parameters together from different sources,


Thanks



--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

Can't post? Read http://www.gromacs.org/mailing_lists/users.php


[gmx-users] High quality movies, bmp frames and videomach=poor quality

2010-04-20 Thread Jennifer Williams

Hi,

I am using tachyon to get high quality still pictures by generating a  
.dat file, changing the render numbers and then typing this command in  
the window


"N:/myhome/VMD1.8.7/tachyon_WIN32.exe" -aasamples 12 plot01.dat  
-format -o plot01.tga


This gives great images but now I am stuck on generating movies as I  
don?t seem to have much choice when using the movie-maker and am  
clueless at scripting. I need a movie file in .mpg or .avi for  
submission to a journal (and not bigger than 10MB).


Whenever I choose tachyon as the renderer within the movie maker I get  
frames generated in .dat which I can?t do anything with and when I use  
tachyon internal the files are automatically written in .bmp. I can?t  
seem to get my frames in targa or any format other than .bmp.


Is there a script command which could give me high quality snapshots in .tga?

If I continue with my frames in .bmp and generate a movie in videomach  
the image quality is very poor (I do this via preview movie and select  
the highest quality (fast 400x400). TO avoid this I have been using  
convert in unix to stick frames together but this won't give me a mpg  
or avi file!


Can anyone give me some advice or a script line similar to the one  
above which would work on a whole trajectory rather than a single  
frame?.


Thanks

Jenny




--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
Please don't post (un)subscribe requests to the list. Use the
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/mailing_lists/users.php


[gmx-users] reference for gromos 87 within PRODRG

2010-04-19 Thread Jennifer Williams



Hi,

I have generated some topology files using the gromos87 forcefield  
within the PRODRG server. Does anyone have a reference for this  
force-field? I can't seem to find it,


Thanks

Jenny



--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

Can't post? Read http://www.gromacs.org/mailing_lists/users.php


Re: [gmx-users] compilation problems orte error

2010-02-10 Thread Jennifer Williams


Sorry for the delay in replying back. I start the job using the  
following script file:


#$ -S /bin/bash
#$ -l h_rt=47:59:00
#$ -j y
#$ -pe mpich2 8
#$ -cwd
cd /home/jwillia4/GRO/gromacs-4.0.7/JJW_003/PH_TORUN
/home/jwillia4/GRO/bin/mpirun -np 8 /home/jwillia4/GRO/bin/mdrun_mpi  
-v -s md.tpr


The strange thing is that sometimes it works and the job runs to  
completion and sometimes it crashes immediately with the orte error so  
I know that it is not the input files causing the problems. It seems  
entirely random.


Has it to do with the -pe mpich2 8 line? I was previously using Open  
MPI installed on the cluster for common use but now have downloaded  
everything into my home directory. The script has been adapted from  
the time when I didn't have my own OpenMPI in my home directory.  
Perhaps it needs further alteration but I don't know what.


How would I do about checking whether MPI is running?

If you spot anything suspicious in the above commands please let me know.

Thanks

Jenny


Quoting Chandan Choudhury :


As Justin said give the command line options for mdrun and also check that
your mpi environment is running.  Better to run a parallel job and check its
output.

Chadnan

--
Chandan kumar Choudhury
NCL, Pune
INDIA


On Mon, Feb 8, 2010 at 8:02 PM, Justin A. Lemkul  wrote:




Jennifer Williams wrote:



Dear All,

I am having problems compiling gromacs 4.0.7 in parallel. I am following
the
Quick and Dirty Installation instructions on the gromacs webpage.
I downloaded the the versions of fftw, OpenMPI and gromacs-4.0.7 following
these instructions.

Everything seems to compile OK and I get all the serial executables
including mdrun written to my bin directory and they seem to run fine.
However when I try to run mdrun_mpi on 6 nodes I get the following:

[vlxbig16:08666] [NO-NAME] ORTE_ERROR_LOG: Not found in file
runtime/orte_init_stage1.c at line 182
[vlxbig16:08667] [NO-NAME] ORTE_ERROR_LOG: Not found in file
runtime/orte_init_stage1.c at line 182
[vlxbig16:08700] [NO-NAME] ORTE_ERROR_LOG: Not found in file
runtime/orte_init_stage1.c at line 182
[vlxbig16:08670] [NO-NAME] ORTE_ERROR_LOG: Not found in file
runtime/orte_init_stage1.c at line 182
[vlxbig16:08681] [NO-NAME] ORTE_ERROR_LOG: Not found in file
runtime/orte_init_stage1.c at line 182
[vlxbig16:08659] [NO-NAME] ORTE_ERROR_LOG: Not found in file
runtime/orte_init_stage1.c at line 182
--
It looks like orte_init failed for some reason; your parallel process is
likely to abort.  There are many reasons that a parallel process can
fail during orte_init; some of which are due to configuration or
environment problems.  This failure appears to be an internal failure;
here's some additional information (which may only be relevant to an
Open MPI developer):

 orte_rml_base_select failed
 --> Returned value -13 instead of ORTE_SUCCESS


Does anyone have any idea what is causing this? Computer support at my
University is not sure.



How are you launching mdrun_mpi (command line)?

-Justin



Thanks





--


Justin A. Lemkul
Ph.D. Candidate
ICTAS Doctoral Scholar
MILES-IGERT Trainee
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
Please don't post (un)subscribe requests to the list. Use the www interface
or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/mailing_lists/users.php







Dr. Jennifer Williams
Institute for Materials and Processes
School of Engineering
University of Edinburgh
Sanderson Building
The King's Buildings
Mayfield Road
Edinburgh, EH9 3JL, United Kingdom
Phone: ++44 (0)131 650 4 861


--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
Please don't post (un)subscribe requests to the list. Use the
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/mailing_lists/users.php


[gmx-users] compilation problems orte error

2010-02-08 Thread Jennifer Williams


Dear All,

I am having problems compiling gromacs 4.0.7 in parallel. I am following the
Quick and Dirty Installation instructions on the gromacs webpage.
I downloaded the the versions of fftw, OpenMPI and gromacs-4.0.7  
following these instructions.


Everything seems to compile OK and I get all the serial executables  
including mdrun written to my bin directory and they seem to run fine.  
However when I try to run mdrun_mpi on 6 nodes I get the following:


[vlxbig16:08666] [NO-NAME] ORTE_ERROR_LOG: Not found in file  
runtime/orte_init_stage1.c at line 182
[vlxbig16:08667] [NO-NAME] ORTE_ERROR_LOG: Not found in file  
runtime/orte_init_stage1.c at line 182
[vlxbig16:08700] [NO-NAME] ORTE_ERROR_LOG: Not found in file  
runtime/orte_init_stage1.c at line 182
[vlxbig16:08670] [NO-NAME] ORTE_ERROR_LOG: Not found in file  
runtime/orte_init_stage1.c at line 182
[vlxbig16:08681] [NO-NAME] ORTE_ERROR_LOG: Not found in file  
runtime/orte_init_stage1.c at line 182
[vlxbig16:08659] [NO-NAME] ORTE_ERROR_LOG: Not found in file  
runtime/orte_init_stage1.c at line 182

--
It looks like orte_init failed for some reason; your parallel process is
likely to abort.  There are many reasons that a parallel process can
fail during orte_init; some of which are due to configuration or
environment problems.  This failure appears to be an internal failure;
here's some additional information (which may only be relevant to an
Open MPI developer):

  orte_rml_base_select failed
  --> Returned value -13 instead of ORTE_SUCCESS


Does anyone have any idea what is causing this? Computer support at my  
University is not sure.



Thanks



--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
Please don't post (un)subscribe requests to the list. Use the
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/mailing_lists/users.php


[gmx-users] g_sdf and visualizing using gOpenMol

2009-12-15 Thread Jennifer Williams

Hi,

I have a mesoporous silica with attached organic surface groups (alkyl  
chains). I am trying to quantify to what extent these surface groups  
interact with the mesoporous silica surface they are attached to (as  
opposed to projecting straight into the pore space).


I thought about using the g_sdf tool (on one surface group at a time  
initially) but am not completely sure what I am doing is correct so I  
am hoping someone can guide me. My index file looks like this


[MCM]
1
[MCM]
2
[MCM]
3
[SI]
4 5 6 7 8 9 10 11 12 13

Where groups 1, 2 and 3 are carbons of the alkyl chain. SI atoms are  
taken to represent the surface of the mesoporous silica. I run g_sdf  
in gromacs and get a  refmol.gro file (which looks very different to  
my structure) and a gom_plt.dat file. I've installed gOpenMol but I am  
having trouble loading the gom_plt.dat file. It is necessary to load  
coordinates in first so I load in refmol.gro. As soon as I try to load  
the .dat file (using plot> contour>import) the gui window closes and I  
get the following


Will apply a physical translation (x,y,z): 12.236523 18.866343 9.419045
Will apply a MT translation (x,y,z): 0.07 0.04 0.02
Minimum value 0.00 maximum value 134.129913
 Signal was caught 
=> : Success
Signal code is: 11, errno

I am not sure whether this is because I have unreasonable data in my  
.dat file (or because my refmol.gro file looks strange) or because I  
haven't got the hang of using gOpenMol. I don't have any idea what the  
contour plot should like. Could someone e-mail me a working   
gom_plt.dat file so I can at least check that I can correctly load it  
and see what a plot looks like?


Any further advice is appreciated,

Thanks

Jenny




--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
Please don't post (un)subscribe requests to the list. Use the
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/mailing_lists/users.php


Re: [gmx-users] using trjconv to get a 2x2x2 unit cell for traj.xtc

2009-12-08 Thread Jennifer Williams


Hi Justin

Thanks for the suggestions. However, I don't really want to stop the  
molecule moving into the next image-I start off with an alkyl chain  
where a large portion of the structure overlaps the pbc.


i.e C1-C2-C3 \\ C4-C5-C6

where || is the unit cell boundary. To see the whole alkyl chain I use  
the graphics/representations/periodic tab in VMD. The problem is that  
VMD will not allow me to draw a bond between C3 and C4 meaning I get  
odd looking movies. VMD prefers to join the C3 to the C4 in the same  
unit cell rather than the C4 sitting right next to it on the other  
side of the pbc. This means long bonds stretching the length of my  
unit cell and missing bonds between atoms sitting next to each other!  
I have spent a while trying to solve this and have posted to the VMD  
forum (as it is actually a problem with VMD and not gromacs).


If I could somehow do the same thing to the traj.xtc as I did to  
confgro.out using genconf -nbox 2 2 2 that would solve my  
visualisation problem. Do you know of a way to do this? Even a  
round-about way?  I suppose if I dumped each frame in the trajectory  
as a pdb, used genbox on them to get a 2x2x2 box and then somehow  
stuck them together into an .xtc this might work (but any suggestions  
which are less messy are greatly appreciated!)


Jenny






Quoting "Justin A. Lemkul" :




Jennifer Williams wrote:



Hello,

I am trying to find a way around a visualisation problem I am   
having in VMD. Some of my molecules go over periodic boundary   
conditions meaning that bonds sometimes appear missing when looking  
 at movies (I am trying to fix this using wrap, unwrap and join in   
VMD but as yet no luck).


I was wondering if there is a way in gromacs to multiply the number  
 of unit cells shown in a trajectory. i.e instead of a 1x1x1 I want  
 the new unit cell to be 2x2x2. This would mean the section of the   
structure I want to zoom in doesn't go over the pbc.


For the confout.gro file I have done this using

genconf -nbox 2 2 2 -f confout.gro -o confbig.out

and this enables me to at least see a static image where all bonds   
are present.


but in order to view a movie, I need to carry out something similar  
 on the traj.xtc file. I have seen that with trjconv there is the   
option


-box  Size for new cubic box

but my unit cell is not cubic, it is a parallelepiped. The cell   
dimensions are :


4.64210   3.77847   1.89596   0.0   0.0  -2.18150   0.0  
   0.0   0.0


I tried using this anyway with the following command:

trjconv -box 9.28414 8.72598 3.79192 -f traj.xtc -o trajbig.xtc

but the resulting .xtc file wouldn't load into VMD so I assume that  
 the .xtc file and the .gro file didn't match.




Correcting periodicity should alleviate some of the problems.  Using
trjconv -pbc mol (or -pbc nojump) is typically the best first step.

-Justin


Any ideas?,

Thanks

Jenny









--


Justin A. Lemkul
Ph.D. Candidate
ICTAS Doctoral Scholar
MILES-IGERT Trainee
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
Please don't post (un)subscribe requests to the list. Use the www
interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/mailing_lists/users.php




Dr. Jennifer Williams
Institute for Materials and Processes
School of Engineering
University of Edinburgh
Sanderson Building
The King's Buildings
Mayfield Road
Edinburgh, EH9 3JL, United Kingdom
Phone: ++44 (0)131 650 4 861


--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
Please don't post (un)subscribe requests to the list. Use the
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/mailing_lists/users.php


[gmx-users] using trjconv to get a 2x2x2 unit cell for traj.xtc

2009-12-08 Thread Jennifer Williams



Hello,

I am trying to find a way around a visualisation problem I am having  
in VMD. Some of my molecules go over periodic boundary conditions  
meaning that bonds sometimes appear missing when looking at movies (I  
am trying to fix this using wrap, unwrap and join in VMD but as yet no  
luck).


I was wondering if there is a way in gromacs to multiply the number of  
unit cells shown in a trajectory. i.e instead of a 1x1x1 I want the  
new unit cell to be 2x2x2. This would mean the section of the  
structure I want to zoom in doesn't go over the pbc.


For the confout.gro file I have done this using

genconf -nbox 2 2 2 -f confout.gro -o confbig.out

and this enables me to at least see a static image where all bonds are  
present.


but in order to view a movie, I need to carry out something similar on  
the traj.xtc file. I have seen that with trjconv there is the option


-box  Size for new cubic box

but my unit cell is not cubic, it is a parallelepiped. The cell  
dimensions are :


4.64210   3.77847   1.89596   0.0   0.0  -2.18150   0.0
0.0   0.0


I tried using this anyway with the following command:

trjconv -box 9.28414 8.72598 3.79192 -f traj.xtc -o trajbig.xtc

but the resulting .xtc file wouldn't load into VMD so I assume that  
the .xtc file and the .gro file didn't match.


Any ideas?,

Thanks

Jenny







--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
Please don't post (un)subscribe requests to the list. Use the
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/mailing_lists/users.php


[gmx-users] functional form of PRODRG dihedrals

2009-12-08 Thread Jennifer Williams

Hello,

I have a general question re the PRODRG server used to create  
topologies for gromacs.


Using PRODRG for a simple alkane I get:

[ dihedrals ]
; ai  aj  ak  al  fuc0, c1, m, ...
   4   3   2   1   1  0.05.9 3  0.05.9 3 ; dih   CAG  CAE  CAC

The dihedrals created are given the function fu = 1 so I assume that  
they are not Ryckaert-Bellemans functions (where fu=3).


What is the functional form of the dihedral given by PRODRG? Does it  
correspond to equation 4.6.1 on page 62 of the gromacs manual (v. 4.0)?


If so which numbers in the above dihedral correspond to the symbols in  
equation 4.6.1?


Thanks

Jenny







Dr. Jennifer Williams
Institute for Materials and Processes
School of Engineering
University of Edinburgh
Sanderson Building
The King's Buildings
Mayfield Road
Edinburgh, EH9 3JL, United Kingdom
Phone: ++44 (0)131 650 4 861


--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

Can't post? Read http://www.gromacs.org/mailing_lists/users.php


[gmx-users] Charge group moving too much

2009-11-10 Thread Jennifer Williams

Hello,

One of my simulations just crashed after 660,000 steps with the  
following message:


Step 660560:
The charge group starting at atom 4100 moved than the distance allowed  
by the domain decomposition (1.50) in direction X

distance out of cell 2.356293
Old coordinates:1.7520.4470.751
New coordinates:3.5931.5320.674
Old cell boundaries in direction X:0.0002.116
New cell boundaries in direction X:0.0002.122

---
Program mdrun, VERSION 4.0.5
Source code file: domdec.c, line: 3651

Fatal error:
A charge group moved too far between two domain decomposition steps
This usually means that your system is not well equilibrated
---

"Why Weren't You at My Funeral ?" (G. Groenhof)

Error on node 0, will try to stop all the nodes
Halting parallel program mdrun on CPU 0 out of 8

I am modelling flexible alkyl chains which are anchored to a frozen  
silica wall. This particular atom (4100) is at the edge of the  
periodic boundary condition so I am guessing that what happens is that  
while the anchored side stays where it is, the other free end wanders  
into the next cell and this perhaps triggers the error message? If  
this is indeed the problem, is there some way to stop the program  
stopping with an error message when this happens?


I did do an energy minimisation before starting and my structure looks  
OK. It seems strange that it would crash only after 66 steps.


Also what is meant by the following? I am using NVT so my simulation  
cell should not change.


Old cell boundaries in direction X:0.0002.116
New cell boundaries in direction X:0.0002.122


My cell is about 4.3 x 4.3 x 7.5 nm
I am using gromacs-4.0.5

Thanks




--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
Please don't post (un)subscribe requests to the list. Use the
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/mailing_lists/users.php


[gmx-users] Freezing a molecule

2009-10-27 Thread Jennifer Williams

Hello,

One quick questions...

I have a structure for which I now want to freeze a portion.
I already have a .top file where the entire structure is flexible  
(bond angles, stretches and torions defined).
When freezing, do I need to delete all those bond stretches, angles  
and torsions associated with the frozen part from my .top file?
(I am not using constraints on the frozen portion as I have seen this  
can cause problems).
Will the bonding parameters all be set to zero when the these atoms  
are frozen or will their contribution be calculated if I leave these  
terms in.



Thanks




--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.


___
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

Can't post? Read http://www.gromacs.org/mailing_lists/users.php


[gmx-users] Not all bonded interactions have been properly assigned to the domain decomposition cells

2009-10-27 Thread Jennifer Williams
serint4 = 0
   userreal1= 0
   userreal2= 0
   userreal3= 0
   userreal4= 0
grpopts:
   nrdf:5392
   ref_t: 300
   tau_t: 0.1
anneal:  No
ann_npoints:   0
   acc:0   0   0
   nfreeze:   Y   Y   Y   N
N   N

   energygrp_flags[  0]: 0
   efield-x:
  n = 0
   efield-xt:
  n = 0
   efield-y:
  n = 0
   efield-yt:
  n = 0
   efield-z:
  n = 0
   efield-zt:
  n = 0
   bQMMM= FALSE
   QMconstraints= 0
   QMMMscheme   = 0
   scalefactor  = 1
qm_opts:
   ngQM = 0

Initializing Domain Decomposition on 8 nodes
Dynamic load balancing: auto
Will sort the charge groups at every domain (re)decomposition

NOTE: Periodic molecules: can not easily determine the required  
minimum bonded cut-off, using half the non-bonded cut-off


Minimum cell size due to bonded interactions: 0.750 nm
Maximum distance for 5 constraints, at 120 deg. angles, all-trans: 0.376 nm
Estimated maximum distance required for P-LINCS: 0.376 nm
Using 0 separate PME nodes
Scaling the initial minimum size with 1/0.8 (option -dds) = 1.25
Optimizing the DD grid for 8 cells with a minimum initial size of 0.938 nm
The maximum allowed number of cells is: X 4 Y 4 Z 8
Domain decomposition grid 2 x 1 x 4, separate PME nodes 0
Domain decomposition nodeid 0, coordinates 0 0 0

Table routines are used for coulomb: TRUE
Table routines are used for vdw: TRUE
Will do PME sum in reciprocal space.

 PLEASE READ AND CITE THE FOLLOWING REFERENCE 
U. Essman, L. Perela, M. L. Berkowitz, T. Darden, H. Lee and L. G. Pedersen
A smooth particle mesh Ewald method
J. Chem. Phys. 103 (1995) pp. 8577-8592
  --- Thank You ---  

Using a Gaussian width (1/beta) of 0.480244 nm for Ewald
Using shifted Lennard-Jones, switch between 0.9 and 1.2 nm
Cut-off's:   NS: 1.5   Coulomb: 1.5   LJ: 1.2
System total charge: 0.000
Generated table with 2000 data points for Ewald.
Tabscale = 500 points/nm
Generated table with 2000 data points for LJ6Shift.
Tabscale = 500 points/nm
Generated table with 2000 data points for LJ12Shift.
Tabscale = 500 points/nm
Generated table with 2000 data points for 1-4 COUL.
Tabscale = 500 points/nm
Generated table with 2000 data points for 1-4 LJ6.
Tabscale = 500 points/nm
Generated table with 2000 data points for 1-4 LJ12.
Tabscale = 500 points/nm
Configuring nonbonded kernels...
Testing x86_64 SSE support... present.

Initializing Parallel LINear Constraint Solver

 PLEASE READ AND CITE THE FOLLOWING REFERENCE 
B. Hess
P-LINCS: A Parallel Linear Constraint Solver for molecular simulation
J. Chem. Theory Comput. 4 (2008) pp. 116-122
  --- Thank You ---  

The number of constraints is 800
There are inter charge-group constraints,
will communicate selected coordinates each lincs iteration

Linking all bonded interactions to atoms
There are 3236 inter charge-group exclusions,
will use an extra communication step for exclusion forces for PME

The initial number of communication pulses is: X 1 Z 1
The initial domain decomposition cell size is: X 2.01 nm Z 1.90 nm

The maximum allowed distance for charge groups involved in interactions is:
 non-bonded interactions   1.500 nm
two-body bonded interactions  (-rdd)   1.500 nm
  multi-body bonded interactions  (-rdd)   1.500 nm
  atoms separated by up to 5 constraints  (-rcon)  1.896 nm

When dynamic load balancing gets turned on, these settings will change to:
The maximum number of communication pulses is: X 1 Z 1
The minimum size for domain decomposition cells is 1.500 nm
The requested allowed shrink of DD cells (option -dds) is: 0.80
The allowed shrink of domain decomposition cells is: X 0.75 Z 0.79
The maximum allowed distance for charge groups involved in interactions is:
 non-bonded interactions   1.500 nm
two-body bonded interactions  (-rdd)   1.500 nm
  multi-body bonded interactions  (-rdd)   1.500 nm
  atoms separated by up to 5 constraints  (-rcon)  1.500 nm

Making 2D domain decomposition grid 2 x 1 x 4, home cell index 0 0 0

There are: 5244 Atoms
There are: 476 VSites
Charge group distribution at step 0: 583 565 583 565 666 684 666 684
Grid: 9 x 6 x 6 cells

Constraining the starting coordinates (step 0)

Constraining the coordinates at t0-dt (step 0)

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


Dr. Jennifer Williams
Institute for Materials and Processes
School of Engineering
University of Edinburgh
Sanderson Building
The King's Buildings
Mayfield Road
Edinburgh, EH9 3JL, United Kingdom
Phone: ++44 (0)131 650 4 861


--
The University of Edinburgh is a charitable body, registered in
Scotland, 

Re: [gmx-users] how to stop duplicate atoms from being deleted

2009-10-26 Thread Jennifer Williams


Hi Justin,

Thanks for the new suggestion.

However, wouldn't this again involve the sorting of my .pdb file?

I have modified the specbond.dat

3
MSI   SI4 MOO 2 0.16  MCM   MCM
MSI   SI4 MOH   OH2 0.16  MCM   MCM
MOH   OH2 MHH 1 0.101 MCM   MCM

BUT the specbond.dat is called after the duplicate atoms are deleted:

Now there are 4 atoms. Deleted 4280 duplicates.
Opening library file  
/home/jwillia4/gromacs-4.0.5/code/share/gromacs/top/specbond.dat

3 out of 3 lines of specbond.dat converted succesfully

So at present my pdb file is shortened down to 4 atoms before specdat  
even gets a chance to work on it.


The only solution I can see to this is to rename each and every atom  
in my .pdb file with a different residue number and name.  This then  
means that I would need a huge specbond.dat where the bonds were  
defined for each of these residues-this would essentially mean  
defining thousands of bonds as I wouldn't be able to define the bonds  
by atom type but by their individual reside number. Is this correct or  
am I missing something?


What I really would like is to stop the program deleting these  
duplicate atoms-then I could use pdb2gmx to generate the entire  
topology file. Do you know a way of doing this (and here I am really  
hoping that I don't have to mess with the source code :) )?


Thanks -I am learning a lot from your advice

Jenny






Quoting "Justin A. Lemkul" :




Justin A. Lemkul wrote:

I can see how this rapidly becomes difficult :)  I don't believe   
that pdb2gmx can handle such "multi-directional" bonding, since the  
 residues that are connected are not necessarily numerically   
sequential.


I should amend this statement (typing quicker than the brain can keep
up!) - It is not that pdb2gmx cannot handle multi-directional bonding,
it is moreso that I don't think it cannot be done using the +/-
convention of the .rtp files.

Using specbond.dat, as I suggested before, should be a viable alternative.

-Justin

--


Justin A. Lemkul
Ph.D. Candidate
ICTAS Doctoral Scholar
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin


___
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
Please don't post (un)subscribe requests to the list. Use the www
interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/mailing_lists/users.php




Dr. Jennifer Williams
Institute for Materials and Processes
School of Engineering
University of Edinburgh
Sanderson Building
The King's Buildings
Mayfield Road
Edinburgh, EH9 3JL, United Kingdom
Phone: ++44 (0)131 650 4 861


--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.


___
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
Please don't post (un)subscribe requests to the list. Use the
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/mailing_lists/users.php


Re: [gmx-users] how to stop duplicate atoms from being deleted

2009-10-23 Thread Jennifer Williams
  
which oxygens belong to which monomer. The only way I can see past  
this is a more elaborate naming system which would introduce yet more  
combinations.


So I?m throwing this out as a last resort before I give up. Has anyone  
had any success using this method for a similar system? Quartz?


Sorry for my rambling

Jenny



Quoting "Justin A. Lemkul" :




Jennifer Williams wrote:


Hi Justin,

Thanks for the reply. I am in fact studying one huge molecule. All   
of my atoms are bonded together in one large structure (kind of   
like a zeolite) so I have necessarily defined them as a single   
residue.




I would argue that you have a polymer, which can certainly be handled
by pdb2gmx.  See below.

There is no way I can split this molecule into smaller subunits and  
 thus define a number of residues-it wouldn't make sense to do so.




If you have a lot of repetition, I would think it would be quite easy
to split it apart.

Yes in my .rtp file I have only defined each atom type once. To   
define each and every atom in my one residue would mean defining   
4284 atoms!




If you have a repeating structure, you have a polymer, so you can just
decompose a repeat unit into a single .rtp entry.  That's the entire
purpose of pdb2gmx, we certainly wouldn't want to create an .rtp entry
for every single possible protein either!

For more information, see here:

http://www.gromacs.org/Documentation/How-tos/Polymers

I am having real trouble in creating topology files for my   
structure. At the moment, the only way I can do this is by using a   
tool in DL_POLY to create a field file and then manually change it   
to a .top file. This is really fiddely and I have a number of   
similar structures to do this for. I was hoping that I could do a   
similar step in Gromacs and get a .top file straight away-even if   
it means a bit more work setting it up.


Is there any hope or is pdb2gmx simply not designed to work for   
this sort of system?




You can certainly use pdb2gmx, it is intended to be versatile so it can
be used with any repeating structure of monomers, homogenous (like a
repeating polymer) or heterogenous (like a protein).  See the link
above.

-Justin


Thanks

Jenny


Quoting "Justin A. Lemkul" :


Quoting Jennifer Williams :



Hello

I am studying a mesoporous silica for which there is no topology in
gromacs-to try to automate the process of generating a topology file
(x2top doesn?t work), I am using pdb2gmx (or rather trying to).

I have parameters for my silica structure and have added a new section
for my molecule to the .rtp file, .atp file, atommass.dat,
atom_nom.dbl, nb.itp and bon.itp files.

The problem is that when I use my .pdb file to generate a topology,
pdb2gmx checks for duplicates and removes almost all of my atoms. It
leaves only one of each type. I should have a few hundred of each atom
type?here is the output from pdb2gmx?

Analyzing pdb file
There are 1 chains and 0 blocks of water and 1 residues with 4284 atoms
  chain  #res #atoms
  1 ' ' 1   4284
All occupancies are one

All ok up to here?and then?.

Processing chain 1 (4284 atoms, 1 residues)
There are 552 donors and 2580 acceptors
There are 1603 hydrogen bonds
Checking for duplicate atoms
Now there are 4 atoms. Deleted 4280 duplicates.

Can anyone explain why this is happening? ?none of my atoms have the
same coordinates. Is there a file that I have forgotten to alter?  Is
there is fix to turn off the checking of duplicate atoms? I don?t want
any of my atoms to be deleted!


You have all of your atoms defined within one residue.  I'm   
assuming  your .rtp
entry contains the definition of a single repeat unit, so each   
monomer should
be a separate residue.  The coordinates don't matter, it's because  
  within each
residue, you have the same atom names, so pdb2gmx removes them   
when it finds

them.



Below I paste an extract of my pdb file?



I'm assuming you'll have to probably reconstruct this file to   
re-organize the
atoms to define continuous residues.  It appears they are grouped   
by  atom name,

which is probably not what you want.

-Justin


CRYST1   46.421   43.630   75.838  90.00  90.00 120.00 P 1   1
ATOM  1  SI   MCM   1 -21.090  -1.951 -29.596  1.00  0.00  
  SI
ATOM  2  SI   MCM   1 -21.090  -1.951 -10.636  1.00  0.00  
  SI

??..
ATOM   1153  OMCM   1  20.602 -18.404 -20.904  1.00  0.00  
   O
ATOM   1154  OMCM   1  20.602 -18.404  -1.945  1.00  0.00  
   O

?
ATOM   3181  OH   MCM   1  -6.620 -18.769 -32.169  1.00  0.00
ATOM   3182  OH   MCM   1  -6.620 -18.769 -13.210  1.00  0.00
.
ATOM   3733  HMCM   1  -6.674 -18.381 -33.035  1.00  0.00  
   H
ATOM   3734  HMCM   1  -6.616 -18.600 -14.144  1.00  0.00  
   H


Any advice appreciated,

Thanks in advance



--
The University of Edinburgh is a char

Re: [gmx-users] how to stop duplicate atoms from being deleted

2009-10-23 Thread Jennifer Williams


Hi Justin,

Thanks for the reply. I am in fact studying one huge molecule. All of  
my atoms are bonded together in one large structure (kind of like a  
zeolite) so I have necessarily defined them as a single residue.


There is no way I can split this molecule into smaller subunits and  
thus define a number of residues-it wouldn't make sense to do so.


Yes in my .rtp file I have only defined each atom type once. To define  
each and every atom in my one residue would mean defining 4284 atoms!


I am having real trouble in creating topology files for my structure.  
At the moment, the only way I can do this is by using a tool in  
DL_POLY to create a field file and then manually change it to a .top  
file. This is really fiddely and I have a number of similar structures  
to do this for. I was hoping that I could do a similar step in Gromacs  
and get a .top file straight away-even if it means a bit more work  
setting it up.


Is there any hope or is pdb2gmx simply not designed to work for this  
sort of system?


Thanks

Jenny


Quoting "Justin A. Lemkul" :


Quoting Jennifer Williams :



Hello

I am studying a mesoporous silica for which there is no topology in
gromacs-to try to automate the process of generating a topology file
(x2top doesn?t work), I am using pdb2gmx (or rather trying to).

I have parameters for my silica structure and have added a new section
for my molecule to the .rtp file, .atp file, atommass.dat,
atom_nom.dbl, nb.itp and bon.itp files.

The problem is that when I use my .pdb file to generate a topology,
pdb2gmx checks for duplicates and removes almost all of my atoms. It
leaves only one of each type. I should have a few hundred of each atom
type?here is the output from pdb2gmx?

Analyzing pdb file
There are 1 chains and 0 blocks of water and 1 residues with 4284 atoms
   chain  #res #atoms
   1 ' ' 1   4284
All occupancies are one

All ok up to here?and then?.

Processing chain 1 (4284 atoms, 1 residues)
There are 552 donors and 2580 acceptors
There are 1603 hydrogen bonds
Checking for duplicate atoms
Now there are 4 atoms. Deleted 4280 duplicates.

Can anyone explain why this is happening? ?none of my atoms have the
same coordinates. Is there a file that I have forgotten to alter?  Is
there is fix to turn off the checking of duplicate atoms? I don?t want
any of my atoms to be deleted!


You have all of your atoms defined within one residue.  I'm assuming  
 your .rtp

entry contains the definition of a single repeat unit, so each monomer should
be a separate residue.  The coordinates don't matter, it's because   
within each

residue, you have the same atom names, so pdb2gmx removes them when it finds
them.



Below I paste an extract of my pdb file?



I'm assuming you'll have to probably reconstruct this file to re-organize the
atoms to define continuous residues.  It appears they are grouped by  
 atom name,

which is probably not what you want.

-Justin


CRYST1   46.421   43.630   75.838  90.00  90.00 120.00 P 1   1
ATOM  1  SI   MCM   1 -21.090  -1.951 -29.596  1.00  0.00
SI
ATOM  2  SI   MCM   1 -21.090  -1.951 -10.636  1.00  0.00
SI

??..
ATOM   1153  OMCM   1  20.602 -18.404 -20.904  1.00  0.00
 O
ATOM   1154  OMCM   1  20.602 -18.404  -1.945  1.00  0.00
 O

?
ATOM   3181  OH   MCM   1  -6.620 -18.769 -32.169  1.00  0.00
ATOM   3182  OH   MCM   1  -6.620 -18.769 -13.210  1.00  0.00
.
ATOM   3733  HMCM   1  -6.674 -18.381 -33.035  1.00  0.00
 H
ATOM   3734  HMCM   1  -6.616 -18.600 -14.144  1.00  0.00
 H


Any advice appreciated,

Thanks in advance



--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.


___
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
Please don't post (un)subscribe requests to the list. Use the
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/mailing_lists/users.php







Justin A. Lemkul
Graduate Research Assistant
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalem...@vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin/


___
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
Please don't post (un)subscribe requests to the list. Use the
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/mailing_lists/users.php






[gmx-users] how to stop duplicate atoms from being deleted

2009-10-23 Thread Jennifer Williams


Hello

I am studying a mesoporous silica for which there is no topology in  
gromacs-to try to automate the process of generating a topology file  
(x2top doesn?t work), I am using pdb2gmx (or rather trying to).


I have parameters for my silica structure and have added a new section  
for my molecule to the .rtp file, .atp file, atommass.dat,  
atom_nom.dbl, nb.itp and bon.itp files.


The problem is that when I use my .pdb file to generate a topology,  
pdb2gmx checks for duplicates and removes almost all of my atoms. It  
leaves only one of each type. I should have a few hundred of each atom  
type?here is the output from pdb2gmx?


Analyzing pdb file
There are 1 chains and 0 blocks of water and 1 residues with 4284 atoms
  chain  #res #atoms
  1 ' ' 1   4284
All occupancies are one

All ok up to here?and then?.

Processing chain 1 (4284 atoms, 1 residues)
There are 552 donors and 2580 acceptors
There are 1603 hydrogen bonds
Checking for duplicate atoms
Now there are 4 atoms. Deleted 4280 duplicates.

Can anyone explain why this is happening? ?none of my atoms have the  
same coordinates. Is there a file that I have forgotten to alter?  Is  
there is fix to turn off the checking of duplicate atoms? I don?t want  
any of my atoms to be deleted!


Below I paste an extract of my pdb file?

CRYST1   46.421   43.630   75.838  90.00  90.00 120.00 P 1   1
ATOM  1  SI   MCM   1 -21.090  -1.951 -29.596  1.00  0.00  SI
ATOM  2  SI   MCM   1 -21.090  -1.951 -10.636  1.00  0.00  SI
??..
ATOM   1153  OMCM   1  20.602 -18.404 -20.904  1.00  0.00   O
ATOM   1154  OMCM   1  20.602 -18.404  -1.945  1.00  0.00   O
?
ATOM   3181  OH   MCM   1  -6.620 -18.769 -32.169  1.00  0.00
ATOM   3182  OH   MCM   1  -6.620 -18.769 -13.210  1.00  0.00
.
ATOM   3733  HMCM   1  -6.674 -18.381 -33.035  1.00  0.00   H
ATOM   3734  HMCM   1  -6.616 -18.600 -14.144  1.00  0.00   H

Any advice appreciated,

Thanks in advance



--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.


___
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
Please don't post (un)subscribe requests to the list. Use the
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/mailing_lists/users.php


[gmx-users] msd not linear and c.o.m removal

2009-09-17 Thread Jennifer Williams

Hi,

I am running simulations of gaseous molecules (CH4, C2H6, CO2, N2) in  
a silica pore (4nm in diameter), cell dimensions, 46 x 44 x 37. This  
in an infinite pore in the z direction.


I am interested in looking at diffusion but I am a bit concerned about  
the shape of my MSD curves, for a 5ns run, they start off fairly  
linear and then at around 2.5ns the slope stops being linear and  
levels out - in some cases the gradient goes to zero. So I am getting  
a curve rather than a straight line. This is the same whether or not I  
use the ?type z option in my command.


At the moment I take the fit the curve between 0.5 - 2.5ns and ignore  
the second curvy half of the graph but I'd like to know why this is  
happening-Is it because the statistics become poorer for longer runs  
or is there something fundamentaly wrong with my system?. Does anyone  
know why this is happening?


Another thing I am confused about is whether or not to remove the  
centre of mass when one is interested in obtaining the diffusion  
coefficient.  I have seen on the gromacs forum that linear centre of  
mass removal is OK but a colleague has advised me not to remove it  
when studying diffusion. What about when applying an acceleration to  
obtain transport diffusion coefficients? Is c.o.m removal then needed  
or not?


Currently, I am not using any centre of mass removal in my mdrun. Then  
I calculate the Ds to be


g_msd ?mol -type z -f traj.xtc -s md.tpr -o msd_z.xvg
D[   CO2] 140.5914 (+/- 75.9583) 1e-5 cm^2/s

However when I include the ?rmcomm command when I run g_msd I get a  
very different Ds.


g_msd -rmcomm -type z -mol -f traj.xtc -s md.tpr -o msd_z.xvg
D[   CO2] 0.0118 (+/- 0.0364) 1e-5 cm^2/s

It is difficult to say which value is correct as we don't have  
experimental data for comparison.


Has anyone got any advice on this? I have pasted my typical input files below

Much appreciated,


C2H6.itp

[ atomtypes ]
;   typemasschargeptype   c6c12
CH315.034 0.00 A0.35121.16236218

[ moleculetype ]
; name  nrexcl
ET  2

[ atoms ]
;   nr  typeresnr   residu  atomcgnrcharge  mass
1   CH3 1   ET CH3  10.000  15.03452
2   CH3 1   ET CH3  10.000  15.03452

[ constraints ]
;  ai  aj funct   c0   c1
1   2   1  0.2353

[ exclusions ]
1 2
2 1

.mdp

; VARIOUS PREPROCESSING OPTIONS
; Preprocessor information: use cpp syntax.
; e.g.: -I/home/joe/doe -I/home/mary/hoe
include  = -I../top
; e.g.: -DI_Want_Cookies -DMe_Too
define   =

; RUN CONTROL PARAMETERS
integrator   = md
; Start time and timestep in ps
tinit= 0
dt   = 0.001
nsteps   = 500
; For exact run continuation or redoing part of a run
; Part index is updated automatically on checkpointing (keeps files separate)
simulation_part  = 1
init_step= 0
; mode for center of mass motion removal
comm-mode= none
; number of steps for center of mass motion removal
nstcomm  = 1
; group(s) for center of mass motion removal
comm-grps=

; LANGEVIN DYNAMICS OPTIONS
; Friction coefficient (amu/ps) and random seed
bd-fric  = 0
ld-seed  = 1993

; OUTPUT CONTROL OPTIONS
; Output frequency for coords (x), velocities (v) and forces (f)
nstxout  = 1000
nstvout  = 1000
nstfout  = 1000
; Output frequency for energies to log file and energy file
nstlog   = 1000
nstenergy= 1000
; Output frequency and precision for xtc file
nstxtcout= 1000
xtc-precision= 1000
; This selects the subset of atoms for the xtc file. You can
; select multiple groups. By default all atoms will be written.
xtc-grps =
; Selection of energy groups
energygrps   =

; NEIGHBORSEARCHING PARAMETERS
; nblist update frequency
nstlist  = 10
; ns algorithm (simple or grid)
ns_type  = grid
; Periodic boundary conditions: xyz, no, xy
pbc  = xyz
periodic_molecules   = yes
; nblist cut-off
rlist= 1.5

; OPTIONS FOR ELECTROSTATICS AND VDW
; Method for doing electrostatics
coulombtype  = PME
rcoulomb-switch  = 0
rcoulomb = 1.5
; Relative dielectric constant for the medium and the reaction field
epsilon-r= 1
epsilon_rf   = 1
; Method for doing Van der Waals
vdw-type = Shift
; cut-off lengths
rvdw-switch  = 0.9
rvdw = 1.2
; Apply long range dispersion corrections for Energy and Pressure
DispCorr = EnerPres
; Extension of the potential lookup tables beyond the cut-off
table-extension  = 1
; Seperate tables betw

[gmx-users] large scaling required to acheive optimal mesh load

2009-09-11 Thread Jennifer Williams

Hello user,

I am simulating a unit cell with dimensions 70x70x38 nm using PME. I  
started out with a cut-off of rvdw = rcoloumb = rlist=0.9 and a  
Spacing for the PME/PPPM FFT grid= 0.12, optimize fft = yes


I get the following output when I compile the .tpr file:

Using a fourier grid of 60x60x33, spacing 0.117 0.117 0.117
Estimate for the relative computational load of the PME mesh part: 0.97

NOTE 1 [file SMO_CO2.top, line 2159]:
  The optimal PME mesh load for parallel simulations is below 0.5
  and for highly parallel simulations between 0.25 and 0.33,
  for higher performance, increase the cut-off and the PME grid spacing

I did a number of test-runs increasing the cut-offs and the grid  
spacing by a factor of themselves. However I had to nearly double the  
cut-off and grid spacing in order to get the PME mesh load below 50.  
From the forum notes on the topic I got the impression that only a  
small scaling factor was needed.


My question is, are the values which I have achieved reasonable?

Cut-off: 1.665 and grid spacing 0.222

This is the output using these values

Checking consistency between energy and charge groups...
Calculating fourier grid dimensions for X Y Z
Using a fourier grid of 32x32x18, spacing 0.219 0.219 0.215
Estimate for the relative computational load of the PME mesh part: 0.38
This run will generate roughly 63 Mb of data
writing run input file...

Does changing these values have any effect on the results of the mdrun  
or only on the speed?


Thanks in advance,

Jenny



--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.


___
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
Please don't post (un)subscribe requests to the list. Use the
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/mailing_lists/users.php


[gmx-users] potential energy NaN and strange dependence on cut-offs

2009-09-08 Thread Jennifer Williams

Hi users,

I am running a very simple simulation of methane inside a pore (v.much  
like a carbon nanotube but in my case the tube is supposed to  
represent silica.) I keep this tube frozen.


I start with an energy minimisation-however this runs to completion  
almost instantly and I keep get NaN for my potential energy:


Steepest Descents converged to machine precision in 18 steps,
but did not reach the requested Fmax < 10.
Potential Energy  =nan
Maximum force =  6.5738518e+01 on atom 2133
Norm of force =  1.5461593e+00

Otherwise the trajectory looks OK (methane moving around inside the  
cylinder). If I go on to use the conf.gro file for an mdrun, it runs  
to completion and generates what looks like a reasonable trajectory,  
however the output again contains NaN i.e:


   Energies (kJ/mol)
LJ (SR)   Coulomb (SR)  PotentialKinetic En.   Total Energy
nan0.0e+00nan3.36749e+01nan
  Conserved En.Temperature Pressure (bar)
nan3.00010e+02nan

and calculating the Diffusion coefficient gives:
D[   CH4] 613.6682 (+/- 97.0563) 1e-5 cm^2/s

If I do the same calculation but reduce the cut-offs to 0.9. I get

   Energies (kJ/mol)
LJ (SR)   Coulomb (SR)  PotentialKinetic En.   Total Energy
nan0.0e+00nan3.36750e+01nan
  Conserved En.Temperature Pressure (bar)
nan3.00011e+02nan

D[   CH4] 237.8712 (+/- 53.5975) 1e-5 cm^2/s

And for a cut-off of 1.3nm I get

   Energies (kJ/mol)
LJ (SR)   Coulomb (SR)  PotentialKinetic En.   Total Energ
y
nan0.0e+00nan3.36737e+01na
n
  Conserved En.Temperature Pressure (bar)
nan2.9e+02nan


D[   CH4] 19.7953 (+/- 154.0168) 1e-5 cm^2/s


For this system, the cut-off shouldn?t need to be larger than 0.8 (I  
have plotted graphs of calculated V vs r) so it is worrying that the  
diffusion coefficient is showing such dependence on the cut-offs when  
they should all give the same result.


Can anyone offer any insight into this? I?ve tried changing the  
timestep making it both larger and smaller and many other things. I?ve  
pasted the relevant parts of my files below:


I?m using gromacs 4.0.5 ?at the moment running in serial.

Thanks for any advice,

Top file

[ defaults ]
; nbfunccomb-rule   gen-pairs   fudgeLJ fudgeQQ
1   2   yes 1.01.0
;
;
[ atomtypes ]
;   typemasschargeptype   c6c12
OSM15.99940.00 A 0.2708   1.538176

;
; Include forcefield parameters
#include "CH4.itp"
;
;
[ moleculetype ]
;   Namenrexcl
MCM 3
[ atoms ]
;   nr  typeresnr   residue atomcgnrcharge  mass
1   OSM 1   MCM OSM 1   0   15.9994
2   OSM 1   MCM OSM 2   0   15.9994
..etc
2127OSM 1   MCM OSM 21270   15.9994
2128OSM 1   MCM OSM 21280   15.9994


[ system ]
; Name
CH4 in MCM

[ molecules ]
; Compound#mols
MCM1
CH410

CH4.itp file

[ atomtypes ]
;   type  masschargeptype   c6c12
CH416.043 0.00 A0.37321.24650457
;
[ moleculetype ]
; name  nrexcl
CH42

[ atoms ]
;   nr  typeresnr   residu  atomcgnrcharge  mass
1   CH4  1   CH4 CH4 10.00  16.043



.mdp file

;
;   File 'mdout.mdp' was generated
;   By user: jwillia4 (353773)
;   On host: vlxhead2
;   At date: Fri Jun 26 15:47:37 2009
;
; VARIOUS PREPROCESSING OPTIONS
; Preprocessor information: use cpp syntax.
; e.g.: -I/home/joe/doe -I/home/mary/hoe
include  = -I../top
; e.g.: -DI_Want_Cookies -DMe_Too
define   =

; RUN CONTROL PARAMETERS
integrator   = steep
; Start time and timestep in ps
tinit= 0
dt   = 0.0001
nsteps   = 10
; For exact run continuation or redoing part of a run
; Part index is updated automatically on checkpointing (keeps files separate)
simulation_part  = 1
init_step= 0
; mode for center of mass motion removal
comm-mode= linear
; number of steps for center of mass motion removal
nstcomm  = 1
; group(s) for center of mass motion removal
comm-grps=

; LANGEVIN DYNAMICS OPTIONS
; Friction coefficient (amu/ps) and random seed
bd-fric  = 0
ld-seed  = 1993

; ENERGY MINIMIZATION OPTIONS
; Force tolerance and initial step-size
emtol=
emstep   = 0.001
; Max number of iterations in relax_shells
niter=
; Step size (ps^2) for minimization of flexibl

RE: [gmx-users] very strange domain composition statistics

2009-07-31 Thread Jennifer Williams

Quoting Jennifer Williams :



Hi,

Thanks for your input. Sorry I should have mentioned that I am using
the latest version of gromacs (4.0.5).

This morning I noticed that the strange domain decomp statistics are
only produced when my simulations run on certain nodes. Below I have
pasted the domain decomp statistics for the SAME .tpr file run on 2
different nodes (each time using 6 nodes)- the first looks ok to me
while the second produces crazy numbers.

I have opened a call to computer support at my Uni to ask if there is
some difference between the nodes and for now I specify that my
simulations should run on certain nodes which I have checked are ok.

I don't know if this is something to do with the way I compiled
gromacs or the architecture. I did a standard installation and it
seemed to run smoothly-no error messages. I have attached my config.log.

If you have any ideas please let me know,

Thanks



D O M A I N   D E C O M P O S I T I O N   S T A T I S T I C S

 av. #atoms communicated per step for force:  2 x 1967.2
 av. #atoms communicated per step for LINCS:  2 x 20.4

 Average load imbalance: 30.3 %
 Part of the total run time spent waiting due to load imbalance: 5.8 %
 Steps where the load balancing was limited by -rdd, -rcon and/or -dds: X 9 %

NOTE: 5.8 % performance was lost due to load imbalance
  in the domain decomposition.

 R E A L   C Y C L E   A N D   T I M E   A C C O U N T I N G

 Computing: Nodes Number G-CyclesSeconds %
---
 Domain decomp. 6 11  852.216  341.9 3.0
 Comm. coord.   6101  594.411  238.5 2.1
 Neighbor search6 11 2271.294  911.2 7.9
 Force  6101 5560.129 2230.519.3
 Wait + Comm. F 6101 1216.171  487.9 4.2
 PME mesh   610115071.432 6046.152.2
 Write traj.6   10012.5771.0 0.0
 Update 6101  264.545  106.1 0.9
 Constraints6101  923.418  370.4 3.2
 Comm. energies 6101  942.036  377.9 3.3
 Rest   61167.866  468.5 4.0
---
 Total  6   28866.09611580.0   100.0
---

Parallel run - timing based on wallclock.

   NODE (s)   Real (s)  (%)
   Time:   1930.000   1930.000100.0
   32:10
   (Mnbf/s)   (GFlops)   (ns/day)  (hour/ns)
Performance: 51.737  6.922 44.767  0.536
Finished mdrun on node 0 Fri Jul 31 13:05:01 2009


D O M A I N   D E C O M P O S I T I O N   S T A T I S T I C S

 av. #atoms communicated per step for force:  2 x 1969.0
 av. #atoms communicated per step for LINCS:  2 x 15.7

 Average load imbalance: 500.0 %
 Part of the total run time spent waiting due to load imbalance:
5822746112.0 %
 Steps where the load balancing was limited by -rdd, -rcon and/or -dds: X 9 %

NOTE: 5822746112.0 % performance was lost due to load imbalance
  in the domain decomposition.

 R E A L   C Y C L E   A N D   T I M E   A C C O U N T I N G

 Computing: Nodes Number G-CyclesSeconds %
---
 Write traj.6   1001 18443320128.89043061.8   100.0
 Update 6101  314.1750.0 0.0
 Rest   69223372036.85521534.950.0
---
 Total  618443397799.33643062.0   100.0
---

NOTE: 306 % of the run time was spent communicating energies,
  you might want to use the -nosum option of mdrun

Parallel run - timing based on wallclock.

   NODE (s)   Real (s)  (%)
   Time:   7177.000   7177.000100.0
   1h59:37
   (Mnbf/s)   (GFlops)   (ns/day)  (hour/ns)
Performance: 13.907  1.861 12.038  1.994
Finished mdrun on node 0 Thu Jul 30 01:47:05 2009










Quoting Berk Hess :






Date: Fri, 31 Jul 2009 07:49:49 +0200
From: sp...@xray.bmc.uu.se
To: gmx-users@gromacs.org
Subject: Re: [gmx-users] very strange domain composition statistics

Mark Abraham wrote:

Jennifer Williams wrote:

Hi ,

I am having some problems when running in parallel. Although my jobs
run to completion I am getting some worrying domain decomposition
statistics in particular the average load imbalance and the
performance loss due to load imbalance see below:


Please report your GROMACS version number. If it&

[gmx-users] very strange domain composition statistics

2009-07-30 Thread Jennifer Williams

Hi ,

I am having some problems when running in parallel. Although my jobs  
run to completion I am getting some worrying domain decomposition  
statistics in particular the average load imbalance and the  
performance loss due to load imbalance see below:


D O M A I N   D E C O M P O S I T I O N   S T A T I S T I C S

 av. #atoms communicated per step for force:  2 x 1974.8
 av. #atoms communicated per step for LINCS:  2 x 15.2

 Average load imbalance: 500.0 %
 Part of the total run time spent waiting due to load imbalance:  
4246403072.0 %

 Steps where the load balancing was limited by -rdd, -rcon and/or -dds: X 9 %

NOTE: 4246403072.0 % performance was lost due to load imbalance
  in the domain decomposition.

 R E A L   C Y C L E   A N D   T I M E   A C C O U N T I N G

 Computing: Nodes Number G-CyclesSeconds %
---
 Write traj.6   1001 18443320139.16442130.9   100.0
 Update 6101 18442922984.49142130.0   100.0
 Rest   69223372036.85521069.450.0
---
 Total  618446422611.66942138.0   100.0
---

NOTE: 305 % of the run time was spent communicating energies,
  you might want to use the -nosum option of mdrun

Parallel run - timing based on wallclock.

   NODE (s)   Real (s)  (%)
   Time:   7023.000   7023.000100.0
   1h57:03
   (Mnbf/s)   (GFlops)   (ns/day)  (hour/ns)
Performance: 14.214  1.902 12.302  1.951
Finished mdrun on node 0 Wed Jul 29 23:47:18 2009



Below is my .mdp file: I am using the PME but not having much of a  
feel for how to set the options under  Spacing for the PME/PPPM FFT  
grid, I left these as the default values. Could this be where the  
trouble lies?


My cut-off cannot be larger than 0.9 as my unit cell is only 18.2A in  
one direction.


How do I choose values for PME/PPPM? Ie what kind of values to use for  
nx, ny and nz ?
I read that they should be divisible by npme to get the best  
performance. Is npme the pme_order in the .mdp file? If not where do I  
set this parameter?


Much appreciated,

Jenny



; VARIOUS PREPROCESSING OPTIONS
; Preprocessor information: use cpp syntax.
; e.g.: -I/home/joe/doe -I/home/mary/hoe
include  = -I../top
; e.g.: -DI_Want_Cookies -DMe_Too
define   =

; RUN CONTROL PARAMETERS
integrator   = md
; Start time and timestep in ps
tinit= 0
dt   = 0.001
nsteps   = 100
; For exact run continuation or redoing part of a run
; Part index is updated automatically on checkpointing (keeps files separate)
simulation_part  = 1
init_step= 0
; mode for center of mass motion removal
comm-mode= linear
; number of steps for center of mass motion removal
nstcomm  = 1
; group(s) for center of mass motion removal
comm-grps=

; LANGEVIN DYNAMICS OPTIONS
; Friction coefficient (amu/ps) and random seed
bd-fric  = 0
ld-seed  = 1993

; ENERGY MINIMIZATION OPTIONS
; Force tolerance and initial step-size
emtol=
emstep   =
; Max number of iterations in relax_shells
niter=
; Step size (ps^2) for minimization of flexible constraints
fcstep   =
; Frequency of steepest descents steps when doing CG
nstcgsteep   =
nbfgscorr=

; TEST PARTICLE INSERTION OPTIONS
rtpi =

; OUTPUT CONTROL OPTIONS
; Output frequency for coords (x), velocities (v) and forces (f)
nstxout  = 1000
nstvout  = 1000
nstfout  = 0
; Output frequency for energies to log file and energy file
nstlog   = 1000
nstenergy= 1000
; Output frequency and precision for xtc file
nstxtcout= 1000
xtc-precision= 1000
; This selects the subset of atoms for the xtc file. You can
; select multiple groups. By default all atoms will be written.
xtc-grps =
; Selection of energy groups
energygrps   =

; NEIGHBORSEARCHING PARAMETERS
; nblist update frequency
nstlist  =
; ns algorithm (simple or grid)
ns_type  = grid
; Periodic boundary conditions: xyz, no, xy
pbc  = xyz
periodic_molecules   = yes
; nblist cut-off
rlist= 0.9

; OPTIONS FOR ELECTROSTATICS AND VDW
; Method for doing electrostatics
coulombtype  = PME
rcoulomb-switch  = 0
rcoulomb = 0.9
; Relative dielectric constant for the medium and the reaction field
epsilon_r 

Re: [gmx-users] Huge acceleration needed to reproduce results!

2009-07-28 Thread Jennifer Williams



Hi Chris,

Thanks for the input-its good to get another person's perspective.

Let me explain a bit more about my system...

MCM stands for MCM-41, a mesoporous silica which consits of 1071 atoms  
so even though it is only 1 molecule, it is the more massive group.


I am purposely accelerating the methane to try and calculate their  
diffusion thourgh the pore of the MCM.


As I understand it, the acceleration I specify should be applied to  
each and every molecule of CH4. I.e if I specify an acceleration of  
0.1 nm ps-2, each of my 340 molecules should experience an  
acceleration of 0.1. I.e the acceleration is not in any way divided up  
between the total number of accelerated molecules?



Jenny


Quoting chris.ne...@utoronto.ca:



Title indicates you think that you have CH4 in MCM:

[ system ]
; Name
CH4 in MCM


Actual system is MCM in CH4:

[ molecules ]
; Compound#mols
MCM1
CH4340


.mdp accelerates the CH4:

acc-grps = CH4


Now I'm not sure if this is the source of any problem, equal and
opposite being true, but you are certainly trying to accelerate the
more massive group and it seems like a strange thing to do. Could this
be it?

Chris.

-- original message --

Quoting Jennifer Williams :



Hi,

Yes I realised that gromacs works in ps. I converted my force in kj
mol-1 A-1 to acceleration in nm/ps2. I also took into account that the
msd.xvg is plotted in nm and ps-2 and the calculated gradient printed
at the top of the msd.xvg file is in cm2/s.

One strange thing that I do get is the message ?There were 228
inconsistent shifts. Check your topology? when I carry out the g_msd
with the ?mol option but not when I don?t use -mol. Why is this?

I also came across a forum post
(http://www.mail-archive.com/gmx-users@gromacs.org/msg5.html) that
said ?If the distance between two atoms is close to half the box, the
force may arbitrarily change sign. This is an ill-defined situation
for which there is no obvious solution.?  Could this somehow be
affecting my simulations?

Below are the relevant parts of my .mdp file and other files. If you
see something suspicious please let me know because I?m stuck,

Thanks


; RUN CONTROL PARAMETERS
integrator   = md
; Start time and timestep in ps
tinit= 0
dt   = 0.001
nsteps   = 100
; For exact run continuation or redoing part of a run
; Part index is updated automatically on checkpointing (keeps files  
 separate)

simulation_part  = 1
init_step= 0
; mode for center of mass motion removal
comm-mode= linear
; number of steps for center of mass motion removal
nstcomm  = 1
; group(s) for center of mass motion removal
comm-grps=


; OUTPUT CONTROL OPTIONS
; Output frequency for coords (x), velocities (v) and forces (f)
nstxout  = 1000
nstvout  = 1000
nstfout  = 0
; Output frequency for energies to log file and energy file
nstlog   = 1000
nstenergy= 1000
; Output frequency and precision for xtc file
nstxtcout= 1000
xtc-precision= 1000
; This selects the subset of atoms for the xtc file. You can
; select multiple groups. By default all atoms will be written.
xtc-grps =
; Selection of energy groups
energygrps   =

; NEIGHBORSEARCHING PARAMETERS
; nblist update frequency
nstlist  =
; ns algorithm (simple or grid)
ns_type  = grid
; Periodic boundary conditions: xyz, no, xy
pbc  = xyz
periodic_molecules   = yes
; nblist cut-off
rlist= 0.9

; OPTIONS FOR ELECTROSTATICS AND VDW
; Method for doing electrostatics
coulombtype  = PME
rcoulomb-switch  = 0
rcoulomb = 0.9
; Relative dielectric constant for the medium and the reaction field
epsilon_r=
epsilon_rf   =

; Method for doing Van der Waals
vdw-type = Cut-off
; cut-off lengths
rvdw-switch  = 0
rvdw = 0.9
; Apply long range dispersion corrections for Energy and Pressure
DispCorr = No
; Extension of the potential lookup tables beyond the cut-off
table-extension  =
; Seperate tables between energy group pairs
energygrp_table  =


; Spacing for the PME/PPPM FFT grid
fourierspacing   = 0.12
; FFT grid size, when a value is 0 fourierspacing will be used
fourier_nx   = 0
fourier_ny   = 0
fourier_nz   = 0
; EWALD/PME/PPPM parameters
pme_order=
ewald_rtol   = 1e-05
ewald_geometry   = 3d
epsilon_surface  = 0
optimize_fft = yes


; OPTIONS FOR WEAK COUPLING ALGORITHMS
; Temperature coupling
tcoupl   = nose-hoover
; Groups to couple separate

Re: [gmx-users] Huge acceleration needed to reproduce results!

2009-07-28 Thread Jennifer Williams


Hi,

Yes I realised that gromacs works in ps. I converted my force in kj  
mol-1 A-1 to acceleration in nm/ps2. I also took into account that the  
msd.xvg is plotted in nm and ps-2 and the calculated gradient printed  
at the top of the msd.xvg file is in cm2/s.


One strange thing that I do get is the message ?There were 228  
inconsistent shifts. Check your topology? when I carry out the g_msd  
with the ?mol option but not when I don?t use -mol. Why is this?


I also came across a forum post  
(http://www.mail-archive.com/gmx-users@gromacs.org/msg5.html) that  
said ?If the distance between two atoms is close to half the box, the  
force may arbitrarily change sign. This is an ill-defined situation  
for which there is no obvious solution.?  Could this somehow be  
affecting my simulations?


Below are the relevant parts of my .mdp file and other files. If you  
see something suspicious please let me know because I?m stuck,


Thanks


; RUN CONTROL PARAMETERS
integrator   = md
; Start time and timestep in ps
tinit= 0
dt   = 0.001
nsteps   = 100
; For exact run continuation or redoing part of a run
; Part index is updated automatically on checkpointing (keeps files separate)
simulation_part  = 1
init_step= 0
; mode for center of mass motion removal
comm-mode= linear
; number of steps for center of mass motion removal
nstcomm  = 1
; group(s) for center of mass motion removal
comm-grps=


; OUTPUT CONTROL OPTIONS
; Output frequency for coords (x), velocities (v) and forces (f)
nstxout  = 1000
nstvout  = 1000
nstfout  = 0
; Output frequency for energies to log file and energy file
nstlog   = 1000
nstenergy= 1000
; Output frequency and precision for xtc file
nstxtcout= 1000
xtc-precision= 1000
; This selects the subset of atoms for the xtc file. You can
; select multiple groups. By default all atoms will be written.
xtc-grps =
; Selection of energy groups
energygrps   =

; NEIGHBORSEARCHING PARAMETERS
; nblist update frequency
nstlist  =
; ns algorithm (simple or grid)
ns_type  = grid
; Periodic boundary conditions: xyz, no, xy
pbc  = xyz
periodic_molecules   = yes
; nblist cut-off
rlist= 0.9

; OPTIONS FOR ELECTROSTATICS AND VDW
; Method for doing electrostatics
coulombtype  = PME
rcoulomb-switch  = 0
rcoulomb = 0.9
; Relative dielectric constant for the medium and the reaction field
epsilon_r=
epsilon_rf   =

; Method for doing Van der Waals
vdw-type = Cut-off
; cut-off lengths
rvdw-switch  = 0
rvdw = 0.9
; Apply long range dispersion corrections for Energy and Pressure
DispCorr = No
; Extension of the potential lookup tables beyond the cut-off
table-extension  =
; Seperate tables between energy group pairs
energygrp_table  =


; Spacing for the PME/PPPM FFT grid
fourierspacing   = 0.12
; FFT grid size, when a value is 0 fourierspacing will be used
fourier_nx   = 0
fourier_ny   = 0
fourier_nz   = 0
; EWALD/PME/PPPM parameters
pme_order=
ewald_rtol   = 1e-05
ewald_geometry   = 3d
epsilon_surface  = 0
optimize_fft = yes


; OPTIONS FOR WEAK COUPLING ALGORITHMS
; Temperature coupling
tcoupl   = nose-hoover
; Groups to couple separately
tc-grps  = System
; Time constant (ps) and reference temperature (K)
tau_t= 0.1
ref_t= 150

; Pressure coupling
Pcoupl   = No
Pcoupltype   =
; Time constant (ps), compressibility (1/bar) and reference P (bar)
tau-p=
compressibility  =
ref-p=
; Scaling of reference coordinates, No, All or COM
refcoord_scaling = no
; Random seed for Andersen thermostat
andersen_seed=


; GENERATE VELOCITIES FOR STARTUP RUN
gen_vel  = yes
gen_temp = 150
gen_seed = 173529

; OPTIONS FOR BONDS
constraints  = none
; Type of constraint algorithm
constraint-algorithm = Lincs
; Do not constrain the start configuration
continuation = no
; Use successive overrelaxation to reduce the number of shake iterations
Shake-SOR= no
; Relative tolerance of shake
shake-tol= 0.0001
; Highest order in the expansion of the constraint coupling matrix
lincs-order  = 4
; Number of iterations in the final step of LINCS. 1 is fine for
; normal simulations, but use 2 to conserve energy in NVE runs.
; For energy minimization with constraints it sho

[gmx-users] Huge acceleration needed to reproduce results!

2009-07-27 Thread Jennifer Williams

Hello users,

I am trying to reproduce a calculation that I carried out in DL_POLY.  
It is to calculate the transport diffusion coefficient for CH4 in a  
frozen mesoporous silica.


In DL_POLY I used an external force of 0.1 kJ mol-1 A-1. (0.1 KJ per  
mole per angstrom). This equates to 10 dl_poly internal units which I  
add in this way at each timestep;


Fsum = Fsum + Fex

In Gromacs, I want to apply the same force as I used in DL_POLY so I  
calculated the required acceleration using F=ma. Where I took the mass  
to be the mass of one molecule of methane (16 a.m.u).


The final value for acceleration that I came up with (which  
corresponds to a force of 0.1kj mol-1 A-1 on each molecule) was 0.0625  
nm ps-2.


The first hiccup was when I used this value, the MSD was negative  
(though linear in the negative region of the graph). I assumed that  
this had something to do with the orientation of the unit cell and  
tried applying 0.0   0.0  -0.0625. The plot then looked much better.


The problem is when I calculate the Mean displacement of the CH4  
molecules. (I do this using a slightly altered version of the g_msd  
code). The Mean displacement  from gromacs is very different to that  
which I calculate using DL_POLY,


Gromacs gives 95.0, dl_poly 21347.0.

The MSD however (where I don?t add an acceleration are similar) so the  
problem lies with the force I am adding.


To test that it wasn?t some bug in my code to calculate the Mean  
displacement, I also looked at how the acceleration/force altered the  
MSD in DL_POLY and gromacs.


In DL_POLY, adding an external force of 0.1kj mol-1 A-1 would change  
the MSD of methane by 3 orders of magnitude compared to a run with no  
force added.


My equivalent acceleration of -0.0625 in gromacs, in comparison,  
barely changes the MSD from that of a run with no acceleration added.  
In fact it takes an acceleration of -200 in the z direction to cause  
such a difference in the Ds coefficients between runs with and without  
acceleration.


Does anyone have any idea what is going on here? An acceleration of  
200 ps nm-2 surely is not reasonable is it?. It seems very large  
compared to the example in the manual of 0.1. This would then imply  
that my back of the envelope calculation for relating force and  
acceleration is wrong. Am I missing something? I?m quite sure that the  
force I am adding in DL_POLY is equivalent to 0.1 kJ mol-1 A-1 so why  
are my methane molecules moving so much less in gromacs in response to  
the equivalent acceleration?


Also I noticed that although in the .mdp file I specify:

; Non-equilibrium MD stuff
acc-grps = CH4
accelerate   = 0.0 0.0 -200.00

In the md.log file I get the following output

acc:0   0-154.549   0   0 45.4514

Can someone clarify what this means?

Any advice/comments appreciated



--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.


___
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
Please don't post (un)subscribe requests to the list. Use the
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/mailing_lists/users.php


[gmx-users] questions on acceleration

2009-07-17 Thread Jennifer Williams

Hello,

I have 3 short questions regarding adding an acceleration.

1. I want to apply a force in the z direction to some of my molecules  
to mimic a chemical potential gradient.
I read an old gmx_users post where someone wanted to add an external  
force to their molecules and the advice given was to add a piece of  
code after the function efield in the code in src/mdlib/sim_util.c .  
That is, to make a function similar to efield which would do Fsum =  
Fsum + Fext


Is there any reason for doing this rather than using the accelerate  
option in the .mdp file?


2.   If I do use the accelerate option in the .mdp file, I need to  
convert the acceleration applied into a force. Just to check?. I  
assume that if I use F=ma, where F and a are the force and  
acceleration applied to a molecule, then m should be the mass of that  
molecule? I.e for a CO2 molecule, m will be 44 a.m.u?


3.  If I use the accelerate option under NEMD, does it still make  
sense to remove the C.O.M?



Thanks


--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.


___
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
Please don't post (un)subscribe requests to the list. Use the
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/mailing_lists/users.php


[gmx-users] Changing gmx_msd.c

2009-07-17 Thread Jennifer Williams

Hi gmx_users,

I am trying to use NEMD to get transport diffusion coefficients.

I am adding an acceleration to my atoms in the z direction and then I  
need to calculate this quantity:


1/t <  ? [r(t) ? r(0)] >

My plan was to alter the code that calculates the MSD (gmx_msd.c  
located in src/tools). This calculates;


1/t < 1/N. ? [r(t) ? r(0)]^2 >

I presumed that all I needed to do was locate and remove the squared  
function and the division by the number of molecules (1/N).


I would also change the factor for the conversion of nm^2/ps to 10e-5  
cm^2/s since I would now be doing nm/ps to 10e-5 cm/s.


However my problem is that I have never worked with C++ so changing  
the code involved some guess-work.


For the removal of the square I did;

r2   += rv[m]*rv[m];   > r2   += rv[m];   (line 252)


r2 += r*r; > r2 += r;(line 270)

to remove the division by the number of molecules (1/N), I removed line

  g/=nx; (line 279)

I realise that these changes are only valid for the first case CASE  
(static_real_calc_norm), which I presumed was for the basic  
calculation invoked by this command;


g_msd ?f traj.xtc ?s md.tpr

i.e. I don?t plan on using ?pbc ?mol etc so I didn?t change the code  
for these cases. (I also realise that unless I change the code, the  
calculated error wont be valid-but I was going to worry about that  
later!)


However, the changes I made do not change the output from that of the  
original MSD code. I am using the correct, changed executable and not  
an old one (i.e I copied the compiled executable to gromacs/bin). So  
it looks like I am not changing the correct part of the code, or  
perhaps the command line g_msd ?f traj.xtc ?s md.tpr is not changing  
the case I altered?


Has anyone else changed the code for this purpose or would anyone be  
able to give me some pointers on what I need to alter in the code.  
Today is the first time I have ever looked at any C++ code (though I  
know some fortran) so please be as specific as possible :)



Much appreciated



--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.


___
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
Please don't post (un)subscribe requests to the list. Use the
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/mailing_lists/users.php


RE: [gmx-users] Re: Best way to handle linear rigid molecules.

2009-07-15 Thread Jennifer Williams
M_CO2.top, line 8230]:
   atom C_CO (Res CO2-1) has mass 0


ERROR 2 [file MCM_CO2.top, line 8230]:
   atom O_CO (Res CO2-1) has mass 0



ERROR 3 [file MCM_CO2.top, line 8230]:
   atom O_CO (Res CO2-1) has mass 0



ERROR 4 [file MCM_CO2.top, line 8230]:
   virtual site DUM (Res CO2-1) has non-zero mass 22.005
Check your topology.



ERROR 5 [file MCM_CO2.top, line 8230]:
   virtual site DUM (Res CO2-1) has non-zero mass 22.005
Check your topology.

So the approach of modelling dummy atoms with a mass and the real
atoms as mass-less generates errors. How do I solve this? Importantly
I need to get the MSD of the CO2 molecules so I am a bit confused with
the approach of using virtual atom sites to model CO2. Which groups
would I then select for the MSD? The virtual site with the mass or the
real atoms site which is mass-less?  Is this approach realistic?

I have also tried to model CO2 2 different approaches

1.  Using bond constraints and setting the bond angle to 180 and using
a large bending constant to keep the molecule effectively rigid.

Here I got the following error message and the run stopped after 2   
timesteps;


Step 1, time 0.001 (ps)  LINCS WARNING
relative constraint deviation after LINCS:
rms 0.969889, max 10.563440 (between atoms 1072 and 1074)
bonds that rotated more than 30 degrees:
  atom 1 atom 2  angle  previous, current, constraint length
1072   1073   89.90.1210   0.6876  0.1209
1072   1074   92.40.1210   1.3978  0.1209
1075   1076   86.80.1209   0.1303  0.1209
1075   1077   86.80.1209   0.1307  0.1209
Wrote pdb files with previous and current coordinates


2.  Using constraints to keep the C-O bond fixed to 0.12 and another
constraint between O and O of 0.24. I hoped this would keep the O-C-O
linear.

Here, the md.log file contained the following;


6 constraints are involved in constraint triangles,
will apply an additional matrix expansion of order 4 for couplings
between constraints inside triangles

The job ran to completion but on viewing the trajectory, the CO2
molecules are only approximately linear.

If anyone has an example file for a triatomic linear molecule like
CO2, I?d be very grateful if you could post it so I could use it as a
template,

I am using gromacs version 4.0.5

Thanks











Hi,

This question has been asked several times before on this list.
The proper way to handle this in Gromacs is to introduce two dummy
mass particles
in such a way the total mass and the moment of inertia is identical to CO2.
Then you can construct the positions of the three mass-less C and O  
 atoms from

the positions of the masses using virtual sites (see the manual).

Berk



Quoting Jennifer Williams :

> Hello users,
>
> I am using gromacs v 4.0.5
>
> What is the best way to treat linear triatomic molecules such as CO2.
> Is there some kind of rigid algorthym as in DL_POLY?
>
> In the asbsence of a ?rigid? algorthym. I initially intended to ?fix?
> the C=O bond length using constraints (as below) and somehow ?fix? the
> O=C=O bond angle to 180.
>
>
> [ atomtypes ]
> ;   type  masschargeptype   c6c12
>C_CO212.011150.5406  A  0.  0.
>O_CO215.9994-0.2703  A  0.29847 1.10765301
>
> [ moleculetype ]
> ; name  nrexcl
> CO2 2
>
> [ atoms ]
> ;   nr  typeresnr   residu  atomcgnrcharge mass
> 1   C_CO2 1   CO2C_CO2  10.5406  12.01115
> 2   O_CO2 1   CO2O_CO2  1   -0.2703  15.9994
> 3   O_CO2 1   CO2O_CO2  1   -0.2703  15.9994
>
> [ constraints ]
> ;  ai  aj funct   c0   c1
> 1   2  1  0.12088
> 1   3  1  0.12088
>
> [ angles ]
> ;  aiajak   funct   c0  c1
> 2 1 3   1   180.00
>
> First question... How do I go about fixing this angle to 180degrees?
> Should I just make c1 extremely large? Is there a more elegant way of
> going about this? Someone also mentioned using a virtual site
> representation for the C. Any comments on this welcome.
>
> I read some posts which imply that the angle bending code apparently can't
> handle angles of exactly 180 degrees.
>
> http://www.mail-archive.com/gmx-users@gromacs.org/msg09979.html
>
> This was a while ago so perhaps this has been dealt with in the current
> version of gromacs?
>
>
> Thanks



--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.


___
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
Please don&

[gmx-users] Re: Best way to handle linear rigid molecules.

2009-07-14 Thread Jennifer Williams

Hi gmx users

I am still having problems modelling small linear molecules. Following  
advice from the forum (see below), I tried to define CO2 using 2 dummy  
atoms (with averaged masses) which I placed in the same positions as  
my oxygen atoms. (This involves having to replicate coordinates for O  
in my pdb file which will be fiddly when I have hundreds of CO2  
molecules).


Here is my .itp file for CO2;


[ atomtypes ]
;   type  masschargeptype   c6c12
   C_CO0.000 0.5406   A  0.  0.
   O_CO0.000-0.2703   A  0.29847 1.10765301
   DUM22.004975  0.   V  0.  0.


[ moleculetype ]
; name  nrexcl
CO2 2

[ atoms ]
;   nr  typeresnr   residu  atomcgnrcharge  mass
1   C_CO 1   CO2C_CO  10.5406  0.000
2   O_CO 1   CO2O_CO  1   -0.2703  0.000
3   O_CO 1   CO2O_CO  1   -0.2703  0.000
4   DUM  1   CO2DUM   10.  22.004975
5   DUM  1   CO2DUM   10.  22.004975


[ constraints ]
;  ai  aj funct   c0   c1
2   3   1  0.24176

[ virtualsites2]
;  aiajak   funct   c0  c1
4 2 3   1   0.
5 2 3   1   0.24176

I am not sure what the final terms in the virtual sites section stand  
for and I have not found a lot of detail in the manual (section 4.7  
and 5.5.2). I assumed that the first term defines the virtual atom  
site, the second and third are the positions of the real atoms which  
the virtual atom is placed in relation to. What does the function ?1?  
stand for?  Here I tried to place the virtual atoms in the same  
positions as my ?real? oxygens- I am not sure if this will be allowed?


However, when I run this I get:


ERROR 1 [file MCM_CO2.top, line 8230]:
  atom C_CO (Res CO2-1) has mass 0


ERROR 2 [file MCM_CO2.top, line 8230]:
  atom O_CO (Res CO2-1) has mass 0



ERROR 3 [file MCM_CO2.top, line 8230]:
  atom O_CO (Res CO2-1) has mass 0



ERROR 4 [file MCM_CO2.top, line 8230]:
  virtual site DUM (Res CO2-1) has non-zero mass 22.005
   Check your topology.



ERROR 5 [file MCM_CO2.top, line 8230]:
  virtual site DUM (Res CO2-1) has non-zero mass 22.005
   Check your topology.

So the approach of modelling dummy atoms with a mass and the real  
atoms as mass-less generates errors. How do I solve this? Importantly  
I need to get the MSD of the CO2 molecules so I am a bit confused with  
the approach of using virtual atom sites to model CO2. Which groups  
would I then select for the MSD? The virtual site with the mass or the  
real atoms site which is mass-less?  Is this approach realistic?


I have also tried to model CO2 2 different approaches

1.	Using bond constraints and setting the bond angle to 180 and using  
a large bending constant to keep the molecule effectively rigid.


Here I got the following error message and the run stopped after 2 timesteps;

Step 1, time 0.001 (ps)  LINCS WARNING
relative constraint deviation after LINCS:
rms 0.969889, max 10.563440 (between atoms 1072 and 1074)
bonds that rotated more than 30 degrees:
 atom 1 atom 2  angle  previous, current, constraint length
   1072   1073   89.90.1210   0.6876  0.1209
   1072   1074   92.40.1210   1.3978  0.1209
   1075   1076   86.80.1209   0.1303  0.1209
   1075   1077   86.80.1209   0.1307  0.1209
Wrote pdb files with previous and current coordinates


2.	Using constraints to keep the C-O bond fixed to 0.12 and another  
constraint between O and O of 0.24. I hoped this would keep the O-C-O  
linear.


Here, the md.log file contained the following;


6 constraints are involved in constraint triangles,
will apply an additional matrix expansion of order 4 for couplings
between constraints inside triangles

The job ran to completion but on viewing the trajectory, the CO2  
molecules are only approximately linear.


If anyone has an example file for a triatomic linear molecule like  
CO2, I?d be very grateful if you could post it so I could use it as a  
template,


I am using gromacs version 4.0.5

Thanks











Hi,

This question has been asked several times before on this list.
The proper way to handle this in Gromacs is to introduce two dummy  
mass particles

in such a way the total mass and the moment of inertia is identical to CO2.
Then you can construct the positions of the three mass-less C and O atoms from
the positions of the masses using virtual sites (see the manual).

Berk



Quoting Jennifer Williams :


Hello users,

I am using gromacs v 4.0.5

What is the best way to treat linear triatomic molecules such as CO2.
Is there some kind of rigid algorthym as in DL_POLY?

In the asbsence of a ?rigid? algorthym. I initially intended to ?fix?
the C=O bond length using constraints (as below) and somehow ?fix? the
O=C=O

[gmx-users] Best way to handle linear rigid molecules.

2009-07-10 Thread Jennifer Williams

Hello users,

I am using gromacs v 4.0.5

What is the best way to treat linear triatomic molecules such as CO2.  
Is there some kind of rigid algorthym as in DL_POLY?


In the asbsence of a ?rigid? algorthym. I initially intended to ?fix?  
the C=O bond length using constraints (as below) and somehow ?fix? the  
O=C=O bond angle to 180.



[ atomtypes ]
;   type  masschargeptype   c6c12
   C_CO212.011150.5406  A   0.  0.
   O_CO215.9994-0.2703  A   0.29847 1.10765301

[ moleculetype ]
; name  nrexcl
CO2 2

[ atoms ]
;   nr  typeresnr   residu  atomcgnrcharge  mass
1   C_CO2 1   CO2C_CO2  10.5406  12.01115
2   O_CO2 1   CO2O_CO2  1   -0.2703  15.9994
3   O_CO2 1   CO2O_CO2  1   -0.2703  15.9994

[ constraints ]
;  ai  aj funct   c0   c1
1   2   1  0.12088
1   3   1  0.12088

[ angles ]
;  aiajak   funct   c0  c1
2 1 3   1   180.00

First question... How do I go about fixing this angle to 180degrees?   
Should I just make c1 extremely large? Is there a more elegant way of  
going about this? Someone also mentioned using a virtual site  
representation for the C. Any comments on this welcome.


I read some posts which imply that the angle bending code apparently can't
handle angles of exactly 180 degrees.

http://www.mail-archive.com/gmx-users@gromacs.org/msg09979.html

This was a while ago so perhaps this has been dealt with in the  
current version of gromacs?



Thanks


--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.


___
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
Please don't post (un)subscribe requests to the list. Use the
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/mailing_lists/users.php


[gmx-users] pdb2gmx adding a residue to the .rtp file

2009-07-08 Thread Jennifer Williams

Hi,

I am trying to use the pdb2gmx to generate a .top file from my .pdb.
I am getting the following error?.

All occupancies are one
Opening library file /home/jwillia4/gromacs/share/gromacs/top/ffoplsaa.atp
Atomtype 1
Reading residue database... (ffoplsaa)
Opening library file /home/jwillia4/gromacs/share/gromacs/top/ffoplsaa.rtp
Residue 56
---
Program pdb2gmx, VERSION 4.0.5
Source code file: resall.c, line: 138

Fatal error:
Atom type opls_23 (residue MCM) not found in atomtype database


As my residue is not already present I have added the following to the  
end of the ffoplsaa.rtp file



[ MCM ]
 [ atoms ]
OHopls_23-0.70  1
HOopls_24 0.435 1
SI   SI   0.832 1
OSopls_62-0.400 1
 [ bonds ]
SIOS
SIOH
OHHO
 [ dihedrals ]
OSSIOHHO
SIOSSIOS
SIOSSIOH
OHSIOHHO


The following is already present in the ffoplsaabon.itp file (standard  
version);


; This bit added by DvdS for quartz simulation 2005-05-07
; These parameters are taken from GROMOS and were also used in
; Wensink et al. Langmuir 16 (2000) pp. 7392-7400
;
[ bondtypes ]
; ij  func   b0  kb
  SIOS  10.16300 251040.  ; From GROMACS
  SIOH  10.16300 251040.  ;

[ angletypes ]
;  ijk  func   th0   cth
   SI   OS   SI1   155.000 397.480
   OS   SI   OS1   109.500 397.480
   OH   SI   OS1   109.500 397.480
   SI   OH   HO1   109.500 397.480

[ dihedraltypes ]
;  ijkl   func coefficients
; Added DvdS for Quartz simulations
   SI   OS1 0.000   3.766  3
   SI   OH1 0.000   3.766  3

And the entries in the .atp file which I am getting the error message  
about are already there (they are in the standard version-I didn?t  
need to add them myself). So I can?t figure out why they are not being  
read.


I also tried adding the name of my residue (MCM) to aminoacids.dat and  
incrementing the count at the top.


Do I need to make any other changes?

Thanks


--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.


___
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
Please don't post (un)subscribe requests to the list. Use the
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/mailing_lists/users.php


[gmx-users] Why is triclinic geometry not retained in confout.gro?

2009-07-07 Thread Jennifer Williams


Hello gmx-users

I?m having problems with my mdrun (probably a newbie question). I?m  
using the latest version of gromacs.


I have my .pdb file (or .gro), mdp and .top files ready and can  
generate the .tpr file. When I run this (with an simple energy  
minimisation) the simulation runs to completion (I?ve only tried a  
very short run)  but I get a strange confout.gro file as output. My  
unit cell is triclinic and my input files (both .pdb or .gro) looks  
fine in VMD.


The confgro.out file is strange in that the box has been converted to  
cubic when I view it in VMD. Is this normal? Why doesn?t it retain the  
triclinic shape I defined in the pdb file?


I?ve been over the topology but my inexperienced eyes can?t see  
anything wrong. One thing I did notice is that when I looked at the  
tpr file, all of my atom numbers were shifted by 1 with regards to my  
topology file i.e


In the tpr file, the first angle listed is for atoms 364, 0 and 413

  Angle:
 nr: 9492
 iatoms:
0 type=11 (ANGLES) 364 0 413

But in the topology file the first angle I have listed is for 365, 1 and 414.

[ angles ]
;   ai  aj  ak  funct   c0  c1
365 1   414 1   109.04  289.095916

In my pdb and top file my atoms are labelled from 1-1071 whereas in  
the .tpr file they are labelled from 0-1070. Is this something I  
should be worried about?


Below I have pasted sections of my top file, my pdb file and my .mdp  
file. I?d appreciate if someone could look over and see that my  
triclinic unit cell is correctly defined (although the input file  
looks OK when viewed in VMD).


If anyone has the time or inclination to try and run my files (if that  
helps spot the error), I would be happy to e-mail them and would be  
very grateful,


If you see anything else that looks odd, please feel free to point it  
out as I am new to gromacs,


Thanks in advance,

Jenny






pdb

CRYST1   46.421   43.630   18.960  90.00  90.00 120.00 P 1   1
ATOM  1  SI X   1 -22.104  -1.646  -1.173  1.00  0.00 SI
ATOM  2  SI X   1   8.325 -18.877   8.329  1.00  0.00  SI
ATOM  3  SI X   1  27.146 -12.854   3.831  1.00  0.00  SI
ATOM  4  SI X   1 -14.415 -11.322   4.375  1.00  0.00  SI
ATOM  5  SI X   1 -10.624 -15.731  -7.960  1.00  0.00  SI
...
ATOM289  O  X   1  19.588 -18.099   7.519  1.00  0.00
ATOM290  O  X   1 -19.450   0.838  -2.667  1.00  0.00
...
ATOM794  O  X   1  22.966 -15.478  -8.908  1.00  0.00
ATOM795  O  X   1  17.234  -5.878  -2.785  1.00  0.00
ATOM796  OH X   1  -7.634 -18.464  -3.746  1.00  0.00
ATOM797  H  X   1  -7.655 -18.213  -4.662  1.00  0.00   H
ATOM798  OH X   1   7.669 -17.509   7.819  1.00  0.00
ATOM799  H  X   1   8.061 -17.122   7.046  1.00  0.00   H

ATOM   1068  OH  X   1   4.808 -14.731   1.210  1.00  0.00
ATOM   1069  H   X   1   3.887 -14.515   1.123  1.00  0.00   H
ATOM   1070  OH  X   1  18.839   0.763  -8.266  1.00  0.00
ATOM   1071  H   X   1  18.283   0.835  -9.032  1.00  0.00   H
END

Topology file

[ defaults ]
; nbfunccomb-rule   gen-pairs   fudgeLJ fudgeQQ
1   2   no  1.01.0
;
;
[ atomtypes ]
;   typemasschargeptype   c6c12
SI 28.08  1.28 A0.000 0.00
O  15.999-0.64 A0.27081.538176
OH 15.999-0.53 A0.30  1.538176
H   1.008 0.21 A0.000 0.000
;
[ moleculetype ]
;   Namenrexcl
MCM 3
[ atoms ]
;   nr  typeresnr   residue atomcgnrcharge  mass
1   SI  1   MCM SI  1   1.2804993   28.086
2   SI  1   MCM SI  2   1.2804993   28.086
3   SI  1   MCM SI  3   1.2804993   28.086
4   SI  1   MCM SI  4   1.2804993   28.086
5   SI  1   MCM SI  5   1.2804993   28.086
...
289 O   1   MCM O   289 -0.64024965 15.9994
290 O   1   MCM O   290 -0.64024965 15.9994
...
794 O   1   MCM O   794 -0.64024965 15.9994
795 O   1   MCM O   795 -0.64024965 15.9994
796 OH  1   MCM OH  796 -0.52612471 15.9994
797 H   1   MCM H   797 0.20599988  1.00797
798 OH  1   MCM OH  798 -0.52612471 15.9994
799 H   1   MCM H   799 0.20599988  1.00797
...

1068   

[gmx-users] Conversion of units to nm in .xyz or pdb

2009-06-29 Thread Jennifer Williams

hello gromacs users,

I am working with a silica structure which is parallepiped in  
structure, i.e the unit cell is

46.4207   37.7846  18.9596 angles =90, 90, 120

I have a structure file in .xyz and .pdb but the units are in  
Angstroms and for gromacs I assume that I need to convert these units  
to nm? This is where I am encountering problems.


I have tried simply dividing the atomic coordinates and cell lengths  
by a factor of 10 but the resulting structure doesnt look right in  
vmd. Do I need to use some sort of symmetry transformation to take the  
non-cubic symmetry into account when scaling the coordinates? Is there  
some program (either within gromacs or external) which will convert  
atomic coordinates from Angstroms to nm while carrying out the  
necessary symmetry operations?



Thanks





--
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.


___
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at http://www.gromacs.org/search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

Can't post? Read http://www.gromacs.org/mailing_lists/users.php