<html>
<head>
<style>
.hmmessage P
{
margin:0px;
padding:0px
}
body.hmmessage
{
font-size: 10pt;
font-family:Verdana
}
</style>
</head>
<body class='hmmessage'>
Hi,<br><br>This is quite wrong.<br><br>Your normal atoms are virtual sites defined from the two masses<br>(you did it the other way around).<br>If you put the masses at the oxygen postitions (which gives an incorrect<br>moment of intertia, so incorrect dynamics, but correct thermodynamics),<br>then your coefficients are 0.5 for C and 0 and 1 for the O's.<br><br>Berk<br><br>> Date: Tue, 14 Jul 2009 14:39:10 +0100<br>> From: Jennifer.Williams@ed.ac.uk<br>> To: gmx-users@gromacs.org<br>> Subject: [gmx-users] Re: Best way to handle linear rigid molecules.<br>> <br>> Hi gmx users<br>> <br>> I am still having problems modelling small linear molecules. Following <br>> advice from the forum (see below), I tried to define CO2 using 2 dummy <br>> atoms (with averaged masses) which I placed in the same positions as <br>> my oxygen atoms. (This involves having to replicate coordinates for O <br>> in my pdb file which will be fiddly when I have hundreds of CO2 <br>> molecules).<br>> <br>> Here is my .itp file for CO2;<br>> <br>> <br>> [ atomtypes ]<br>> ; type mass charge ptype c6 c12<br>> C_CO 0.0000000 0.5406 A 0.0000 0.00000000<br>> O_CO 0.0000000 -0.2703 A 0.29847 1.10765301<br>> DUM 22.004975 0.0000 V 0.0000 0.00000000<br>> <br>> <br>> [ moleculetype ]<br>> ; name nrexcl<br>> CO2 2<br>> <br>> [ atoms ]<br>> ; nr type resnr residu atom cgnr charge mass<br>> 1 C_CO 1 CO2 C_CO 1 0.5406 0.0000000<br>> 2 O_CO 1 CO2 O_CO 1 -0.2703 0.0000000<br>> 3 O_CO 1 CO2 O_CO 1 -0.2703 0.0000000<br>> 4 DUM 1 CO2 DUM 1 0.0000 22.004975<br>> 5 DUM 1 CO2 DUM 1 0.0000 22.004975<br>> <br>> <br>> [ constraints ]<br>> ; ai aj funct c0 c1<br>> 2 3 1 0.24176<br>> <br>> [ virtualsites2]<br>> ; ai aj ak funct c0 c1<br>> 4 2 3 1 0.0000<br>> 5 2 3 1 0.24176<br>> <br>> I am not sure what the final terms in the virtual sites section stand <br>> for and I have not found a lot of detail in the manual (section 4.7 <br>> and 5.5.2). I assumed that the first term defines the virtual atom <br>> site, the second and third are the positions of the real atoms which <br>> the virtual atom is placed in relation to. What does the function ?1? <br>> stand for? Here I tried to place the virtual atoms in the same <br>> positions as my ?real? oxygens- I am not sure if this will be allowed?<br>> <br>> However, when I run this I get:<br>> <br>> <br>> ERROR 1 [file MCM_CO2.top, line 8230]:<br>> atom C_CO (Res CO2-1) has mass 0<br>> <br>> <br>> ERROR 2 [file MCM_CO2.top, line 8230]:<br>> atom O_CO (Res CO2-1) has mass 0<br>> <br>> <br>> <br>> ERROR 3 [file MCM_CO2.top, line 8230]:<br>> atom O_CO (Res CO2-1) has mass 0<br>> <br>> <br>> <br>> ERROR 4 [file MCM_CO2.top, line 8230]:<br>> virtual site DUM (Res CO2-1) has non-zero mass 22.005<br>> Check your topology.<br>> <br>> <br>> <br>> ERROR 5 [file MCM_CO2.top, line 8230]:<br>> virtual site DUM (Res CO2-1) has non-zero mass 22.005<br>> Check your topology.<br>> <br>> So the approach of modelling dummy atoms with a mass and the real <br>> atoms as mass-less generates errors. How do I solve this? Importantly <br>> I need to get the MSD of the CO2 molecules so I am a bit confused with <br>> the approach of using virtual atom sites to model CO2. Which groups <br>> would I then select for the MSD? The virtual site with the mass or the <br>> real atoms site which is mass-less? Is this approach realistic?<br>> <br>> I have also tried to model CO2 2 different approaches<br>> <br>> 1. Using bond constraints and setting the bond angle to 180 and using <br>> a large bending constant to keep the molecule effectively rigid.<br>> <br>> Here I got the following error message and the run stopped after 2 timesteps;<br>> <br>> Step 1, time 0.001 (ps) LINCS WARNING<br>> relative constraint deviation after LINCS:<br>> rms 0.969889, max 10.563440 (between atoms 1072 and 1074)<br>> bonds that rotated more than 30 degrees:<br>> atom 1 atom 2 angle previous, current, constraint length<br>> 1072 1073 89.9 0.1210 0.6876 0.1209<br>> 1072 1074 92.4 0.1210 1.3978 0.1209<br>> 1075 1076 86.8 0.1209 0.1303 0.1209<br>> 1075 1077 86.8 0.1209 0.1307 0.1209<br>> Wrote pdb files with previous and current coordinates<br>> <br>> <br>> 2. Using constraints to keep the C-O bond fixed to 0.12 and another <br>> constraint between O and O of 0.24. I hoped this would keep the O-C-O <br>> linear.<br>> <br>> Here, the md.log file contained the following;<br>> <br>> <br>> 6 constraints are involved in constraint triangles,<br>> will apply an additional matrix expansion of order 4 for couplings<br>> between constraints inside triangles<br>> <br>> The job ran to completion but on viewing the trajectory, the CO2 <br>> molecules are only approximately linear.<br>> <br>> If anyone has an example file for a triatomic linear molecule like <br>> CO2, I?d be very grateful if you could post it so I could use it as a <br>> template,<br>> <br>> I am using gromacs version 4.0.5<br>> <br>> Thanks<br>> <br>> <br>> <br>> <br>> <br>> <br>> <br>> <br>> <br>> <br>> <br>> Hi,<br>> <br>> This question has been asked several times before on this list.<br>> The proper way to handle this in Gromacs is to introduce two dummy <br>> mass particles<br>> in such a way the total mass and the moment of inertia is identical to CO2.<br>> Then you can construct the positions of the three mass-less C and O atoms from<br>> the positions of the masses using virtual sites (see the manual).<br>> <br>> Berk<br>> <br>> <br>> <br>> Quoting Jennifer Williams <Jennifer.Williams@ed.ac.uk>:<br>> <br>> > Hello users,<br>> ><br>> > I am using gromacs v 4.0.5<br>> ><br>> > What is the best way to treat linear triatomic molecules such as CO2.<br>> > Is there some kind of rigid algorthym as in DL_POLY?<br>> ><br>> > In the asbsence of a ?rigid? algorthym. I initially intended to ?fix?<br>> > the C=O bond length using constraints (as below) and somehow ?fix? the<br>> > O=C=O bond angle to 180.<br>> ><br>> ><br>> > [ atomtypes ]<br>> > ; type mass charge ptype c6 c12<br>> > C_CO2 12.01115 0.5406 A 0.0000 0.00000000<br>> > O_CO2 15.9994 -0.2703 A 0.29847 1.10765301<br>> ><br>> > [ moleculetype ]<br>> > ; name nrexcl<br>> > CO2 2<br>> ><br>> > [ atoms ]<br>> > ; nr type resnr residu atom cgnr charge mass<br>> > 1 C_CO2 1 CO2 C_CO2 1 0.5406 12.01115<br>> > 2 O_CO2 1 CO2 O_CO2 1 -0.2703 15.9994<br>> > 3 O_CO2 1 CO2 O_CO2 1 -0.2703 15.9994<br>> ><br>> > [ constraints ]<br>> > ; ai aj funct c0 c1<br>> > 1 2 1 0.12088<br>> > 1 3 1 0.12088<br>> ><br>> > [ angles ]<br>> > ; ai aj ak funct c0 c1<br>> > 2 1 3 1 180.00<br>> ><br>> > First question... How do I go about fixing this angle to 180degrees?<br>> > Should I just make c1 extremely large? Is there a more elegant way of<br>> > going about this? Someone also mentioned using a virtual site<br>> > representation for the C. Any comments on this welcome.<br>> ><br>> > I read some posts which imply that the angle bending code apparently can't<br>> > handle angles of exactly 180 degrees.<br>> ><br>> > http://www.mail-archive.com/gmx-users@gromacs.org/msg09979.html<br>> ><br>> > This was a while ago so perhaps this has been dealt with in the current<br>> > version of gromacs?<br>> ><br>> ><br>> > Thanks<br>> <br>> <br>> <br>> -- <br>> The University of Edinburgh is a charitable body, registered in<br>> Scotland, with registration number SC005336.<br>> <br>> <br>> _______________________________________________<br>> gmx-users mailing list gmx-users@gromacs.org<br>> http://lists.gromacs.org/mailman/listinfo/gmx-users<br>> Please search the archive at http://www.gromacs.org/search before posting!<br>> Please don't post (un)subscribe requests to the list. Use the <br>> www interface or send it to gmx-users-request@gromacs.org.<br>> Can't post? Read http://www.gromacs.org/mailing_lists/users.php<br><br /><hr />See all the ways you can stay connected <a href='http://www.microsoft.com/windows/windowslive/default.aspx' target='_new'>to friends and family</a></body>
</html>