ADiheRestraint Class Reference

#include <FreeEnergyRestrain.h>

Inheritance diagram for ADiheRestraint:

Public Member Functions
	ADiheRestraint ()
void	PrintInfo ()
virtual Bool_t	TwoTargets ()=0
virtual double	GetDiheTarget1 ()=0
virtual double	GetDiheTarget2 ()=0
Public Member Functions inherited from ARestraint
	ARestraint ()
virtual	~ARestraint ()
int	GetNumGroups ()
void	SetKf (double Kf)
double	GetKf ()
void	SetLambdaKf (double LambdaKf)
void	SetLambdaRef (double LambdaRef)
double	GetLambdaKf ()
double	GetLambdaRef ()
void	SetGroup (AGroup &Group, int GroupIndex)
void	SetGroups (AGroup &Group1)
void	SetGroups (AGroup &Group1, AGroup &Group2)
void	SetGroups (AGroup &Group1, AGroup &Group2, AGroup &Group3)
void	SetGroups (AGroup &Group1, AGroup &Group2, AGroup &Group3, AGroup &Group4)
void	UpdateCOMs (GlobalMasterFreeEnergy &CFE)
void	DistributeForce (int WhichGroup, AVector Force, GlobalMasterFreeEnergy &CFE)
virtual AVector	GetGradient (int WhichGroup)=0
virtual double	GetEnergy ()=0
virtual void	GetStr (char *Str)=0
virtual Bool_t	IsForcing ()
virtual double	Get_dU_dLambda ()
virtual void	SetRefPos (AVector)
virtual void	SetRefDist (double)
virtual void	SetRefAngle (double)
virtual void	SetBound (Bound_t)
virtual void	SetLowerAngle (double)
virtual void	SetUpperAngle (double)
virtual void	SetIntervalAngle (double)
virtual void	SetStartPos (AVector)
virtual void	SetStopPos (AVector)
virtual void	SetStartDist (double)
virtual void	SetStopDist (double)
virtual void	SetStartAngle (double)
virtual void	SetStopAngle (double)

Protected Member Functions
double	GetE (double RefDihe, double Const)
AVector	GetGrad (int WhichGroup, double RefDihe, double Const)
AVector	gradU (AVector &P1P2P3, AVector &P4P5P6, AVector &dP1, AVector &dP2, AVector &dP3, AVector &dP4, AVector &dP5, AVector &dP6)
Protected Member Functions inherited from ARestraint
double	GetAngle (AVector &A, AVector &B, AVector &C)
double	GetDihe (AVector &A, AVector &B, AVector &C, AVector &D)
void	EarlyExit (const char *Str, int AtomID)

Additional Inherited Members
Protected Attributes inherited from ARestraint
double	m_Kf
int	m_NumGroups
AGroup *	m_pGroups
AVector *	m_pCOMs
Static Protected Attributes inherited from ARestraint
static double	m_LambdaKf = 1.0
static double	m_LambdaRef = 0.0

Detailed Description

Definition at line 146 of file FreeEnergyRestrain.h.

Constructor & Destructor Documentation

ADiheRestraint()

ADiheRestraint::ADiheRestraint

(

)

Definition at line 537 of file FreeEnergyRestrain.C.

References ARestraint::m_NumGroups, ARestraint::m_pCOMs, and ARestraint::m_pGroups.

                               {
//----------------------------------------------------------------------------
// each ADiheRestraint restrains the dihedral angle between 4 groups of atoms
//----------------------------------------------------------------------------
  m_NumGroups = 4;
  m_pGroups = new AGroup[m_NumGroups];
  m_pCOMs = new AVector[m_NumGroups];
}

Member Function Documentation

GetDiheTarget1()

virtual double ADiheRestraint::GetDiheTarget1 ( )

pure virtual

Implemented in AForcingDiheRestraint, ABoundDiheRestraint, and AFixedDiheRestraint.

Referenced by PrintInfo().

GetDiheTarget2()

virtual double ADiheRestraint::GetDiheTarget2 ( )

pure virtual

Implemented in AForcingDiheRestraint, ABoundDiheRestraint, and AFixedDiheRestraint.

Referenced by PrintInfo().

GetE()

double ADiheRestraint::GetE ( double  RefDihe, double  Const  )

protected

Definition at line 587 of file FreeEnergyRestrain.C.

References ARestraint::GetDihe(), and ARestraint::m_pCOMs.

