[gmx-users] Simulation of polarizable Carbon nanotubes

Fri Jun 26 00:21:55 CEST 2015

Hi Justin,
Thanks a lot for your reply and suggestion. I have pulled the master branch and performed "git checkout drude". I think, the file named md.cpp in /src/program/mdrun is the main file to run the code. Please let me know if I am wrong.
I am just wondering how can I execute the code in my current version of gromacs (4.6.1), which is installed in my computer and cluster. Do I need to re-install the the clone? I've tried to explore the information from Git Tutorial in Gromacs website, but could not understand much. Please accept my apology.
Regarding anisotropic polarization: I see that [water_polarization] directive was used in the topology file in the past to model anisotropic polarization of water molecule, while [polarization] was used for isotropic polarization. I am assuming by mentioning "Build the manual for the syntax of the new topology directives" in your previous email, you meant to create a new directive, for instance [cnt_polarization], in my topology file. I am not sure how to do that. I do not see any information regarding polarization directives in the gromacs manual as well. I would highly appreciate if you please provide some information how can I proceed with this.
Thanks once again. I always appreciate your help.

Regards,Golam

----------------------------

Message: 3
Date: Wed, 24 Jun 2015 20:31:10 -0400
From: Justin Lemkul <jalemkul at vt.edu>
To: Discussion list for GROMACS users <gmx-users at gromacs.org>
Subject: Re: [gmx-users] Simulation of polarizable Carbon nanotubes
Message-ID: <558B4BCE.9070801 at vt.edu>
Content-Type: text/plain; charset=utf-8; format=flowed

