package name.mjw.jamber;
import javax.vecmath.Vector3d;
/**
* Angle object representing a harmonic interaction between three atoms
*
* @author mjw
*
*/
public class Angle extends ForceFieldTerm {
/**
* Atom i in angle term.
*/
private Atom i;
/**
* Atom j in angle term.
*/
private Atom j;
/**
* Atom k in angle term.
*/
private Atom k;
/**
* Force constant in kcal/(mol rad**2)
*/
private double forceConstant;
/**
* Equilibrium angle parameter in degrees.
*/
private double equilibriumAngle;
/**
* Harmonic angle
*
* @param i
* Atom i in angle term
* @param j
* Atom j in angle term
* @param k
* Atom k in angle term
* @param equilibriumAngle
* Equilibrium angle in degrees
* @param forceConstant
* Force constant in kcal/(mol rad**2)
*/
public Angle(Atom i, Atom j, Atom k, double equilibriumAngle,
double forceConstant) {
this.setAtomI(i);
this.setAtomJ(j);
this.setAtomK(k);
this.setEquilibriumAngle(equilibriumAngle);
this.setForceConstant(forceConstant);
}
/**
* Set atom i of the angle term
*/
private void setAtomI(Atom i) {
this.i = i;
}
/**
* Gets atom i of the angle term
*/
public Atom getI() {
return i;
}
/**
* Set atom i of the angle term
*/
private void setAtomJ(Atom j) {
this.j = j;
}
/**
* Gets atom j of the angle term
*/
public Atom getJ() {
/**
* Gets atom j of the angle term
*/
return j;
}
/**
* Sets atom k in the angle.
*/
private void setAtomK(Atom k) {
this.k = k;
}
/**
* Gets atom k in the angle.
*/
public Atom getK() {
return k;
}
/**
* Set the angle force constant parameter, in kcal/mol/(rad**2).
*/
private void setForceConstant(double forceConstant) {
this.forceConstant = forceConstant;
}
/**
* Return the force constant parameter, in kcal/mol/(rad**2).
*/
public double getForceConstant() {
return forceConstant;
}
/**
* Set the equilibrium angle parameter value, in degrees.
*/
private void setEquilibriumAngle(double equilibriumAngle) {
this.equilibriumAngle = equilibriumAngle;
}
/**
* Return the equilibrium angle parameter value, in degrees.
*/
public double getEquilibriumAngle() {
return equilibriumAngle;
}
/**
* Returns the current potential energy of the angle in kcal/mol.
*/
@Override
public double getPotentialEnergy() {
return forceConstant
* square(Math.toRadians(getCurrentAngle() - equilibriumAngle));
}
/**
* Returns the current (scalar) gradient of the potential energy term. This
* is the analytical derivative of the potential energy at this point
*/
@Override
public double getAnalyticalGradient() {
return 2 * forceConstant
* Math.toRadians(getCurrentAngle() - equilibriumAngle);
}
/**
* Returns the current angle in degrees
*
* See <a
* href="http://cbio.bmt.tue.nl/pumma/index.php/Theory/Potentials">eq.
* (4)</a>.
*
*/
public double getCurrentAngle() {
Vector3d Rij = new Vector3d(i.getPosition());
Vector3d Rj = new Vector3d(j.getPosition());
Vector3d Rkj = new Vector3d(k.getPosition());
// Both bond vectors must be defined towards
// the atom that they share (here atom j).
Rij.sub(Rj);
Rkj.sub(Rj);
// Normalise both vectors.
Rij.normalize();
Rkj.normalize();
// Dot product; angle between two vectors
return Math.toDegrees(Math.acos(Rij.dot(Rkj)));
}
/**
* Calculates and imparts the force of the angle to its member atoms
*/
@Override
public void evaluateForce() {
Vector3d temp3d = new Vector3d();
temp3d.sub(i.getPosition(), j.getPosition());
double unitForce = getAnalyticalGradient() / temp3d.length();
temp3d.scale(unitForce);
i.addForce(temp3d);
k.subForce(temp3d);
}
@Override
public String toString() {
StringBuilder result = new StringBuilder();
String NEW_LINE = System.getProperty("line.separator");
result.append(this.getClass().getName()).append(" Object {").append(NEW_LINE);
result.append(String.format(" Atoms:\t\t\t\t\t\t%-2s-%-2s-%-2s",
i.getName(),
j.getName(),
k.getName())).append(NEW_LINE);
result.append(String.format(" AtomTypes:\t\t\t\t\t%-2s-%-2s-%-2s",
i.getAMBERAtomType(),
j.getAMBERAtomType(),
k.getAMBERAtomType())).append(NEW_LINE);
result.append(" forceConstant ( kcal/mol/(rad**2) ) : \t\t").append(forceConstant).append(NEW_LINE);
result.append(" equilibriumAngle (degrees): \t\t\t").append(equilibriumAngle).append(NEW_LINE);
result.append(" current angle (degrees): \t\t\t").append(this.getCurrentAngle()).append(NEW_LINE);
result.append(" current potential contribution (kcal/mol): \t").append(this.getPotentialEnergy()).append(NEW_LINE);
result.append("}");
return result.toString();
}
}