IMP  2.1.0
The Integrative Modeling Platform
WormLikeChain.h
Go to the documentation of this file.
1 /**
2  * \file IMP/misc/WormLikeChain.h
3  * \brief Worm-like-chain score for polymer chains.
4  *
5  * Copyright 2007-2013 IMP Inventors. All rights reserved.
6  */
7 
8 #ifndef IMPMISC_WORM_LIKE_CHAIN_H
9 #define IMPMISC_WORM_LIKE_CHAIN_H
10 
11 #include <IMP/misc/misc_config.h>
12 #include <IMP/UnaryFunction.h>
13 #include <IMP/kernel/internal/constants.h>
14 #include <IMP/kernel/internal/units.h>
15 
16 IMPMISC_BEGIN_NAMESPACE
17 
18 //! Worm-like-chain energy for polymer chains
19 /** This function implements one polymer force/extension curve. It
20  is a physical energy and all the inputs are in angstroms
21  and the outputs in kcal/mol (/angstrom).
22 
23  \note The actual worm like chain force blows up as the extension approaches
24  the maximum. Since that makes optimization problematic, we
25  just linearly extend the force after 99% of maximum.
26  */
27 class WormLikeChain : public UnaryFunction {
28  public:
29  //! Define the energy term
30  /** \param[in] l_max maximum length of the chain in angstroms
31  \param[in] lp persistence length in angstroms
32  */
33  WormLikeChain(Float l_max, Float lp) : lmax_(l_max), lp_(lp) {
34  IMP_USAGE_CHECK(l_max > lp, "The persistence length should be less "
35  << "than the total length for this model");
36  }
37 
38  virtual DerivativePair evaluate_with_derivative(double feature) const;
39 
40  virtual double evaluate(double feature) const;
41 
43 
44  void do_show(std::ostream &out) const;
45 
46  private:
47  unit::Piconewton cderiv(unit::Angstrom l) const {
48  unit::Piconewton pn = IMP::internal::KB *
49  IMP::internal::DEFAULT_TEMPERATURE / lp_ *
50  (.25 / square(1.0 - (l / lmax_)) - .25 + (l / lmax_));
51  return pn;
52  }
53 
54  unit::Picojoule eval(unit::Angstrom m) const {
55  unit::Picojoule J =
56  IMP::internal::KB * IMP::internal::DEFAULT_TEMPERATURE / lp_ *
57  (.25 * square(lmax_) / (lmax_ - m) - m * .25 + .5 * square(m) / lmax_);
58  return J;
59  }
60 
61  unit::Angstrom cutoff() const { return .99 * lmax_; }
62 
63  unit::Angstrom lmax_, lp_;
64 };
65 
66 #ifndef IMP_DOXYGEN
67 inline double WormLikeChain::evaluate(double v) const {
68  return evaluate_with_derivative(v).first;
69 }
70 
71 inline DerivativePair WormLikeChain::evaluate_with_derivative(double v) const {
72  static const unit::Picojoule zero = eval(unit::Angstrom(0));
73  unit::Angstrom l(v);
74  if (l < unit::Angstrom(0)) l = unit::Angstrom(0);
75  unit::Picojoule ret;
76 
77  unit::Piconewton doubled;
78  if (l < cutoff()) {
79  ret = (eval(l) - zero);
80  doubled = cderiv(l);
81  } else {
82  unit::Picojoule springterm = (l - cutoff()) * cderiv(cutoff());
83  ret = (eval(cutoff()) + springterm - zero);
84  doubled = cderiv(cutoff());
85  IMP_LOG_VERBOSE("Overstretched " << cderiv(cutoff()) << " " << doubled
86  << " " << l << " " << lmax_ << " "
87  << cutoff() << std::endl);
88  }
89  unit::YoctoKilocalorie zc = convert_J_to_Cal(ret);
90  double value = (zc * unit::ATOMS_PER_MOL).get_value();
91  // std::cout << "Force is " << doubled << " for length " << l << std::endl;
92  // convert from picoNewton
93  unit::YoctoKilocaloriePerAngstrom du = unit::convert_J_to_Cal(doubled);
94 
95  double deriv = (du * unit::ATOMS_PER_MOL).get_value();
96  // std::cout << "Which converts to " << d << std::endl;
97  return std::make_pair(value, deriv);
98 }
99 
100 inline void WormLikeChain::do_show(std::ostream &out) const {
101  out << "params " << lmax_ << " " << lp_ << std::endl;
102 }
103 #endif
104 
105 IMPMISC_END_NAMESPACE
106 
107 #endif /* IMPMISC_WORM_LIKE_CHAIN_H */
Import IMP/kernel/UnaryFunction.h in the namespace.
Abstract single variable functor class for score functions.
WormLikeChain(Float l_max, Float lp)
Define the energy term.
Definition: WormLikeChain.h:33
#define IMP_USAGE_CHECK(expr, message)
A runtime test for incorrect usage of a class or method.
#define IMP_OBJECT_METHODS(Name)
Define the basic things needed by any Object.
Worm-like-chain energy for polymer chains.
Definition: WormLikeChain.h:27
std::pair< double, double > DerivativePair
A pair representing a function value with its first derivative.
Definition: base/types.h:23
double Float
Basic floating-point value (could be float, double...)
Definition: base/types.h:20
#define IMP_LOG_VERBOSE(expr)