Re: [ESPResSo-devel] error in dawaanr interaction

[Rudolf Weeber]

Hi Markus,
On Fri, Sep 26, 2014 at 03:59:49PM +0300, Markus Gusenbauer wrote:
> thanks for your fast answer. When i searched for dprefac*, I
> couldn't find it in the computation of dawaanr.

>Just with ROTATION,
> which I don't use. Now I am curious again, where does l_b act in the
> code in case of dawaanr? Also when I search for bjerrum, it is also
> not acting for dawaanr interactions.
src/core/magnetic_non_p3m_methods.cpp:
[...]
// Calculates dipolar energy and/or force between two particles
double calc_dipole_dipole_ia(Particle* p1, Particle *p2, int force_flag)
{
[...]
  u= coulomb.Dprefactor* ( pe1/r3 -   pe4*pe2*pe3);

[...]
  // Force, if requested
  if(force_flag) { 
[...]
    // Add the force to the particles
    p1->f.f[0] +=coulomb.Dprefactor*ffx;
    p1->f.f[1] +=coulomb.Dprefactor*ffy;
    p1->f.f[2] +=coulomb.Dprefactor*ffz;
    p2->f.f[0] -=coulomb.Dprefactor*ffx;
    p2->f.f[1] -=coulomb.Dprefactor*ffy;
    p2->f.f[2] -=coulomb.Dprefactor*ffz;

[...]
      // Torques
#ifdef ROTATION
[...]
    p1->f.torque[0]+=coulomb.Dprefactor*(-ax/r3+bx*cc);
    p1->f.torque[1]+=coulomb.Dprefactor *(-ay/r3+by*cc);
    p1->f.torque[2]+=coulomb.Dprefactor*(-az/r3+bz*cc);
[...] 
    p2->f.torque[0]+=coulomb.Dprefactor* (ax/r3+bx*d);
    p2->f.torque[1]+=coulomb.Dprefactor*(ay/r3+by*d);
    p2->f.torque[2]+=coulomb.Dprefactor* (az/r3+bz*d);
#endif
[...]
}    
As far as I can see, the pre-factor is there, everywhere.


The Bjerum length is set here:
core/interaction_data.cpp:
[...]
void recalc_coulomb_prefactor()
{
[...]
#ifdef DIPOLES
  if(temperature > 0.0)
    coulomb.Dprefactor = coulomb.Dbjerrum * temperature; 
  else
    coulomb.Dprefactor = coulomb.Dbjerrum;
#endif
}

Hope that helps.
Regards, Rudolf