On 27/01/2011 12:11 PM, Makoto Yoneya wrote:
Dear Berk and all:

I'd tried to rewrite the routine

src/gmxlib/nonbonded/nb_generic.c
to modify the LJ potential to the shifted and truncated one.

First, I'd add the new switch(ivdw) as case 4, but when I'd tried;
[ defaults ]
; nbfunc comb-rule   gen-pairs   fudgeLJ  fudgeQQ
   4      1     no    1.0   1.0,
I'd got,
ERROR 1 [file topol.top, line 8]:
   Invalid nonbond function selector '4' using LJ.

So you have to update the machinery that parses the .top to recognise that the value of 4 is now legal.

Then, I'd rewrite the standard switch(ivdw) case 1 as follows
(with real tc12c5 and int factor).

#ifdef TRUNCATED_LJ
        tc12c6 = 2.0*c12/c6;
        factor = ceil(floor(tc12c6*rinvsix)/(tc12c6*rinvsix));
        fscal += factor*(12.0*Vvdw_rep-6.0*Vvdw_disp)*rinvsq;
        Vvdwtot = Vvdwtot+factor*(Vvdw_rep-Vvdw_disp+0.5*c6/tc12c6);
#else
        fscal += (12.0*Vvdw_rep-6.0*Vvdw_disp)*rinvsq;
        Vvdwtot = Vvdwtot+Vvdw_rep-Vvdw_disp;
#endif

However with this modification (recompile with #define TRUNCATED_LJ),
the simulation result was exactly the same.

I'd appreciate to tell me what kind of mistake I'd done?

Did you set the environment variable to actually call the generic nonbonded lists?

Mark

Makoto Yoneya, Dr.
AIST
http://staff.aist.go.jp/makoto-yoneya/

Yes.

The pow function is expensive though. The code will run much faster
if you can use rinvsix, such as check for  2*rinvsix>  c6/c12.
(I forgot the factor 2 in my previous mail).

Berk

From: makoto-yoneya at aist.go.jp
To: gmx-users at gromacs.org
Date: Tue, 11 Jan 2011 10:10:56 +0900
Subject: [gmx-users] truncated LJ potential

Dear Berk:

Thanks again for the further reply.

The LJ potential and force code in the above looks like in the c6-c12
form
not in epsilon-sigma one.
The LJ potential modification I'd like to try is based on the epsilon-
sigma form and the mixing rule is the Lorents-Bertelot's one.
Could you kindly tell me the LJ potential and force routine in the
epsilon-
Sigma form.
There is no such code.

You can simply check for rinvsix>  c6/c12
Is it means that the LJ potential in epsion-sigma form (with the
Lorents-Bertelot
mixing rule) is evaluated after the coversion into c6-c12 form in GROMACS?
Then, may I evaluate the modified LJ potential:

V(r)
    = 4*epsilon*{ (sigma/r)^(12) - (sigma/r)^6 + (1/4) } for r<=
2^(1/6)*sigma
    = 0 for r>   2^(1/6)*sigma
with translated into the equivalent c6-c12 form:

V(r)
    = (c12/r)^(12) - (c6/r)^6 + (c6/2)*(c6/2*c12) for r<= (2*c12/c6)^(1/6)
    = 0 for r>  (2*c12/c6)^(1/6)

in the following routines.
You can also set the environment variable nb_generic.c and modify
src/gmxlib/nonbonded/nb_generic.c, but might lead to somewhat
slower simulations.
If my understanding in the above would correct, I'll try that.

--
gmx-users mailing list    gmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

Reply via email to