[gmx-developers] Verlet Algorithm

[Elena Della valle]

I used the gromacs 4.6.5
Best regards 
Elena della Valle


> Il giorno 04 ott 2016, alle ore 16:49, David van der Spoel <spoel at xray.bmc.uu.se> ha scritto:
> 
>> On 04/10/16 15:15, Elena della Valle wrote:
>> 
>> Hi to all,
>> 
>> I'm Elena della Valle, a Ph.D. student coming from la Sapienza
>> University of Rome.
>> 
>> I'm writing because I have some questions about some modifications that
>> I have done to the update of velocities and positions in the verlet
>> algorithm. My aim is to implement the magnetic field in gromacs by
>> introducing the therm of the larmor frequency to the velocities and
>> positions. By literature I read that this has been done by the update of
>> the verlet velocities and positions with the frequency larmor therm.
>> These are  the equations in order to be more clear on what i wanted to do:
>> 
>> I did this by modifying in the update.c file in gromacs the
>> update_do_vv_vel and update_do_vv_pos ad follows:
>> 
>> if ((ptype[n] != eptVSite) && (ptype[n] != eptShell) && !nFreeze[gf][d])
>>            {
>>                v[n][0] = mv1*(mv1*v[n][0] +
>> 0.5*(w_dt*mv2*f[n][0]))+0.5*accel[ga][0]*dt +
>> w_dt*charge[n]*campoB*v[n][1]*mv1*mv1 +
>> 0.25*dt*invmass[n]*charge[n]*campoB*((w_dt*mv2*mv1*f[n][1]) +
>> accel[ga][1]*dt -2*w_dt*charge[n]*campoB*mv1*mv1*v[n][0]);
>>                v[n][1] = mv1*(mv1*v[n][1] +
>> 0.5*(w_dt*mv2*f[n][1]))+0.5*accel[ga][1]*dt +
>> w_dt*charge[n]*campoB*v[n][0]*mv1*mv1 -
>> 0.25*dt*invmass[n]*charge[n]*campoB*((w_dt*mv2*mv1*f[n][0]) +
>> accel[ga][0]*dt +2*w_dt*charge[n]*campoB*mv1*mv1*v[n][1]);
>>                v[n][2] = mv1*(mv1*v[n][2] +
>> 0.5*(w_dt*mv2*f[n][2]))+0.5*accel[ga][2]*dt;
>> 
>> 
>> 
>>            //printf("frame %lf \n", mv1);
>>           // printf("frame %d: %lf\t%lf%lf\n", n, v[n][0], v[n][1],
>> v[n][2]);
>>            }
>> 
>> in the do_update_vv_pos:
>> 
>> if ((ptype[n] != eptVSite) && (ptype[n] != eptShell) && !nFreeze[gf][d])
>>            {
>>                xprime[n][0]   = mr1*(mr1*x[n][0]+mr2*dt*v[n][0]) +
>> mr1*0.5*dt*(w_dt*mr2*f[n][0] + w_dt*charge[n]*campoB*v[n][1]*mr1);
>>                xprime[n][1]   = mr1*(mr1*x[n][1]+mr2*dt*v[n][1]) +
>> mr1*0.5*dt*(w_dt*mr2*f[n][1] - w_dt*charge[n]*campoB*v[n][0]*mr1);
>>                xprime[n][2]   = mr1*(mr1*x[n][2]+mr2*dt*v[n][2]) +
>> mr1*0.5*dt*(w_dt*mr2*f[n][2]);
>>            }
>> 
>> 
>> After that in the grompp file I modified the integrator in md-vv and I
>> run the simulation and it seems to use the position and velocities
>> modifications.
>> 
>> I want know if you think that this way is correct, if I implemented in
>> the right way
>> 
>> And Also I would like to know if in the grompp file I have to modify
>> other parameters (I use a beredensed coupling)
>> 
>> Thanks in advance
>> Sorry for bothering you
>> Best Regards
>> Elena della Valle
>> 
>> --
>> Elena della Valle
>> Ph.D. Student in Electronic Engineering
>> 
>> Department of Information Engineering, Electronics and Telecommunications
>> Sapienza, University of Rome
>> via Eudossiana, 18 00184 Rome
>> 
>> 
>> 
> Just a quick question: which version did you use? All development should be in the master version and this looks like uses something older.
> 
> 
> 
> -- 
> David van der Spoel, Ph.D., Professor of Biology
> Dept. of Cell & Molec. Biol., Uppsala University.
> Box 596, 75124 Uppsala, Sweden. Phone:    +46184714205.
> spoel at xray.bmc.uu.se    http://folding.bmc.uu.se
> -- 
> Gromacs Developers mailing list
> 
> * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-developers_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-developers or send a mail to gmx-developers-request at gromacs.org.