[gmx-users] Simulation of polarizable Carbon nanotubes

Peter Kroon p.c.kroon at rug.nl
Mon Feb 23 09:47:33 CET 2015


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.
>
> [ 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)
>
> 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.
>
> 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.
>>




More information about the gromacs.org_gmx-users mailing list