[gmx-users] Simulation of polarizable Carbon nanotubes
Justin Lemkul
jalemkul at vt.edu
Fri Feb 20 19:20:59 CET 2015
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
==================================================
More information about the gromacs.org_gmx-users
mailing list