Package csb :: Package statistics :: Package samplers :: Package mc

[frames] | no frames]

# Source Code for Package csb.statistics.samplers.mc

      1  """
      2  Abstract Monte Carlo samplers.
      3  """
      4
      5  import numpy.random
      6
      7  import csb.numeric
      8  import csb.core
      9
     10  from abc import ABCMeta, abstractmethod, abstractproperty
     11  from csb.statistics.samplers import AbstractSampler, AbstractState, State, EnsembleState
     12
     13 -class AbstractMC(AbstractSampler):
     14      """
     15      Abstract Monte Carlo sampler class. Subclasses implement various
     16      Monte carlo equilibrium sampling schemes.
     17
     18      @param state: Initial state
     19      @type state: L{AbstractState}
     20      """
     21
     22      __metaclass__ = ABCMeta
     23
     24 -    def __init__(self, state):
     25
     26          self._state = None
     27          self.state = state
     28
     29 -    def _checkstate(self, state):
     30
     31          if not isinstance(state, AbstractState):
     32              raise TypeError(state)
     33
     34      @abstractproperty
     35 -    def energy(self):
     36          """
     37          Energy of the current state.
     38          """
     39          pass
     40
     41      @property
     42 -    def state(self):
     43          """
     44          Current state.
     45          """
     46          return self._state
     47      @state.setter
     48 -    def state(self, value):
     49          self._checkstate(value)
     50          self._state = value
     51
     52      @abstractmethod
     53 -    def sample(self):
     54          """
     55          Draw a sample.
     56          @rtype: L{AbstractState}
     57          """
     58          pass
     59
     60 -class AbstractPropagationResult(object):
     61      """
     62      Abstract class providing the interface for the result
     63      of a deterministic or stochastic propagation of a state.
     64      """
     65
     66      __metaclass__ = ABCMeta
     67
     68      @abstractproperty
     69 -    def initial(self):
     70          """
     71          Initial state
     72          """
     73          pass
     74
     75      @abstractproperty
     76 -    def final(self):
     77          """
     78          Final state
     79          """
     80          pass
     81
     82      @abstractproperty
     83 -    def heat(self):
     84          """
     85          Heat produced during propagation
     86          @rtype: float
     87          """
     88          pass
     89
     90 -class PropagationResult(AbstractPropagationResult):
     91      """
     92      Describes the result of a deterministic or stochastic
     93      propagation of a state.
     94
     95      @param initial: Initial state from which the
     96                      propagation started
     97      @type initial: L{State}
     98
     99      @param final: Final state in which the propagation
    100                    resulted
    101      @type final: L{State}
    102
    103      @param heat: Heat produced during propagation
    104      @type heat: float
    105      """
    106
    107
    108 -    def __init__(self, initial, final, heat=0.0):
    109
    110          if not isinstance(initial, AbstractState):
    111              raise TypeError(initial)
    112
    113          if not isinstance(final, AbstractState):
    114              raise TypeError(final)
    115
    116          self._initial = initial
    117          self._final = final
    118          self._heat = None
    119
    120          self.heat = heat
    121
    122 -    def __iter__(self):
    123
    124          return iter([self._initial, self.final])
    125
    126      @property
    127 -    def initial(self):
    128          return self._initial
    129
    130      @property
    131 -    def final(self):
    132          return self._final
    133
    134      @property
    135 -    def heat(self):
    136          return self._heat
    137      @heat.setter
    138 -    def heat(self, value):
    139          self._heat = float(value)
    140
    141 -class Trajectory(csb.core.CollectionContainer, AbstractPropagationResult):
    142      """
    143      Ordered collection of states, representing a phase-space trajectory.
    144
    145      @param items: list of states defining a phase-space trajectory
    146      @type items: list of L{AbstractState}
    147      @param heat: heat produced during the trajectory
    148      @type heat: float
    149      @param work: work produced during the trajectory
    150      @type work: float
    151      """
    152
    153 -    def __init__(self, items, heat=0.0, work=0.0):
    154
    155          super(Trajectory, self).__init__(items, type=AbstractState)
    156
    157          self._heat = heat
    158          self._work = work
    159
    160      @property
    161 -    def initial(self):
    162          return self[0]
    163
    164      @property
    165 -    def final(self):
    166          return self[self.last_index]
    167
    168      @property
    169 -    def heat(self):
    170          return self._heat
    171      @heat.setter
    172 -    def heat(self, value):
    173          self._heat = float(value)
    174
    175      @property
    176 -    def work(self):
    177          return self._work
    178      @work.setter
    179 -    def work(self, value):
    180          self._work = float(value)
    181
    182 -class TrajectoryBuilder(object):
    183      """
    184      Allows to  build a Trajectory object step by step.
    185
    186      @param heat: heat produced over the trajectory
    187      @type heat: float
    188      @param work: work produced during the trajectory
    189      @type work: float
    190      """
    191
    192 -    def __init__(self, heat=0.0, work=0.0):
    193          self._heat = heat
    194          self._work = work
    195          self._states = []
    196
    197      @staticmethod
    198 -    def create(full=True):
    199          """
    200          Trajectory builder factory.
    201
    202          @param full: if True, a TrajectoryBuilder instance designed
    203                       to build a full trajectory with initial state,
    204                       intermediate states and a final state. If False,
    205                       a ShortTrajectoryBuilder instance designed to
    206                       hold only the initial and the final state is
    207                       returned
    208          @type full: boolean
    209          """
    210
    211          if full:
    212              return TrajectoryBuilder()
    213          else:
    214              return ShortTrajectoryBuilder()
    215
    216      @property
    217 -    def product(self):
    218          """
    219          The L{Trajectory} instance build by a specific instance of
    220          this class
    221          """
    222          return Trajectory(self._states, heat=self._heat, work=self._work)
    223
    224 -    def add_initial_state(self, state):
    225          """
    226          Inserts a state at the beginning of the trajectory
    227
    228          @param state: state to be added
    229          @type state: L{State}
    230          """
    231          self._states.insert(0, state.clone())
    232
    233 -    def add_intermediate_state(self, state):
    234          """
    235          Adds a state to the end of the trajectory
    236
    237          @param state: state to be added
    238          @type state: L{State}
    239          """
    240          self._states.append(state.clone())
    241
    242 -    def add_final_state(self, state):
    243          """
    244          Adds a state to the end of the trajectory
    245
    246          @param state: state to be added
    247          @type state: L{State}
    248          """
    249          self._states.append(state.clone())
    250
    251 -class ShortTrajectoryBuilder(TrajectoryBuilder):
    252
    253 -    def add_intermediate_state(self, state):
    254          pass
    255
    256      @property
    257 -    def product(self):
    258          """
    259          The L{PropagationResult} instance built by a specific instance of
    260          this class
    261          """
    262
    263          if len(self._states) != 2:
    264              raise ValueError("Can't create a product, two states required")
    265
    266          initial, final = self._states
    267          return PropagationResult(initial, final, heat=self._heat)
    268
    269
    270 -class MCCollection(csb.core.BaseCollectionContainer):
    271      """
    272      Collection of single-chain samplers.
    273
    274      @param items: samplers
    275      @type items: list of L{AbstractSingleChainMC}
    276      """
    277
    278 -    def __init__(self, items):
    279
    280          from csb.statistics.samplers.mc.singlechain import AbstractSingleChainMC
    281
    282          super(MCCollection, self).__init__(items, type=AbstractSingleChainMC)
    283
    284
    285 -def augment_state(state, temperature=1.0, mass_matrix=None):
    286      """
    287      Augments a state with only positions given by momenta drawn
    288      from the Maxwell-Boltzmann distribution.
    289
    290      @param state: State to be augmented
    291      @type state: L{State}
    292
    293      @param temperature: Temperature of the desired Maxwell-Boltzmann
    294                          distribution
    295      @type temperature: float
    296
    297      @param mass_matrix: Mass matrix to be used in the Maxwell-Boltzmann
    298                          distribution; None defaults to a unity matrix
    299      @type mass_matrix: L{InvertibleMatrix}
    300
    301      @return: The initial state augmented with momenta
    302      @rtype: L{State}
    303      """
    304
    305      d = len(state.position)
    306      mm_unity = None
    307
    308      if mass_matrix is None:
    309          mm_unity = True
    310
    311      if mm_unity == None:
    312          mm_unity = mass_matrix.is_unity_multiple
    313
    314      if mm_unity == True:
    315          momentum = numpy.random.normal(scale=numpy.sqrt(temperature),
    316                                         size=d)
    317      else:
    318          covariance_matrix = temperature * mass_matrix
    319          momentum = numpy.random.multivariate_normal(mean=numpy.zeros(d),
    320                                                      cov=covariance_matrix)
    321
    322      state.momentum = momentum
    323
    324      return state
    325