Referenced by AFixedDiheRestraint::GetEnergy(), ABoundDiheRestraint::GetEnergy(), and AForcingDiheRestraint::GetEnergy().

                                                         {
//---------------------------------------------------------------------------
// calculate and return the Energy for this angle restraint.
//
//    E = (E0/2) * (1 - cos(Chi - ChiRef))
//
// where Chi is the dihedral angle between 4 centers-of-mass of restrained atoms,
// m_pCOMs[0] -- m_pCOMs[1] -- m_pCOMs[2] -- m_pCOMs[3].
// ChiRef is the reference angle.
//
// Note:  COM's are calculated before this routine is called.
//---------------------------------------------------------------------------
  double  Angle;

  Angle = GetDihe(m_pCOMs[0], m_pCOMs[1], m_pCOMs[2], m_pCOMs[3]);
  return( (Const/2.0) * (1.0 - cos(Angle-RefAngle)) );
}

GetGrad()

AVector ADiheRestraint::GetGrad ( int  WhichGroup, double  RefDihe, double  Const  )

protected

Definition at line 606 of file FreeEnergyRestrain.C.

References ASSERT, AVector::cross(), AVector::Dist(), AVector::dot(), ARestraint::GetDihe(), gradU(), kAlmostMinusOne, kAlmostOne, ARestraint::m_pCOMs, and AVector::Set().

Referenced by AFixedDiheRestraint::GetGradient(), ABoundDiheRestraint::GetGradient(), and AForcingDiheRestraint::GetGradient().

                                                               {
//---------------------------------------------------------------------------
// calculate and return the gradient for this dihedral angle restraint.
//
//    E = (E0/2) * (1 - cos(Chi - ChiRef))
//
// return:  grad(E)
//
// Notes: COM's are calculated before this routine is called.
//---------------------------------------------------------------------------
  AVector A, B, C, D;

  ASSERT((WhichGroup==0)||(WhichGroup==1)||(WhichGroup==2)||(WhichGroup==3));

  if ((WhichGroup==0) || (WhichGroup==1)) {
    A = m_pCOMs[0];
    B = m_pCOMs[1];
    C = m_pCOMs[2];
    D = m_pCOMs[3];
  }
  // re-state the problem so the gradient is solved for either atoms 0 or 1
  else {
    A = m_pCOMs[3];
    B = m_pCOMs[2];
    C = m_pCOMs[1];
    D = m_pCOMs[0];
    if (WhichGroup==3) {WhichGroup=0;}
    if (WhichGroup==2) {WhichGroup=1;}
  }

  AVector CD(D - C);
  AVector CB(B - C);
  AVector BC(C - B);
  AVector BA(A - B);
  AVector AC(C - A);
  AVector CDxCB, BCxBA;
  AVector Vec;
  double  phi;
  double  top, bot, u;

  CDxCB = CD.cross(CB);
  BCxBA = BC.cross(BA);

  top = CDxCB.dot(BCxBA);
  bot = CDxCB.Dist() * BCxBA.Dist();

  u = top/bot;
  // protect against acos(<-1.0), acos(>1.0), sqrt(<0), and divide by 0
  if (u < kAlmostMinusOne) {u = kAlmostMinusOne;}
  if (u > kAlmostOne)      {u = kAlmostOne;}

  // get dihedral using atan
  phi = GetDihe(A,B,C,D);

  AVector dP1, dP2, dP3;
  AVector dP4, dP5, dP6;
  ASSERT((WhichGroup==0) || (WhichGroup==1));
  if (WhichGroup==0) {
    dP1.Set( 0,      0,      0    );
    dP2.Set( 0,      0,      0    );
    dP3.Set( 0,      0,      0    );
    dP4.Set( 0,     -BC[2],  BC[1]);
    dP5.Set( BC[2],  0,     -BC[0]);
    dP6.Set(-BC[1],  BC[0],  0    );
  }
  else {
    dP1.Set( 0,     -CD[2],  CD[1]);
    dP2.Set( CD[2],  0,     -CD[0]);
    dP3.Set(-CD[1],  CD[0],  0    );
    dP4.Set( 0,      AC[2], -AC[1]);
    dP5.Set(-AC[2],  0,      AC[0]);
    dP6.Set( AC[1], -AC[0],  0    );
  }

  Vec = gradU(CDxCB, BCxBA, dP1, dP2, dP3, dP4, dP5, dP6);
  Vec *= (Const/2.0) * sin(phi-RefAngle) * (-1.0/sqrt(1.0 - u*u));

  // flip gradient for negative angles
  if (phi < 0) {
    Vec *= -1.0;
  }
  
  return(Vec);
}