(Please don't CC me separately, just post to the list)

On 6/24/15 7:32 PM, S.M. Golam Mortuza wrote:
> Dear Gromacs Users,
>
> I am trying to run a simulation using polarizable force field for carbon
> nanotube (CNT), which has 304 carbon atoms. As a test purpose, I am trying to
> simulate a system containing a pair of KCl ion in an infinite carbon nanotube.
> However, my ultimately goal is to simulate a system comprising water, KCl inside
> the CNT.
>
> I have been trying to mimic the paper
> /http://www.sciencedirect.com/science/article/pii/S0927025607000456/
>
> Also, I have been following the thread on this topic:
> https://www.mail-archive.com/gromacs.org_gmx-users@maillist.sys.kth.se/msg10520.html
>
> Again, as a test purpose, I have put only one dummy atom (atom id 306), which is
> harmonically connected with a shell atom (atom id 305), inside a hexagonal ring
> (atom id: 1 20 38 39 57 58). I have used position restrain for the nanotube and
> the dummy atom. I think I've followed every procedure that Justin and Peter
> suggested in the thread.  Some parts of my itp file of the CNT looks like:
>
> /[ atoms ]
> ;  nr      type  resnr residue  atom  cgnr    charge      mass  typeB
> chargeB      massB
> 1    CJ    1    CNT    C    1    0.25    12.011    ;    qtot    0/
> /2    CJ    1    CNT    C    2    0    12.011    ;    qtot    0
> 3    CJ    1    CNT    C    3    0    12.011    ;    qtot    0/
> /...................................................................................../
> /.....................................................................................
> /
> /305    SL    1    CNT    SL    1    -1.5    1    ;    qtot    0 ; Shell
> 306    VS    1    CNT    VS    300    0    0    ;    qtot    0 ; Dummy/
> /
> /
> /[ polarization ]
> ; virtual_site    shell    functiontype    alpha (in nm^3)
>      306        305    1        0.1296
> /
>
>
> //
> /[bonds]/
> /......../
>
> //
> /[angles]/
> /....../
>
> //
> /[dihedrals]/
> /....../
> /
> /
> /[ virtual_sitesn ]
> ; Site funct from
> 306 1 1 20 38 39 57 58/
>
> /[ exclusions ]
> ; iatom excluded from interaction with i
> 305 1    2    3    4    5    6    7    8    9  ..... 304/  ; I have used
> exclusions from the shell atom with all carbon atoms of the nanotube. I am not
> sure whether I should exclude more exclusions.
>
> Sometimes, the simulation gives error as:
>
> /EM did not converge in 20 iterations/
>
> /Segmentation fault (core dumped)/
>
>
> I tried to play with different parameters, such as polarizablity, exclusion
> list, etc. However, sometimes, I can see following messages in my window, and
> the simulation does not move further.
>
> /starting mdrun 'KCl electrolyte in CNT'
> 100000 steps,    100.0 ps.
> MDStep=    0/ 0 EPot: -4.83364532e+02, rmsF: 2.17e+01
>
> step 0MDStep=    1/ 0 EPot: -4.34285614e+02, rmsF: 2.19e+01
> MDStep=    2/ 0 EPot: -3.03536743e+02, rmsF: 2.22e+01
> MDStep=    3/ 0 EPot: -1.36375397e+02, rmsF: 2.24e+01
> MDStep=    4/ 0 EPot: 1.27241879e+01, rmsF: 2.25e+01
>
> NOTE: Turning on dynamic load balancing
>
> MDStep=    5/ 0 EPot: 1.01205765e+02, rmsF: 2.27e+01
> MDStep=    6/ 0 EPot: 1.15893883e+02, rmsF: 2.28e+01
> MDStep=    7/ 0 EPot: 7.45831146e+01, rmsF: 2.28e+01
> MDStep=    8/ 0 EPot: 1.53622913e+01, rmsF: 2.29e+01
> MDStep=    9/ 0 EPot: -2.12497330e+01, rmsF: 2.29e+01
> MDStep=  10/ 0 EPot: -1.05190954e+01, rmsF: 2.29e+01
> MDStep=  11/ 0 EPot: 4.68694420e+01, rmsF: 2.29e+01
> MDStep=  12/ 0 EPot: 1.29762100e+02, rmsF: 2.29e+01
> MDStep=  13/ 0 EPot: 2.02292053e+02, rmsF: 2.29e+01
> MDStep=  14/ 0 EPot: 2.34559998e+02, rmsF: 2.29e+01
> MDStep=  15/ 0 EPot: 2.15615555e+02, rmsF: 2.29e+01
> MDStep=  16/ 0 EPot: 1.57610504e+02, rmsF: 2.30e+01
> MDStep=  17/ 0 EPot: 8.87061996e+01, rmsF: 2.30e+01
> MDStep=  18/ 0 EPot: 4.16818733e+01, rmsF: 2.30e+01
> MDStep=  19/ 0 EPot: 3.34845123e+01, rmsF: 2.31e+01
> MDStep=  20/ 0 EPot: 6.07384720e+01, rmsF: 2.31e+01
> MDStep=  21/ 0 EPot:          inf, rmsF: 2.30e+01/
>
> I have checked that my shell and dummy atoms remain close to each other.
> However, sometimes the CNT blows up. I would really appreciate if anyone can
> help me with this.
>
> Justin,
> I have seen that your paper on polarizable force field has been published
> already. I am using gromacs version 4.6.1. Can I use the codes for polarizable
> force field in this version. Any sort of help would be highly appreciated.
>
> Please let me know, if you need any further information. Thanks in advance.
>

The only benefit from our code is the anisotropic polarization function.  You 
have defined an isotropic polarization in this model, which is not what the 
paper specifies.  Pull the latest master branch and do a 'git checkout drude' to 
get the code.  Build the manual for the syntax of the new topology directives.

Note that the extended Lagrangian formalism is incompatible with this model 
because the "parent atom" of the shell/Drude is actually a virtual site.  The 
Lagrangian will collapse spectacularly.  Stick with SCF.

De-complicate the system with some model system like benzene or something. 
Don't put in ions or anything yet.  The attachment of a shell to a dummy is a 
bit strange and may need to be considered differently in the code, but that's 
not something I can say with any degree of certainty.

-Justin

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

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.

End of gromacs.org_gmx-users Digest, Vol 134, Issue 122
*******************************************************