[gmx-users] Simulation of polarizable Carbon nanotubes
Justin Lemkul
jalemkul at vt.edu
Mon Feb 23 13:39:55 CET 2015
On 2/23/15 3:47 AM, Peter Kroon wrote:
>
> On 02/22/2015 05:33 PM, Jashimuddin Ashraf wrote:
>> Thanks for your reply Dr. Lemkul
>>
>> I changed the bond length to zero in my topology file, forcefield file and
>> in the n2t file. I also added the [ exclusion] section-
> If memory serves, you don't need a bond between a shell particle and
> it's "parent"; the polarization directive takes care of it's position.
> If you want, you can define a non-interacting connection (bond type 5)
> between them.
This is correct. The [polarization] directive is just another bonded
interaction. The value of the force constant is back-calculated from the atomic
polarizability, alpha.
>>
>> [ exclusions ]
>> ; iatom excluded from interaction with i
>> 1 2 3 4 5 6 7 8
>> 2 1 3 4 5 6 7 8
>> 3 1 2 4 5 6 7 8
>> 4 1 2 3 5 6 7 8
>> 5 1 2 3 4 6 7 8
>> 6 1 2 3 4 5 7 8
>> 7 1 2 3 4 5 6
>> 8 1 2 3 4 5 6
>>
>>
>> in my topology. But the same error keeps showing up. What did I possibly do
>> wrong?
> AFAIK this will exclude all atoms from eachother. Just the last line
> should be enough if you just want to exclude the shell-carbon interactions
>>
>> In your mail, you were kind to reply with-
>>
>> "I only looked at the paper briefly, but it seems they are working with a
>> model that makes use of anisotropic polarization. In GROMACS, this is
>> currently only available for water, so the model would not be supported."
>>
>> Are you suggesting that, this model is something we should not work with in
>> GROMACS right now? (I am really sorry, I could not understand this part
>> properly)
From my quick read of the paper (focusing only a narrow section of the
methods), the authors list alpha_zz only, implying that the shell is only
polarizable along the normal to the ring. Therefore the other elements of the
polarization tensor are zero.
>>
>> Also, it would be very nice, if you could give us an idea about the time
>> required to publish the codes and adding them in the new versions of
>> GROMACS.
>>
That depends on reviewers :)
The paper should be submitted this week, after which I will be working to submit
the code to the master development branch. The code will not be included in 5.1
(I missed the deadline, due to some debugging that took longer than expected),
but the code will be available as a patch to the master branch in the coming weeks.
-Justin
>> Thanks in advance,
> Lastly, for defining the virtual site you could look at the
> virtual_sitesn directives and define a center of geometry or a center of
> mass. Im not sure if it has any advantages over a virtual site as it is
> defined now though.
> If your simulation blows up, try analysing the distance between the
> shell and it's parent; it should give you a hint about what's going wrong.
>
> Good luck!
>
> Regards,
> Peter Kroon
>>
>> On Sat, Feb 21, 2015 at 12:20 AM, Justin Lemkul <jalemkul at vt.edu> wrote:
>>
>>>
>>> On 2/20/15 12:57 PM, Jashimuddin Ashraf wrote:
>>>
>>>> Dear users,
>>>>
>>>> I want to perform a molecular dynamic simulation of polarizable carbon
>>>> nanotubes. I intend to implement this paper-
>>>>
>>>> http://www.sciencedirect.com/science/article/pii/S0927025607000456
>>>>
>>>> I digged up the manual but could not find much help from it. I went
>>>> through
>>>> some mails in the gromacs user maillist, studied some .itp files and
>>>> learned some elementary stuffs regarding the addition of virtual sites and
>>>> shell atoms.I understand that I have to add both shell and virtual sites
>>>> in
>>>> this case.
>>>>
>>>> Now, before jumping right into a big nanotube molecule, I was trying to
>>>> perform a simulation with a single benzene ring with a virtual site placed
>>>> at the center and a shell attached to the virtual site.
>>>>
>>>> In my forcefield.itp file have the virtual site and the shell declared
>>>> like
>>>> this-
>>>>
>>>> ------------------------------------------------------------
>>>> -------------------------------------------------------------------------
>>>> [ atomtypes ]
>>>> ; name mass charge ptype sigma eps
>>>> CJ1 1 12.01100 0.25 A 3.40000e-01 3.60100e-01
>>>> VS 1 0 0 D 0 0
>>>> SL 1 0 -1.5 S 0 0
>>>>
>>>> [ bondtypes ]
>>>> `; i j func b0 kb gamma
>>>> CJ1 CJ1 3 0.1418 478.9000 21.867
>>>> VS SL 1 0.06 2409
>>>> ------------------------------------------------------------
>>>> -------------------------------------------------------------------------
>>>>
>>>> Here, I have considered a bond between the virtual site and the shell (the
>>>> paper mentions something like it but does not provide with the bond
>>>> length). Is it a mistake?
>>>>
>>>>
>>> The equilibrium bond length should be zero (i.e. no induced polarization).
>>>
>>> and inside my topol.top file, I have-
>>>>
>>>> ------------------------------------------------------------
>>>> -------------------------------------------------------------------------
>>>> [ atoms ]
>>>> ; nr type resnr residue atom cgnr charge mass
>>>> typeB chargeB massB
>>>> 1 CJ1 1 C C 1 0.25 12.011 ;
>>>> qtot 0.25
>>>> 2 CJ1 1 C C 1 0.25 12.011 ;
>>>> qtot 0.5
>>>> 3 CJ1 1 C C 1 0.25 12.011 ;
>>>> qtot 0.75
>>>> 4 CJ1 1 C C 1 0.25 12.011 ;
>>>> qtot 1.00
>>>> 5 CJ1 1 C C 1 0.25 12.011 ;
>>>> qtot 1.25
>>>> 6 CJ1 1 C C 1 0.25 12.011 ;
>>>> qtot 1.50
>>>> 7 VS 1 C VS 1 0 0 ;
>>>> qtot 1.50
>>>> 8 SL 1 C S 1 -1.5 0 ;
>>>> qtot 0
>>>>
>>>>
>>>> [ polarization ]
>>>> ; virtual_site shell functiontype alpha (in nm^3)
>>>> 7 8 1 0.1
>>>>
>>>>
>>> I only looked at the paper briefly, but it seems they are working with a
>>> model that makes use of anisotropic polarization. In GROMACS, this is
>>> currently only available for water, so the model would not be supported.
>>> It will be soon (I know I've been saying that for a while, but our paper
>>> regarding Drude simulations in GROMACS is just about done, after which I
>>> can provide the code).
>>>
>>>
>>>> [ bonds ]
>>>> ; ai aj funct c0 c1 c2 c3
>>>> 1 2 3
>>>> 1 5 3
>>>> 2 4 3
>>>> 3 4 3
>>>> 3 6 3
>>>> 5 6 3
>>>> 7 8 1
>>>>
>>>> [ dihedrals ]
>>>> ; ai aj ak al funct c0 c1
>>>> c2 c3 c4 c5
>>>> 5 1 2 4 3
>>>> 2 1 5 6 3
>>>> 1 2 4 3 3
>>>> 6 3 4 2 3
>>>> 4 3 6 5 3
>>>> 1 5 6 3 3
>>>>
>>>> [ virtual_sites3 ]
>>>> ; detailed calculation not shown here
>>>> ; Dummy from funct a b
>>>> 7 4 5 6 1 0.5 0
>>>>
>>>>
>>> You're missing a critical element here; the paper says that the shell does
>>> not interact with the carbon atoms of the ring, so you need to define
>>> exclusions manually.
>>>
>>> -Justin
>>>
>>>
>>> ------------------------------------------------------------
>>>> -------------------------------------------------------------------------
>>>>
>>>> I ran an energy minimization for an emtol of 100 but it gives me a result
>>>> like this-
>>>>
>>>> Steepest Descents converged to machine precision in 141 steps,
>>>> but did not reach the requested Fmax < 100.
>>>> Potential Energy = -6.1810545e+06
>>>> Maximum force = 9.1656689e+12 on atom 4
>>>> Norm of force = 4.5828344e+12
>>>>
>>>> ------------------------------------------------------------
>>>> -------------------------------------------------------------------------
>>>> If I run a production MD with this, the simulation blows up with this-
>>>> MDStep= 40/18 EPot: nan, rmsF: nan
>>>> Warning: Only triclinic boxes with the first vector parallel to the x-axis
>>>> and the second vector in the xy-plane are supported.
>>>> Box (3x3):
>>>> Box[ 0]={ nan, nan, nan}
>>>> Box[ 1]={ nan, nan, nan}
>>>> Box[ 2]={ nan, nan, nan}
>>>> Can not fix pbc.
>>>> MDStep= 40/19 EPot: nan, rmsF: nan
>>>> step 40: EM did not converge in 20 iterations, RMS force nan
>>>>
>>>> ------------------------------------------------------------
>>>> -------------------------------------------------------------------------
>>>> Is something wrong with my system itself? or is there anything wrong with
>>>> my methods?
>>>>
>>>> I have been stuck at this for a very long time now and anything- any
>>>> comment or hint would be very much helpful for me.
>>>>
>>>> Thanks in advance,
>>>> Jashimuddin Ashraf
>>>>
>>>>
>>> --
>>> ==================================================
>>>
>>> Justin A. Lemkul, Ph.D.
>>> Ruth L. Kirschstein NRSA Postdoctoral Fellow
>>>
>>> Department of Pharmaceutical Sciences
>>> School of Pharmacy
>>> Health Sciences Facility II, Room 629
>>> University of Maryland, Baltimore
>>> 20 Penn St.
>>> Baltimore, MD 21201
>>>
>>> jalemkul at outerbanks.umaryland.edu | (410) 706-7441
>>> http://mackerell.umaryland.edu/~jalemkul
>>>
>>> ==================================================
>>> --
>>> Gromacs Users mailing list
>>>
>>> * Please search the archive at http://www.gromacs.org/
>>> Support/Mailing_Lists/GMX-Users_List before posting!
>>>
>>> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>>>
>>> * For (un)subscribe requests visit
>>> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
>>> send a mail to gmx-users-request at gromacs.org.
>>>
>
>
>
>
--
==================================================
Justin A. Lemkul, Ph.D.
Ruth L. Kirschstein NRSA Postdoctoral Fellow
Department of Pharmaceutical Sciences
School of Pharmacy
Health Sciences Facility II, Room 629
University of Maryland, Baltimore
20 Penn St.
Baltimore, MD 21201
jalemkul at outerbanks.umaryland.edu | (410) 706-7441
http://mackerell.umaryland.edu/~jalemkul
==================================================
More information about the gromacs.org_gmx-users
mailing list