[gmx-users] Segmentation Fault (Address not mapped)

[Justin A. Lemkul]

Darrell Koskinen wrote:
> Hi Mark,
> Thanks for your comments.
> 
> I selected ffoplsaa as I was told in a response to one of my posts that
> OPLS would work fine for air type simulations:
> **************************************************************************************** 
> 
>>> Also, I thought that ffoplsaa was to be used for liquids and since my
>  >> simulation is in an air environment, I thought that I should be using
>  >> the encad force field in a vacuum and should therefore be modifying
>  >> an encad n2t file. Please correct me if I am wrong.
> Encad is not made for vacuum. You have to verify your force field
> anyway, for publication purposes. OPLS will work just fine though.
> 
> (The above was extracted from
> http://oldwww.gromacs.org/pipermail/gmx-users/2009-April/041070.html)
> **************************************************************************************** 
> 

OPLS may work just fine, but as a number of us have been saying (including that 
post) "verify your force field."  You have to prove that a condensed-phase force 
field for use with biological molecules accurately represents the species in 
your system.  It may be that it works just fine, but you have to convince a 
scientific audience of that fact.

> 
> You indicate that I need to read the primary literature for the force
> fields. Could you please tell me where this literature can be found? I
> have searched using Google for the forcefields available for use in
> GROMACS (G43b1, G43a1, G43a2, G45a3, G53a5, G53a6, gmx), but have not
> found any primary literature for these force fields.
> 

The literature is out there, perhaps you're not using the right search terms. 
The first hit from searching "Gromos96 53a6" in Google is a relevant reference. 
The G43a1 reference is the GROMOS manual, which is not publicly available, but 
many (if not all) of the later force fields (45a3, 53a5, and 53a6) have had 
their parameters published.  53a5/6 were published in JCC in 2004.  In fact, 
pdb2gmx prints out a few references :)

> With regard to your comment about demonstrating that this force field is
> doing a good job at representing ammonia in the gas phase, how can I
> determine the diffusion of the pure ammonia gas? And what exactly do you

g_msd calculates diffusion constants.

> mean by a density constant? Do you mean the density of the gas? If so,
> then the density is defined by the number of molecules and box size I
> define in my .gro file, is it not? Also, the graphene is in an ammonia

In a system at equilibrium, yes.  You'd have to run an NPT simulation, achieve 
equilibrium, and determine the density.

> vapour, not in a vacuum, and the graphene structure is not infinite but
> has dimensions that fit well inside the simulation box. How would I
> determine if my graphene structure is stable?
> 

What are the structural features of this graphene lattice?  You'd probably have 
to demonstrate that your parameters accurately represent its true structure.

-Justin

