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