gradU()

AVector ADiheRestraint::gradU ( AVector &  P1P2P3, AVector &  P4P5P6, AVector &  dP1, AVector &  dP2, AVector &  dP3, AVector &  dP4, AVector &  dP5, AVector &  dP6  )

protected

Definition at line 693 of file FreeEnergyRestrain.C.

References AVector::Dist(), and AVector::dot().

Referenced by GetGrad().

                                                                              {
//----------------------------------------------------------------
// calculate the gradient for ((P1P2P3.dot.P4P5P6)/(mag(P1P2P3)*mag(P4P5P6)))
// P1P2P3 = (P1)i + (P2)j + (P3)k
// P4P5P6 = (P4)i + (P5)j + (P6)k
// dP1 = (d(P1)/dx)i + (d(P1)/dy)j +(d(P1)/dz)k
// dP2 = (d(P2)/dx)i + (d(P2)/dy)j +(d(P2)/dz)k
// dP3 = (d(P3)/dx)i + (d(P3)/dy)j +(d(P3)/dz)k
// dP4 = (d(P4)/dx)i + (d(P4)/dy)j +(d(P4)/dz)k
// dP5 = (d(P5)/dx)i + (d(P5)/dy)j +(d(P5)/dz)k
// dP6 = (d(P6)/dx)i + (d(P6)/dy)j +(d(P6)/dz)k
//----------------------------------------------------------------
  double  Mag123, Mag456, Dot;
  double  Const1, Const2, Const3;
  double  P1, P2, P3, P4, P5, P6;
  AVector RetVec;

  P1 = P1P2P3[0];  P2 = P1P2P3[1];  P3 = P1P2P3[2];
  P4 = P4P5P6[0];  P5 = P4P5P6[1];  P6 = P4P5P6[2];

  Mag123 = P1P2P3.Dist();
  Mag456 = P4P5P6.Dist();
  Dot    = P1P2P3.dot(P4P5P6);

  Const1 =         1.0 / (Mag123*Mag456);
  Const2 = -Dot * (1.0 / (Mag123*Mag456*Mag456*Mag456));
  Const3 = -Dot * (1.0 / (Mag456*Mag123*Mag123*Mag123));

  RetVec = (dP4*P1 + dP1*P4 + dP5*P2 + dP2*P5 + dP6*P3 + dP3*P6) * Const1 +
           (dP4*P4 + dP5*P5 + dP6*P6)                            * Const2 +
           (dP1*P1 + dP2*P2 + dP3*P3)                            * Const3;

  return(RetVec);
}

PrintInfo()

void ADiheRestraint::PrintInfo ( )

virtual

Implements ARestraint.

Definition at line 547 of file FreeEnergyRestrain.C.

References endi(), ARestraint::GetDihe(), GetDiheTarget1(), GetDiheTarget2(), iout, kPi, ARestraint::m_pCOMs, and TwoTargets().

                               {
//--------------------------------------------------------------------
// print the dihedral angle for this dihedral restraint
//--------------------------------------------------------------------
  double  Dihedral;
  char    Str1[20], Str2[20], Str3[20];

  Dihedral = GetDihe(m_pCOMs[0], m_pCOMs[1], m_pCOMs[2], m_pCOMs[3]) * (180/kPi);
  sprintf(Str1, "%8.3f", Dihedral);
  Dihedral = GetDiheTarget1() * (180/kPi);
  sprintf(Str2, "%8.3f", Dihedral);
  Dihedral = GetDiheTarget2() * (180/kPi);
  sprintf(Str3, "%8.3f", Dihedral);

#if defined(_VERBOSE_PMF)
  iout << "Dihedral = ";
  iout << Str1;
  iout << " degrees";
  iout << "  Target = ";
  iout << Str2;
  iout << " degrees";
  if (TwoTargets()) {
    iout << " to ";
    iout << Str3;
    iout << " degrees";
  }
  iout << std::endl << endi;
#else
  iout << Str1;
  iout << "  ";
  iout << Str2;
  if (TwoTargets()) {
    iout << ", ";
    iout << Str3;
  }
  iout << " | ";
#endif
}

TwoTargets()

virtual Bool_t ADiheRestraint::TwoTargets ( )

pure virtual

Implemented in AForcingDiheRestraint, ABoundDiheRestraint, and AFixedDiheRestraint.

Referenced by PrintInfo().

---

The documentation for this class was generated from the following files:

• FreeEnergyRestrain.h
• FreeEnergyRestrain.C