CompositeBond.h

#ifndef SIMP_COMPOSITE_BOND_H
#define SIMP_COMPOSITE_BOND_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 <util/param/ParamComposite.h>
#include <util/param/Begin.h>
#include <util/param/End.h>
#include <util/global.h>

#include <math.h>

namespace Util 
{  class Random; }

namespace Simp
{

   using namespace Util;

   template <class BareBond, class BarePair>
   class CompositeBond : public ParamComposite 
   {
   
   public:
   
      CompositeBond();

      // Mutators

      void setNBondType(int nBondType);

      void readParameters(std::istream &in);
      
      void set(std::string name, int type, double value);

      // Accessors
      
      double energy(double rsq, int type) const;
   
      double forceOverR(double rsq, int bondTypeId) const;

      double randomBondLength(Random* randomPtr, double beta, int type) 
             const;

      double get(std::string name, int type) const;

   private:

      // Bare bond potential object.
      BareBond bond_;

      // Bare pair potential object.
      BarePair pair_;
      
      CompositeBond(const CompositeBond<BareBond, BarePair>& other);
 
   };
     
 
   /* 
   * Constructor.
   */
   template <class BareBond, class BarePair>
   CompositeBond<BareBond, BarePair>::CompositeBond() 
    : bond_(),
      pair_()
   {
      std::string name("CompositeBond");
      name += "<";
      name += bond_.className();
      name += ",";
      name += pair_.className();
      name += ">";
      setClassName(name.c_str());
   }
   
   
   /* 
   * Copy constructor.
   */   
   template <class BareBond, class BarePair>
   CompositeBond<BareBond, BarePair>::CompositeBond(const CompositeBond<BareBond, BarePair>& other) 
    : bond_(other.bond_),
      pair_(other.pair_)
   {}
     
   /* 
   * Read potential parameters from file.
   */  
   template <class BareBond, class BarePair>
   void CompositeBond<BareBond, BarePair>::readParameters(std::istream &in)
   {

      readParamComposite(in, bond_);
      readParamComposite(in, pair_);

      #if 0
      std::string label;
      Begin* beginPtr = 0;
      End*   endPtr = 0;
      // Read bare bond potential
      beginPtr = &addBegin(bond_.className().c_str());
      beginPtr->setIndent(*this);
      beginPtr->readParameters(in);
      readParamComposite(in, bond_);
      endPtr = &addEnd();
      endPtr->setIndent(*this);
      endPtr->readParameters(in);

      // Read bare pair potential
      beginPtr = &addBegin(pair_.className().c_str());
      beginPtr->setIndent(*this);
      beginPtr->readParameters(in);
      readParamComposite(in, pair_);
      endPtr = &addEnd();
      endPtr->setIndent(*this);
      endPtr->readParameters(in);
      #endif

   }
   
   /* 
   * Set nAtomType
   */
   template <class BareBond, class BarePair>
   void CompositeBond<BareBond, BarePair>::setNBondType(int nBondType) 
   {  
      bond_.setNBondType(nBondType);
      pair_.setNAtomType(1);
   }   

   /* 
   * Calculate interaction energy for a bond.
   */
   template <class BareBond, class BarePair>
   inline double 
   CompositeBond<BareBond, BarePair>::energy(double rsq, int type)
   const 
   {
      double total = bond_.energy(rsq, type);
      total += pair_.energy(rsq, 0, 0);
      return total;
   }
  
   /* 
   * Compute force/distance for a bond.
   */
   template <class BareBond, class BarePair>
   inline double 
   CompositeBond<BareBond, BarePair>::forceOverR(double rsq, int type)
   const 
   {
      double total = bond_.forceOverR(rsq, type);
      total += pair_.forceOverR(rsq, 0, 0);
      return total;
   }

   /*
   * Modify a parameter, identified by a string.
   */
   template <class BareBond, class BarePair>
   void CompositeBond<BareBond, BarePair>
        ::set(std::string name, int type, double value)
   {
      UTIL_THROW("Unrecognized parameter name");
   }

   /*
   * Get a parameter value, identified by a string.
   */
   template <class BareBond, class BarePair>
   double CompositeBond<BareBond, BarePair>::
          get(std::string name, int type) const
   {
      UTIL_THROW("Unrecognized parameter name");
      return 0.0;
   }

   /*
   * Throws exception if called.
   */  
   template <class BareBond, class BarePair>
   double CompositeBond<BareBond, BarePair>
          ::randomBondLength(Random* randomPtr, double beta, int type) const
   {  
      UTIL_THROW("Unimplemented function"); 
      return 0.0; // To avoid compiler warnings.
   }

}
#endif