pysawsim/sawsim.py

   1 # Copyright (C) 2010  W. Trevor King <wking@drexel.edu>
   2 #
   3 # This program is free software: you can redistribute it and/or modify
   4 # it under the terms of the GNU General Public License as published by
   5 # the Free Software Foundation, either version 3 of the License, or
   6 # (at your option) any later version.
   7 #
   8 # This program is distributed in the hope that it will be useful,
   9 # but WITHOUT ANY WARRANTY; without even the implied warranty of
  10 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  11 # GNU General Public License for more details.
  12 #
  13 # You should have received a copy of the GNU General Public License
  14 # along with this program.  If not, see <http://www.gnu.org/licenses/>.
  15 #
  16 # The author may be contacted at <wking@drexel.edu> on the Internet, or
  17 # write to Trevor King, Drudge's University, Physics Dept., 3141 Chestnut St.,
  18 # Philadelphia PA 19104, USA.
  19
  20
  21 """`Seminar` for running `sawsim` and parsing the results.
  22 """
  23
  24 from __future__ import with_statement
  25
  26 try:
  27     from collections import namedtuple
  28 except ImportError:  # work around Python < 2.6
  29     from ._collections import namedtuple
  30 import hashlib
  31 from optparse import Option
  32 import os
  33 import os.path
  34 from random import shuffle
  35 import shutil
  36 from uuid import uuid4
  37
  38 from . import __version__
  39 from .manager import MANAGERS, get_manager, InvokeJob
  40
  41
  42 SAWSIM = 'sawsim'  # os.path.expand(os.path.join('~', 'bin', 'sawsim'))
  43 CACHE_DIR = os.path.expanduser(os.path.join('~', '.sawsim-cache'))
  44 DEFAULT_PARAM_STRING = (
  45     '-s cantilever,hooke,0.05 -N1 '
  46     '-s folded,null -N8 '
  47     "-s 'unfolded,wlc,{0.39e-9,28e-9}' "
  48     "-k 'folded,unfolded,bell,{3.3e-4,0.25e-9}' "
  49     '-q folded -v 1e-6')
  50
  51
  52 # `Event` instances represent domain state transitions.
  53 Event = namedtuple(
  54     typename='Event',
  55     field_names=['force', 'initial_state', 'final_state'])
  56
  57
  58 class SawsimRunner (object):
  59     """
  60     >>> from .manager.thread import ThreadManager
  61     >>> m = ThreadManager()
  62     >>> sr = SawsimRunner(sawsim='bin/sawsim', manager=m)
  63     >>> for run in sr(param_string=DEFAULT_PARAM_STRING, N=2):
  64     ...     print 'New run'
  65     ...     for i,event in enumerate(run):
  66     ...         print i, event  # doctest: +ELLIPSIS, +NORMALIZE_WHITESPACE
  67     New run
  68     0 Event(force=..., initial_state='folded', final_state='unfolded')
  69     1 Event(force=..., initial_state='folded', final_state='unfolded')
  70     2 Event(force=..., initial_state='folded', final_state='unfolded')
  71     3 Event(force=..., initial_state='folded', final_state='unfolded')
  72     4 Event(force=..., initial_state='folded', final_state='unfolded')
  73     5 Event(force=..., initial_state='folded', final_state='unfolded')
  74     6 Event(force=..., initial_state='folded', final_state='unfolded')
  75     7 Event(force=..., initial_state='folded', final_state='unfolded')
  76     New run
  77     0 Event(force=..., initial_state='folded', final_state='unfolded')
  78     1 Event(force=..., initial_state='folded', final_state='unfolded')
  79     2 Event(force=..., initial_state='folded', final_state='unfolded')
  80     3 Event(force=..., initial_state='folded', final_state='unfolded')
  81     4 Event(force=..., initial_state='folded', final_state='unfolded')
  82     5 Event(force=..., initial_state='folded', final_state='unfolded')
  83     6 Event(force=..., initial_state='folded', final_state='unfolded')
  84     7 Event(force=..., initial_state='folded', final_state='unfolded')
  85     >>> m.teardown()
  86     """
  87
  88     optparse_options = [
  89         Option('-s','--sawsim', dest='sawsim',
  90                metavar='PATH',
  91                help='Set sawsim binary (%default).',
  92                default=SAWSIM),
  93         Option('-p','--params', dest='param_string',
  94                metavar='PARAMS',
  95                help='Initial params for fitting (%default).',
  96                default=DEFAULT_PARAM_STRING),
  97         Option('-N', '--number-of-runs', dest='N',
  98                metavar='INT', type='int',
  99                help='Number of sawsim runs at each point in parameter space (%default).',
 100                default=400),
 101         Option('-m', '--manager', dest='manager',
 102                metavar='STRING',
 103                help='Job manager name (one of %s) (%%default).'
 104                % (', '.join(MANAGERS)),
 105                default=MANAGERS[0]),
 106         Option('-C','--use-cache', dest='use_cache',
 107                help='Use cached simulations if they exist (vs. running new simulations) (%default)',
 108                default=False, action='store_true'),
 109         Option('--clean-cache', dest='clean_cache',
 110                help='Remove previously cached simulations if they exist (%default)',
 111                default=False, action='store_true'),
 112         Option('-d','--cache-dir', dest='cache_dir',
 113                metavar='STRING',
 114                help='Cache directory for sawsim unfolding forces (%default).',
 115                default=CACHE_DIR),
 116     ]
 117
 118     def __init__(self, sawsim=None, cache_dir=None,
 119                  use_cache=False, clean_cache=False,
 120                  manager=None):
 121         if sawsim == None:
 122             sawsim = SAWSIM
 123         self._sawsim = sawsim
 124         if cache_dir == None:
 125             cache_dir = CACHE_DIR
 126         self._cache_dir = cache_dir
 127         self._use_cache = use_cache
 128         self._clean_cache = clean_cache
 129         self._manager = manager
 130         self._local_manager = False
 131         self._headline = None
 132
 133     def initialize_from_options(self, options):
 134         self._sawsim = options.sawsim
 135         self._cache_dir = options.cache_dir
 136         self._use_cache = options.use_cache
 137         self._clean_cache = options.clean_cache
 138         self._manager = get_manager(options.manager)()
 139         self._local_manager = True
 140         call_params = {}
 141         for param in ['param_string', 'N']:
 142             try:
 143                 call_params[param] = getattr(options, param)
 144             except AttributeError:
 145                 pass
 146         return call_params
 147
 148     def teardown(self):
 149         if self._local_manager == True:
 150             self._manager.teardown()
 151
 152     def __call__(self, param_string, N):
 153         """Run `N` simulations and yield `Event` generators for each run.
 154
 155         Use the `JobManager` instance `manager` for asynchronous job
 156         execution.
 157
 158         If `_use_cache` is `True`, store an array of unfolding forces
 159         in `cache_dir` for each simulated pull.  If the cached forces
 160         are already present for `param_string`, do not redo the
 161         simulation unless `_clean_cache` is `True`.
 162         """
 163         count = N
 164         if self._use_cache == True:
 165             d = self._param_cache_dir(param_string)
 166             if os.path.exists(d):
 167                 if self._clean_cache == True:
 168                     shutil.rmtree(d)
 169                     self._make_cache(param_string)
 170                 else:
 171                     for data in self._load_cached_data(param_string):
 172                         yield data
 173                         count -= 1
 174                         if count == 0:
 175                             return
 176             else:
 177                 self._make_cache(param_string)
 178
 179         jobs = {}
 180         for i in range(count):
 181             jobs[i] = self._manager.async_invoke(InvokeJob(
 182                     target='%s %s' % (self._sawsim, param_string)))
 183         complete_jobs = self._manager.wait(
 184             [job.id for job in jobs.itervalues()])
 185         for i,job in jobs.iteritems():
 186             j = complete_jobs[job.id]
 187             assert j.status == 0, j.data['error']
 188             if self._use_cache == True:
 189                 self._cache_run(d, j.data['stdout'])
 190             yield self.parse(j.data['stdout'])
 191         del(jobs)
 192         del(complete_jobs)
 193
 194     def _param_cache_dir(self, param_string):
 195         """
 196         >>> s = SawsimRunner()
 197         >>> s._param_cache_dir(DEFAULT_PARAM_STRING)  # doctest: +ELLIPSIS
 198         '/.../.sawsim-cache/...'
 199         """
 200         return os.path.join(
 201             self._cache_dir, hashlib.sha256(param_string).hexdigest())
 202
 203     def _make_cache(self, param_string):
 204         cache_dir = self._param_cache_dir(param_string)
 205         os.makedirs(cache_dir)
 206         with open(os.path.join(cache_dir, 'param_string'), 'w') as f:
 207             f.write('# version: %s\n%s\n' % (__version__, param_string))
 208
 209     def _load_cached_data(self, param_string):
 210         pcd = self._param_cache_dir(param_string)
 211         filenames = os.listdir(pcd)
 212         shuffle(filenames)
 213         for filename in filenames:
 214             if not filename.endswith('.dat'):
 215                 continue
 216             with open(os.path.join(pcd, filename), 'r') as f:
 217                 yield self.parse(f.read())
 218
 219     def _cache_run(self, cache_dir, stdout):
 220         simulation_path = os.path.join(cache_dir, '%s.dat' % uuid4())
 221         with open(simulation_path, 'w') as f:
 222             f.write(stdout)
 223
 224     def parse(self, text):
 225         """Parse the output of a `sawsim` run.
 226
 227         >>> text = '''#Force (N)\\tinitial state\\tFinal state
 228         ... 2.90301e-10\\tfolded\\tunfolded
 229         ... 2.83948e-10\\tfolded\\tunfolded
 230         ... 2.83674e-10\\tfolded\\tunfolded
 231         ... 2.48384e-10\\tfolded\\tunfolded
 232         ... 2.43033e-10\\tfolded\\tunfolded
 233         ... 2.77589e-10\\tfolded\\tunfolded
 234         ... 2.85343e-10\\tfolded\\tunfolded
 235         ... 2.67796e-10\\tfolded\\tunfolded
 236         ... '''
 237         >>> sr = SawsimRunner()
 238         >>> events = list(sr.parse(text))
 239         >>> len(events)
 240         8
 241         >>> events[0]  # doctest: +ELLIPSIS
 242         Event(force=2.9030...e-10, initial_state='folded', final_state='unfolded')
 243         >>> sr._headline
 244         ['Force (N)', 'initial state', 'Final state']
 245         """
 246         for line in text.splitlines():
 247             line = line.strip()
 248             if len(line) == 0:
 249                 continue
 250             elif line.startswith('#'):
 251                 if self._headline == None:
 252                     self._headline = line[len('#'):].split('\t')
 253                 continue
 254             fields = line.split('\t')
 255             if len(fields) != 3:
 256                 raise ValueError(fields)
 257             force,initial_state,final_state = fields
 258             yield Event(float(force), initial_state, final_state)
 259
 260
 261 def main(argv=None):
 262     """
 263     >>> try:
 264     ...     main(['--help'])
 265     ... except SystemExit, e:
 266     ...     pass  # doctest: +ELLIPSIS, +REPORT_UDIFF
 267     Usage: ... [options]
 268     <BLANKLINE>
 269     Options:
 270       -h, --help            show this help message and exit
 271       -s PATH, --sawsim=PATH
 272                             Set sawsim binary (sawsim).
 273       ...
 274     >>> print e
 275     0
 276     >>> main(['--sawsim', 'bin/sawsim', '-N', '2'])
 277     ... # doctest: +ELLIPSIS, +NORMALIZE_WHITESPACE
 278     #Force (N)  Initial state  Final state
 279     ...         folded         unfolded
 280     ...         folded         unfolded
 281     ...         folded         unfolded
 282     ...         folded         unfolded
 283     ...
 284     ...         folded         unfolded
 285     """
 286     from optparse import OptionParser
 287     import sys
 288
 289     if argv == None:
 290         argv = sys.argv[1:]
 291
 292     sr = SawsimRunner()
 293
 294     usage = '%prog [options]'
 295     epilog = '\n'.join([
 296             'Python wrapper around `sawsim`.  Distribute `N` runs using',
 297             'one of the possible job "managers".  Also supports caching',
 298             'results to speed future runs.'
 299             ])
 300     parser = OptionParser(usage, epilog=epilog)
 301     for option in sr.optparse_options:
 302         parser.add_option(option)
 303
 304     options,args = parser.parse_args(argv)
 305
 306     try:
 307         sr_call_params = sr.initialize_from_options(options)
 308
 309         first_run = True
 310         for run in sr(**sr_call_params):
 311             if first_run == True:
 312                 first_run = False
 313                 run = list(run)  # force iterator evaluation
 314                 if sr._headline != None:
 315                     print '#%s' % '\t'.join(sr._headline)
 316             for event in run:
 317                 print '\t'.join([str(x) for x in event])
 318     finally:
 319         sr.teardown()