> Thanks again.
> 
> Darrell
> 
>> Date: Tue, 04 Aug 2009 12:14:07 +1000
>> From: Mark Abraham <Mark.Abraham at anu.edu.au>
>> Subject: Re: [gmx-users] Segmentation Fault (Address not mapped)
>> To: Discussion list for GROMACS users <gmx-users at gromacs.org>
>> Message-ID: <4A77996F.9070106 at anu.edu.au>
>> Content-Type: text/plain; charset=windows-1252; format=flowed
>>
>> Darrell Koskinen wrote:
>>  
>>> Hi Justin,
>>> I just want to use the most appropriate force field and associated 
>>> parameters for my simulation of graphene in ammonia gas. Is OPLS an 
>>> appropriate force field to use?
>>>     
>>
>> It seems fairly unlikely that it is appropriate. Force fields are 
>> parameterized to reproduce a range of behaviour of a certain classes 
>> of chemicals. To understand what these are, you need to read the 
>> primary literature for that force field. There is no reason to suppose 
>> that a force field parameterized to reproduce (say) condensed-phase 
>> properties of peptides should be effective at modelling either 
>> infinite graphene, or gas-phase ammonia. It might work, and similar 
>> things might have been shown to work, but that's your job to prove 
>> and/or read about.
>>
>>  
>>> I have searched for force fields for gas simulations but have not yet 
>>> been able to find one specific to graphene and ammonia gas. Since 
>>> graphene is in a condensed phase and ammonia is in a gas phase and I 
>>> understand that I need to use the same force field for both graphene 
>>> and ammonia, I will have to use a force field for either gas-phase MD 
>>> or condensed-phase MD. Since I am finding it difficult to locate 
>>> gas-phase MD force fields, is it acceptable to use a condensed-phase 
>>> MD force field, which I believe is the case for OPLS (please confirm)?
>>>     
>>
>> Only if you can demonstrate that it's doing a good job. For starters, 
>> does a simulation of gas-phase ammonia using these parameters 
>> reproduce density and diffusion constants? Is an infinite graphene 
>> lattice in vacuum stable? Otherwise any paper you wrote might be 
>> rejected out-of-hand.
>>
>>  
>>> If it is acceptable to use OPLS, then I am assuming that I am to use 
>>> the parameters in the ffoplsaabon.itp and ffoplsaanb.itp as is and 
>>> not modify them (i.e. not create frankenstein force fields), correct? 
>>> Am I correct in my understanding that the only item that requires 
>>> modification is the .n2t file?
>>>     
>>
>> I can't say.
>>
>>  
>>> And, by the way, I have resolved the segmentation fault problem. It 
>>> was being caused by the freezing of the graphene lattice. When I 
>>> removed the freezing, I saw the graphene structure curl upward and am 
>>> thinking that by freezing the atoms in place it was creating great 
>>> forces between the atoms as they were kept from moving into their 
>>> preferred orientation.
>>>     
>>
>> Yep, that's right. The segfault was probably indirectly a consequence 
>> of integrating large opposing forces on these atoms. A small 
>> perturbation in a position might now result in a large change in net 
>> force, a large displacement on a frozen atom, and the whole thing goes 
>> up in silicaceous smoke.
>>
>>  
>>> I subsequently tried to restrain the atoms via position restraints, 
>>> but I still see the graphene structure moving. Could you please let 
>>> me know why the position restraints do not appear to be working? Do I 
>>> need to use greater restraint parameters? I created a position 
>>> restaint file for     
>>
>> Yes and no. Position restraints are good for giving things a gentle 
>> nudge to stay where you've put them. The textbook usage is that while 
>> equilibrating a starting structure into an ensemble, you don't want it 
>> to move much. Then once you're happy with your equilibrium ensemble, 
>> then you remove the restraints and see what really happens. However, 
>> if you're using them to smash something flat that really doesn't want 
>> to be flat, then you run the risk of integration problems as above. 
>> The fact of not wanting to be flat indicates that the model physics is 
>> severely broken, as Dallas has observed. You need to work out how and 
>> why it is broken, and once you fix the real problem, restraints or 
>> freezing might not prove to be needed.
>>
>> Mark
>>
>>  
>>> the graphene structure posre.itp using "genpr -f graphene.gro -n 
>>> index.ndx -o posre.itp -fc 0 0 1000" (note that genrestr is not 
>>> recognized in my version of GROMACS) to restrain the graphene lattice 
>>> to the x-y plane and submitted 2 runs with the following in my 
>>> run.mdp file:
>>>
>>> Run 1:
>>> define =-DPOSRES
>>> constraints =none
>>>
>>> Run 2:
>>> define =-DPOSRES
>>> constraints =all-bonds
>>>
>>> Run 1 with constraints=none appears to restrain the graphene lattice 
>>> better than the Run 2 with constraints=all-bonds.
>>>
>>> Here is an extract from posre.itp for the first few atoms:
>>>
>>> ******************************************
>>> ; position restraints for Grph of Ãzú^K^Q+
>>>
>>> [ position_restraints ]
>>> ; i funct fcx fcy fcz
>>> 1 1 0 0 1000
>>> 2 1 0 0 1000
>>> 3 1 0 0 1000
>>> ******************************************
>>>
>>> Again, I thank you very much for your assistance.
>>>
>>> Darrell
>>>
>>>    
>>>> Date: Mon, 27 Jul 2009 17:48:02 -0400
>>>> From: "Justin A. Lemkul" <jalemkul at vt.edu>
>>>> Subject: Re: [gmx-users] Segmentation Fault (Address not mapped)
>>>> To: Discussion list for GROMACS users <gmx-users at gromacs.org>
>>>> Message-ID: <4A6E2092.402 at vt.edu>
>>>> Content-Type: text/plain; charset=windows-1252; format=flowed
>>>>
>>>>
>>>>
>>>> Darrell Koskinen wrote:
>>>>  
>>>>      
>>>>> Hi Justin,
>>>>> With regard to your comment about using other force fields for my 
>>>>> simulation of graphene surrounded by ammonia gas, are referring to 
>>>>> the force fields for both graphene and ammonia or only to the force 
>>>>> field for graphene?
>>>>>
>>>>>             
>>>> The "force field" should be one in the same.  The parameters should 
>>>> be suitably derived using the same scheme as the original work.
>>>>
>>>>  
>>>>      
>>>>> I reviewed my selection of the force field parameters for graphene 
>>>>> and see that I selected the parameters from the paper by Cornell et 
>>>>> al. since a paper on deformation of carbon nanotubes ("A structural 
>>>>> mechanics approach for the analysis of carbon nanotubes" by Chunyu 
>>>>> Li, Tsu-Wei Chou in International Journal of Solids and Structures 
>>>>> 40 (2003) 2487–2499) used parameters from the paper by Cornell et 
>>>>> al. Is this sufficient to justify the use of these parameters? I 
>>>>> also thought that graphene and ammonia would be considered organic 
>>>>> since they are comprised of carbon, nitrogen, and hydrogen, which 
>>>>> are all common elements found in organic matter.
>>>>>
>>>>>             
>>>> If you feel that precedent is sufficient, then I guess go ahead.  
>>>> But realize that the paper by Cornell et al. refers to parameters 
>>>> suitable for simulations of proteins and nucleic acids, as well as a 
>>>> few organic functional groups. Parameterization was based on peptide 
>>>> backbone geometry, as well as other parameters, likely none of which 
>>>> involved graphene and gaseous ammonia (they were doing liquid 
>>>> simulations).
>>>>
>>>>  
>>>>      
>>>>> With regard to the parameters for ammonia, is it acceptable to use 
>>>>> the parameters from the paper by Cornell et al. or do I need to 
>>>>> find papers specifically related to ammonia gas molecular dynamics 
>>>>> simulations?
>>>>>
>>>>>             
>>>> I would seriously consider finding parameters (if they exist) that 
>>>> have been derived for use with gas-phase simulations.
>>>>
>>>>  
>>>>      
>>>>> With regard to dihedral selection, I looked at the .top file and 
>>>>> see that the function type is listed as 3 in the dihedrals section, 
>>>>> which I believe indicates that my simulation is using a 
>>>>> Ryckaert-Bellemans function. It appears to me that if I am using 
>>>>> the OPLS force field that the Ryckaert-Bellemans dihedral type is 
>>>>> automatically selected. I do see in the ffoplsaabon.itp file that 
>>>>> there is an improper dihedral definition which think I could 
>>>>> potentially use in my simulation "improper_Z_CA_X_Y", but how do I 
>>>>> cause my simulation to use this definition? Do I need to add the 
>>>>> following line to my .top file?
>>>>> "#define improper_Z_CA_X_Y"
>>>>>
>>>>>             
>>>> No, you would have to use a special dihedrals section, that 
>>>> specifies the four atoms involved (function type 1), followed by the 
>>>> specification of that particular improper ("improper_Z_CA_X_Y").  
>>>> See the manual for more details, and/or generate a topology for a 
>>>> protein with pdb2gmx with OPLS-AA and see how they're defined.
>>>>
>>>>  
>>>>      
>>>>> Further, I thought that I could simply modify the appropriate lines 
>>>>> of the ffoplsaanb.itp and ffoplsaabon.itp to enter in the 
>>>>> parameters found in the paper by Cornell et al. so that they will 
>>>>> be used in my simulation. Am I correct in my assumption?
>>>>>
>>>>>             
>>>> Modifying ffoplsaa files to parameters that are not OPLS-AA is 
>>>> probably a bad choice.  If you're dead-set on using these 
>>>> parameters, download the ffamber ports and use the ffamber94 force 
>>>> field, not some Frankensteined-OPLS.  Besides, all of the parameters 
>>>> in the Cornell paper are in kcal/mol/A^2, which will require 
>>>> conversion to Gromacs standard units.  If you've printed them 
>>>> directly in your previous attempts, I can guarantee you you're not 
>>>> getting the values you think you are.
>>>>
>>>> -Justin
>>>>
>>>>  
>>>>      
>>>>> Thanks again for your help.
>>>>>
>>>>> Darrell
>>>>>
>>>>>
>>>>>  >Date: Fri, 17 Jul 2009 22:11:07 -0400
>>>>>  >From: "Justin A. Lemkul" <jalemkul at vt.edu>
>>>>>  >Subject: Re: [gmx-users] Segmentation Fault (Address not mapped)
>>>>>  >To: Discussion list for GROMACS users <gmx-users at gromacs.org>
>>>>>  >Message-ID: <4A612F3B.7060503 at vt.edu>
>>>>>  >Content-Type: text/plain; charset=ISO-8859-1; format=flowed
>>>>>  >
>>>>>  >
>>>>>  >
>>>>>  >darrellk at ece.ubc.ca wrote:
>>>>>  >> Hi Justin,
>>>>>  >> I am using the ffoplsaabon.itp parameters for:
>>>>>  >> bonds (CA-CA, CA-CT, CA-HA, CT-HC, H-N3)
>>>>>  >> angles (CA-CA-CA, CA-CA-HA, CA-CT-HC, H-N3-H)
>>>>>  >> dihedrals (X-CA-CA-X, X-CA-CT-X)
>>>>>  >>
>>>>>  >> I am using the ffolpsaanb.itp VdW parameters for H & N3
>>>>>  >>
>>>>>  >> And I am using the VdW parameters found in the paper "A Second
>>>>>  >> Generation Force Field for the Simulation of Proteins, Nucleic 
>>>>> Acids,
>>>>>  >> and Organic Molecules" for CA, HA, CT, & HC doing the 
>>>>> conversion from
>>>>>  >> Angstroms (paper units) to nm (ffoplsaanb.itp units) for sigma 
>>>>> and from
>>>>>  >> kcal/mol (paper units) to kJ/mol (ffoplsaanb.itp units) for 
>>>>> epsilon. The
>>>>>  >> converted parameters for sigma were different by only +/-2% to 
>>>>> +/-8%
>>>>>  >> from the values in ffoplsaanb.itp, but the converted parameters 
>>>>> for
>>>>>  >> epsilon were different by +/-23% to +/-66% from the values in
>>>>>  >> ffoplsaanb.itp. Is this a valid paper to use for selection of 
>>>>> parameters?
>>>>>  >>
>>>>>  >
>>>>>  >If memory serves, that paper is the derivation by Cornell, et al. 
>>>>> for the
>>>>>  >AMBER94 parameter set, so no, you are not going to find those 
>>>>> same parameters in
>>>>>  >ffoplsaanb.itp.
>>>>>  >
>>>>>  >Furthermore, the title of the paper indicates that the parameters 
>>>>> are for
>>>>>  >proteins, nucleic acids, and organic molecules. Has anyone else 
>>>>> in the
>>>>>  >literature used these parameters for graphene lattices? I would 
>>>>> suggest using a
>>>>>  >force field that others have used for such solid materials, like 
>>>>> nanotubes or
>>>>>  >other graphene surfaces. Otherwise, you will have to somehow 
>>>>> demonstrate that a
>>>>>  >force field designed for use with condensed phase biomolecules is 
>>>>> applicable to
>>>>>  >your gas-phase graphene lattice.
>>>>>  >
>>>>>  >> The graphene lattice is only vibrating slightly and looks like 
>>>>> a bunch of
>>>>>  >> travelling waves.
>>>>>  >>
>>>>>  >> I believe I am using proper dihedrals and not improper 
>>>>> dihedrals in my
>>>>>  >> model.
>>>>>  >>
>>>>>  >
>>>>>  >It would seem. Consider if impropers might be necessary to keep 
>>>>> your planar
>>>>>  >groups planar. That's what they are there for.
>>>>>  >
>>>>>  >> Since the graphene structure will be connected to electrodes at 
>>>>> both ends
>>>>>  >> and will be mounted above a substrate, I assume the atoms in the
>>>>>  >> structure are not completely free and thus I need to model them 
>>>>> being
>>>>>  >> restricted in their movement. Maybe position restraints are more
>>>>>  >> appropriate than freezing to model such a situation.
>>>>>  >>
>>>>>  >> I do not see anything in the trajectory that provides clues as 
>>>>> to what
>>>>>  >> might be causing the segmentation fault.
>>>>>  >>
>>>>>  >
>>>>>  >Then it will be very difficult to get any more remote help ;) If you
>>>>>  >consistently get the crash between step 20,000 and 30,000, then 
>>>>> perhaps split
>>>>>  >your simulation into shorter sections, and during the time when 
>>>>> you expect the
>>>>>  >crash to occur, set nstxtcout = (some small value, like 1 or 10) 
>>>>> to obtain the
>>>>>  >most detail possible. At some point, something is going to go 
>>>>> careening off
>>>>>  >into infinity.
>>>>>  >
>>>>>  >-Justin
>>>>>  >
>>>>>  >> Thanks again for your help.
>>>>>  >>
>>>>>  >> Darrell
>>>>>  >>
>>>>>  >>> Date: Fri, 17 Jul 2009 19:08:38 -0400
>>>>>  >>> From: "Justin A. Lemkul" <jalemkul at vt.edu>
>>>>>  >>> Subject: Re: [gmx-users] Segmentation Fault (Address not mapped)
>>>>>  >>> To: Discussion list for GROMACS users <gmx-users at gromacs.org>
>>>>>  >>> Message-ID: <4A610476.2090709 at vt.edu>
>>>>>  >>> Content-Type: text/plain; charset=ISO-8859-1; format=flowed
>>>>>  >>>
>>>>>  >>>
>>>>>  >>>
>>>>>  >>> darrellk at ece.ubc.ca wrote:
>>>>>  >>>> Hi Justin,
>>>>>  >>>> I froze the graphene sheet because, prior to freezing it, I 
>>>>> noticed that
>>>>>  >>>> it was vibrating and thought that maybe its vibration was not 
>>>>> allowing
>>>>>  >>>> the NH3 molecules to adsorb (bond) to it. But after freezing the
>>>>>  >>>> graphene sheet, I see that see that NH3 molecules are still 
>>>>> not bonding
>>>>>  >>>> to it. Physical experiments of NH3 and a graphene lattice 
>>>>> connected to
>>>>>  >>>> electrodes have shown that NH3 does adsorb to graphene, but 
>>>>> all I see
>>>>>  >>>> are NH3 molecules coming close to the graphene surface and 
>>>>> then bouncing
>>>>>  >>>> away which I am assuming is a result of repulsion between the 
>>>>> negatively
>>>>>  >>>> charged N atom in the ammonia molecule and the pi electrons 
>>>>> in the
>>>>>  >>>> graphene lattice. So I am not sure why experiments have shown 
>>>>> adsorption
>>>>>  >>>> unless adsoption is occurring as a result of a current 
>>>>> flowing through
>>>>>  >>>> the graphene structure or as a result of edge effects at the 
>>>>> interface
>>>>>  >>>> between the electrodes and the graphene lattice.
>>>>>  >>>>
>>>>>  >>> What parameters are you using for the species involved? Could 
>>>>> be that there's
>>>>>  >>> something wrong with your model.
>>>>>  >>>
>>>>>  >>>> Could you tell me how freezing is different that position 
>>>>> restraining as
>>>>>  >>>> this is not completely clear to me?
>>>>>  >>>>
>>>>>  >>> Freezing means positions are absolutely fixed and never 
>>>>> updated. Position
>>>>>  >>> restraining means there is an energy penalty to movement, but 
>>>>> positions can
>>>>>  >>> adjust slightly. I was thinking that if you have some strange 
>>>>> interaction
>>>>>  >>> happening, and freezing was preventing any reaction to the 
>>>>> force, then you'd see
>>>>>  >>> an explosion. I don't necessarily see why you need to apply 
>>>>> either, but that
>>>>>  >>> will depend on the extent of "vibration" that you see in the 
>>>>> graphene sheet.
>>>>>  >>> Are you applying improper dihedrals appropriately? Are the 
>>>>> rings puckering, or
>>>>>  >>> just vibrating slightly (which could be normal)?
>>>>>  >>>
>>>>>  >>>> I will try position restraining the graphene structure and 
>>>>> see if that
>>>>>  >>>> resolves my problem.
>>>>>  >>>>
>>>>>  >>>> I have been able to view a trajectory for simulations of 
>>>>> fewer than
>>>>>  >>>> 20,000 time steps and see the frozen graphene lattice and the 
>>>>> NH3
>>>>>  >>>> molecules floating through space.
>>>>>  >>>>
>>>>>  >>>> Note that the segmentation fault only occurs sometime between 
>>>>> 20,000 and
>>>>>  >>>> 30,000 time steps. Could it be that the "funky" behaviour 
>>>>> associated
>>>>>  >>>> with freezing would take 20,000+ time steps to cause a 
>>>>> segmentation
>>>>>  >>>> fault?
>>>>>  >>>>
>>>>>  >>> Depends entirely upon what's causing the explosion. Any clues 
>>>>> from the trajectory?
>>>>>  >>>
>>>>>  >>> -Justin
>>>>>  >>>
>>>>>  >>>> Thanks again for your help.
>>>>>  >>>>
>>>>>  >>>> Darrell
> _______________________________________________
> gmx-users mailing list    gmx-users at 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-request at gromacs.org.
> Can't post? Read http://www.gromacs.org/mailing_lists/users.php
> 

-- 
========================================

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

========================================