[Top][All Lists]

[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
## [ESPResSo] BOND_ANGLE_HARMONIC

**From**: |
Mehmet Sayar |

**Subject**: |
[ESPResSo] BOND_ANGLE_HARMONIC |

**Date**: |
Wed, 09 Jul 2008 11:36:38 +0300 |

**User-agent**: |
Thunderbird 2.0.0.14 (X11/20080421) |

Hi all,
In the bond_angle_harmonic routine shouldn't there be a check for
cosine > 1.0. In other words, I think the code should be corrected as:
#ifdef BOND_ANGLE_HARMONIC
{
double phi,sinphi;
if ( cosine > TINY_COS_VALUE ) cosine = TINY_COS_VALUE;
if ( cosine < -TINY_COS_VALUE) cosine = -TINY_COS_VALUE;
phi = acos(-cosine);
sinphi = sin(phi);
if ( sinphi < TINY_SIN_VALUE ) sinphi = TINY_SIN_VALUE;
fac *= (phi - iaparams->p.angle.phi0)/sinphi;
}
#endif
If everybody agrees, I will add this to the cvs.

`P.S. I originally send this message to address@hidden But since
``that email did not show up in the email-thread, I am re-sending the message.
`Is address@hidden still in use?

`And I already got an answer from Hanjo to the email sent to address@hidden,
``which was:
`
> So far we have the checks only for functions with a singularity when
> cosine is aproaching 1 (and sinphi is aproaching 0). Can the scalar
> product yield values >1 ( <-1)? (which would yield an error in the acos
> function). If yes then add
> if ( cosine > 1 ) cosine = 1;
> if ( cosine < -1) cosine = -1;
> To the corresponding force and energy routines.
>
> Greetings Hanjo
Does everyone agree? Any objections?
Greetings from Istanbul,
Mehmet
--
________________________________________________________
Mehmet Sayar - Assist. Prof. Koc University
Email:address@hidden Dept. Mechanical Eng.
Phone:+90-212-338-1840 Rumelifeneri Yolu, Sariyer
Fax: +90-212-338-1548 34450 Istanbul, Turkey

**[ESPResSo] BOND_ANGLE_HARMONIC**,
*Mehmet Sayar* **<=**