CfbRebridgeBase.h

#ifndef MCMD_CFB_REBRIDGE_BASE_H
#define MCMD_CFB_REBRIDGE_BASE_H

/*
* Simpatico - Simulation Package for Polymeric and Molecular Liquids
*
* Copyright 2010 - 2017, The Regents of the University of Minnesota
* Distributed under the terms of the GNU General Public License.
*/

#include "CfbEndBase.h"        // base class
#include <util/space/Vector.h>  // math utility class

namespace McMd
{

   using namespace Util;

   class McSystem;
   class Atom;

   class CfbRebridgeBase : public CfbEndBase 
   {
   
   public:
   
      CfbRebridgeBase(McSystem& system);
   
      virtual ~CfbRebridgeBase();

      virtual void readParameters(std::istream& in);
   
      virtual void loadParameters(Serializable::IArchive &ar);
   
      virtual void save(Serializable::OArchive &ar);

      virtual void setup();
   
   protected:

      double length21_;

      double length10_;

      double kappa10_;

      void deleteMiddleAtom(Atom* partPtr, Atom* prevPtr, Atom* nextPtr,
             int prevBType, int nextBType, double &rosenbluth, double &energy);

      void addMiddleAtom(Atom* partPtr, Atom* prevPtr, Atom* nextPtr,
             int prevBType, int nextBType, double &rosenbluth, double &energy);

      void deleteSequence(int nActive, int sign, Atom* endPtr, int *bonds,
                          double &rosenbluth, double &energy);
   
      void addSequence(int nActive, int sign, Atom* beginPtr, int *bonds,
                          double &rosenbluth, double &energy);

   private:
   
      static const int MaxBin_ = 102;

      double l20Table[MaxBin_];

      double angTable[MaxBin_];

      double kappaTable[MaxBin_];

      double normalTable[MaxBin_];
   
      inline void getAngKappaNorm(double l_20, double &preAng,
                                  double &kappa, double &normConst);
  
      inline void orientationBias(Vector u_21, Vector u_20,
          double prefAng, double kappaAng, double &bias);   

      void orientationBiasTable
         (double l_20, double &prefAng, double &kappaAng, double &normConst);
   
      void computeNormalization
         (double prefAng, double kappaAng, double &normConst);
 
   };

   /* 
   * Inline method: get the preferential angle and effective spring constant for
   * given l_20
   */
   void CfbRebridgeBase::getAngKappaNorm(double l_20, double &prefAng,
           double &kappa, double &normConst)
   {
      if (l_20 <= l20Table[0]) {
         prefAng   = angTable[0];
         kappa     = kappaTable[0];
         normConst = normalTable[0];
      } else if (l_20 >= l20Table[MaxBin_-1]) {
         prefAng   = angTable[MaxBin_-1];
         kappa     = kappaTable[MaxBin_-1];
         normConst = normalTable[MaxBin_-1];
      } else {
         int i;
         i = int ( (l_20 - l20Table[0]) / (l20Table[2] - l20Table[1]) );
         prefAng   = angTable[i+1];
         kappa     = kappaTable[i+1];
         normConst = normalTable[i+1];
      }
   }
   
   /* 
   * Inline method: compute the actual angle and the Boltzmann weight.
   */
   void CfbRebridgeBase::orientationBias(
      Vector u_21, Vector u_20, double prefAng, double kappaAng, double &bias)
   {
      double actualAng;
   
      actualAng = acos(u_21.dot(u_20));
      bias = boltzmann(0.5*kappaAng*(actualAng - prefAng)*(actualAng - prefAng));
   }
   
}
#endif