proteins are complex molecules of fundamental 
 importance in biological 
 processes . numerical simulation using molecular dynamics ( md )
has 
 proven to be a powerful tool to predict many important properties 
 such as the native state of the protein or its free energy . in this paper , we will focus on methods , based on md , to calculate 
 reaction rates , which are defined as transition rate between metastable 
 states or conformations of the protein . as a byproduct of our analysis
, 
 we will also calculate the main mechanisms of the reaction , i.e. , 
 the transition pathways .
for example , a cartoon model of the free 
 energy for a biomolecule is shown in figure 1 .
the left region represents the reactant states ( r ) , and the right the product states ( p ) .
trajectories 
 spend most of their time in r or p with 
 infrequent transitions due to the energy barrier that separates the 
 two states .
molecular dynamics sampling of a one - dimensional free - energy surface 
 between reactant ( r ) and product ( p ) states . 
 by periodically recording the molecular conformation , energy , and 
 other observables , various relevant properties of the protein can 
 be computed . in 1977
mccammon et al . 
 applied md to the bovine pancreatic trypsin 
 inhibitor ( bpti ) . while the system was 
 simple ( vacuum with a crude force field ) , the simulation nonetheless 
 contributed to shifting the view of proteins as rigid 
 since this initial simulation , md has been used to 
 study a wide variety of topics , such as identification of integral 
 motions such as hinge bending modes , trna 
 flexibility , and the study of e. coli chaperone groel .
the application of md to larger and more complicated molecules impacts 
 development of force - field parameters such as the charmm and amber families , and solvent models 
 ( generaly categorized as implicit or 
 explicitly defined ) .
recently , md has been used to study an hiv capsid and a complete satellite tobacco mosaic virus , as well as screening designed proteins .
in fact , a large portion of the major supercomputers 
 time is currently dedicated to this type of simulation .
the main difficulty 
 is that the time scales of interest are typically in the millisecond 
 ( 10 s ) while a typical time step in md is on the 
 order of the femtosecond ( 10 s ) .
therefore , a 
 brute force simulation would require on the order of 1000 billion 
 time steps , which is impractical given the current hardware .
high - performance computing infrastructures allow development of 
 efficient parallel algorithms to speed up simulations .
specialized hardware , such as mdgrape and anton , provide orders - of - magnitude speedup over traditional high - power 
 computing ( hpc ) simulations .
similarly , work with graphical processing 
 units ( gpus ) have shown the that gpus achieved greater performance 
 over a single cntral processing unit ( cpu ) , allowing a single gpu to simulate 
 biological systems with comparable performance to a cluster .
most approaches to date attempted to accelerate a single ( or 
 a 
 few ) very long simulations .
due to the sequential nature of md simulations , 
 this is clearly a challenging task , which limited scalability .
, much work has gone into the 
 development of fast and efficient algorithms to dedicate large resources 
 to long simulations of molecules .
for instance , namd and amber were 
 among the first to achieve scaling to hundreds of nodes and microsecond - time 
 scales using parallel implementations .
improvements 
 to constraint algorithms , particle interactions , and memory access patterns have resulted in significant 
 performance improvements over the decades . at the other end
, 
 a recent class of methods is attempting to predict 
 equilibrium properties , such as free energy , reaction rate and reaction 
 pathways , from a large number of short trajectories .
in fact , 
 in this paper , we will demonstrate large - scale simulations using parallel 
 computer grids , showing that scalability and performance can be achieved 
 even on slow networks or when using geographically remote compute 
 nodes .
this class of approaches have been applied with great 
 success by 
 several groups , using a method called the markov state model ( msm ) .
the entry ( i , j ) in the matrix stores the probability 
 of a trajectory to go from cell i to cell j after some lag time . from the eigenvalues , eigenvectors , 
 and other properties of the matrix , we can extract all relevant time 
 scales and kinetics of interest .
this type of approach requires 
 running a large number of trajectories 
 of length  and as such is easily parallelizable .
the accuracy 
 of this method is independent of the presence of energy barriers or 
 of the metastability of regions in conformational space .
instead it 
 relies on efficient local sampling to predict time scales that can 
 be orders of magnitude longer than the aggregate simulation time .
one limitation of msm is that its accuracy depends on the choice 
 of . in fact , getting converged rates is very difficult . at 
 short lag times
 , a significant bias is present in the prediction , 
 that is even in the presence of infinite sampling the predicted rate 
 is incorrect .
this can be addressed by increasing the lag time but 
 this has two limitations .
, it is difficult to obtain converged results with respect 
 to the lag time and sampling size .
once a msm has been computed , 
 one can obtain the implicit time 
 scales from the diagonalization of the stochastic transition matrix . 

these implicit time scales are highly sensitive to the time lag , or 
 the time between samples used to construct the transition matrix : 
 at short time lags implicit time scales have a bias that makes them 
 appear much faster than they really are . at large time
lags the statistical 
 error dominates , and it is hard to find a compromise that balances 
 these two errors . in this paper , we considered a variant method , 
 called the weighted 
 ensemble ( we ) approach which has a higher computational cost but convergence 
 is more robust and easier to monitor .
the we method was 
 introduced by huber and kim to address 
 the problem of capturing infrequent reaction events under brownian 
 dyanamics .
further studies have applied it to simulated annealing 
 ( si ) obtaining higher success rate ( order of magnitude ) for finding 
 the global optimum as compared to traditional si , as well as in the investigation of super oxide dismutase 
 and monoclonal antibody nc6.8 .
have extended we and 
 describe several key properties : ( i ) 
 we produces unbiased results for both markovian and non - markovian 
 dynamics , indicating that the method is much more general than previously 
 understood and ( ii ) the method depends on the number of cells , which 
 may be updated dynamically , allowing the system to discover unknown 
 cells without loss of accuracy in the calculation of reaction rates . 

applications of we to toy models and alanine dipeptid1e , coarse - grained models of adenylate cyclase , and the sodium symport mhp2 transporter protein have yielded promising results . in we 
 a large number of parallel walkers ( md trajectories ) are 
 used .
later , we will explain in more details the algorithm , but the 
 basic process is to start a large number of trajectories from regions 
 in a partitioning of the free - energy surface .
after a small number 
 of md integration steps we consider the population of walkers in each 
 cell . then , using a statistically unbiased procedure , we are able 
 to either remove walkers from overpopulated cells or , on the contrary , 
 repopulate cells that are becoming depleted .
this is done through 
 an appropriate subselection ( killing walkers ) or duplication process . 

this ensures that when the simulation restarts the number of walkers 
 in each cell is the same , and is equal to some target number .
then , 
 we resume the simulation and integrate forward in time each walker 
 using md
. this process is repeated until the quantities of interest , 
 such as the reaction rates , which are computed through an appropriate 
 postprocessing procedure , are converged .
this method avoids 
 the slow time scales found in md by ensuring 
 that all cells are equally populated .
this means for example that 
 regions near the transition barrier between two metastable regions 
 are sampled as often as minimum free energy basins .
this method can 
 be interpreted as a way of enhancing the sampling of conformational 
 space .
we note that its efficiency still depends on an appropriate 
 choice of macro - states or cells , although the fact that these cells 
 form a simple partition and can be constructed in a variety of ways 
 ( e.g. , they can be simple voronoi cells ) gives a lot of flexibility 
 to the method .
some of the advantages of we is that it is easy to 
 setup , is inherently massively parallel and scalable , and is unbiased . 

for example , unlike msm , we is guaranteed to converge to the exact 
 result with enough samples .
notwithstanding these properties , 
 awe , like msm and other related 
 methods , remains computationally expensive due to the large number 
 of walkers that need to be run .
this paper will present a software 
 infrastructure we developed to facilitate running such calculations 
 with accelerated convergence .
this is achieved by an appropriate selection 
 of the initial conditions as described in hence the method is known as accelerated weighted ensemble ( awe ) , 
 an extension of we where the initial weights have been approximated 
 to the steady state weights to accelerate the convergence .
we compute 
 them from a transition matrix obtained from the simulation , based 
 on cells as they come from clustering , perhaps using msmbuilder .
awe 
 also uses a simpler and more accurate method for resampling and generating 
 new walker populations through splitting and merging of walkers .
awe allows one to eliminate the bias introduced 
 in msm by generating 
 the statistically exact distribution of walkers inside each cell . 

msm uses the equilibrium distribution for the walkers in each cell , 
 which causes a bias .
instead , awe uses an out - of - equilibrium distribution 
 corresponding to walkers flowing from reactants to product states 
 ( e.g. , folding / unfolding ) and vice versa , thereby producing exact 
 statistics .
we focused on the following 
 main goals when designing the software:1.allow using any md code in a black - box 
 manner , without having to make intrusive changes inside the code . 

this is important since many md codes are available and people are 
 often required to use a particular md code because of some special 
 required capability .
for instance , the following md codes are usable 
 with our framework : amber , charmm , gromacs , lammps , namd , mdynamix , orac , gromos , tinker , desmond , dl_poly.2.dynamically scale to available 
 resources 
 as they are added and removed as availability or funding allows and 
 robustly handle resource failure .
details are discussed later , but 
 the main idea is to support , for example , starting awe using a dedicated 
 cluster locally then adding cloud infrustructure such as amazon ec2.3.use a scripting language 
 for quick 
 prototyping , interactive simulations , and user - friendly codes .
the use of this simple 
 scripting interface also means that it is not difficult to reuse our 
 framework for other methods that require running a large number of 
 short trajectories with some appropriate postprocessing , e.g. , msm 
 and other related methods .
allow using any md code in a black - box 
 manner , without having to make intrusive changes inside the code . 

this is important since many md codes are available and people are 
 often required to use a particular md code because of some special 
 required capability .
for instance , the following md codes are usable 
 with our framework : amber , charmm , gromacs , lammps , namd , mdynamix , orac , gromos , tinker , desmond , dl_poly .
dynamically scale to available 
 resources 
 as they are added and removed as availability or funding allows and 
 robustly handle resource failure .
details are discussed later , but 
 the main idea is to support , for example , starting awe using a dedicated 
 cluster locally then adding cloud infrustructure such as amazon ec2 . use a scripting language 
 for quick 
 prototyping , interactive simulations , and user - friendly codes .
the use of this simple 
 scripting interface also means that it is not difficult to reuse our 
 framework for other methods that require running a large number of 
 short trajectories with some appropriate postprocessing , e.g. , msm 
 and other related methods . in this paper , 

we present a distributed computing implementation 
 of awe , based on workqueue ( wq ) , that is capable of using computing 
 resources with different architectural properties ( e.g. , gpu , cpu ) , 
 as well as different job execution semantics such as dedicated ( hpc 
 grid ) or cycle - scavenging ( condor ) resources . at peak performance , 
 for our validation , awe - wq simulated an aggregate 1.5 ms , with a peak 
 performance of 1000 ns / h , showing that a large amount of resources 
 can be efficiently managed by the system .
similar software packages 
 include the copernicus framework developed 
 by pronk et al . , adaptive markov state 
 models ( adaptive msms ) by bowman et al . , and the weighted ensemble simulation toolkit with parallelization 
 and analysis ( westpa ) developed by zwier et al .
copernicus is a framework for running ensemble molecular 
 dynamics that allows multiple resources to be used , supports multiple 
 project types , and adaptive sampling .
adaptive msms iteratively build 
 msm models to determine the contribution of each state to the slowest 
 kinetic process . based on this information further simulations are 
 run from states based proportionally to their contribution to the 
 uncertainty .
westpa is a software framework for running we simulations . 
 currently under development by lillian chong and dan zuckermann
, westpa 
 has been used to study several systems , such as molecular association and the sodium symporter mhp1 .
the fault - tolerance model , intelligent task scheduling , 
 and caching 
 distinguish awe - wq .
these features allow the program to dynamically 
 handle resource addition and removal , automatically handle worker 
 failures , recover from machine failures , support clusters of heterogeneous 
 computers , and minimize data transfer .
we first describe the 
 algorithm , design challenges , and implementation 
 in the design and implementation section and 
 provide a brief overview of usage .
second , in the results section , we demonstrate that awe - wq meets the criteria 
 described in the design and implementation section as well as provide a validation of the implementation on 
 a nontrivial protein .
finally , we provide public access to the software 
 and the data sets used for the results as well as describe the current 
 limitations and future directions in availability 
 and future directions section .
the we approach 
 proceeds as follows and illustrated in figure 2.1.partition 
 state - space into c cells.2.determine w number 
 of simulations ( walkers ) to maintain in each cell.3.parallel step ( walk ) : run the n = w  c short simulations , 
 assign each of the n final conformations to a cell.4.synchronization barrier 
 ( resample ) : 
 once all the n simulations are assigned , merge and 
 split the walkers to maintain w walkers in each cell 
 and update the associated weights.5.go to step 3 unless weights converged .
parallel step ( walk ) : run the n = w  c short simulations , 
 assign each of the n final conformations to a cell .
synchronization barrier 
 ( resample ) : 
 once all the n simulations are assigned , merge and 
 split the walkers to maintain w walkers in each cell 
 and update the associated weights . go to step 3 unless weights converged .
given a population 
 of walkers with weights wi , we first calculate a mean weight wm = i(wi / ntargetwalkers ) , where ntargetwalkers is the desired number of walkers that we are trying to achieve . 

then , if a walker has a weight greater than wm , it is split into an integer number of walkers of weight wm , with a remainder that is reinserted into 
 the list of walkers to process .
walkers that have a weight less than wm are merged in a statistically exact way to 
 create walkers with a weight greater or equal to wm . by repeating this process of splitting and merging 
 walkers
, we are able to generate a population of walkers with weights 
 exactly equal to wm , in a statistically 
 exact manner .
this is different from the procedure in the work of 
 huber and kim that generates walkers 
 with weights between wm/2 and 2wm but not exactly equal to the target weight wm .
the choice of method for defining cells 
 is orthogonal to we . as 
 long as a cell definition exists , in cells.dat currently , awe - wq is 
 able to use it .
this allows one to define cells arbitrarily to examine 
 the effect of a cell definition method .
we based the cells on the 
 msm construction in order to accelerate the convergence of cell weights . 
 as explained previously ,
the main reason to use we is to remove the 
 bias present in msm models due to the lag time and the ( incorrect ) 
 distribution of walkers inside each cell . upon reaching steady - state 

the weight of each cell converges to 
 its probability , allowing this procedure to calculate free energy . 

each region represents a cell on 
 the energy surface , the circles the walking simulations , the colors 
 the association with the reactant ( r ) or product ( p ) regions .
calculation of transition 
 rates requires a further modification : 
 define two sets of states , r and p , associated with the reactant and product regions ( figure 2 ) .
each walker is assigned a color , such as blue 
 for r and red for p , corresponding 
 to the previously visited region . at each step in the simulation whenever 
 a blue walker enters p , its color changes to red , 
 and vice versa . the rate from r to p is then directly obtained by computing the flux of blue particles 
 that change color , and similarly for the p to r rate . if extended to multiple colors and sets multiple 
 kinetic rates can be determined . in practical terms
, the r and p regions may correspond to the unfolded and 
 folding regions and the fluxes to the folding and unfolding rates . 
 by providing an initial approximation of the free energy the convergence 
 to steady
the method requires the following 
 three steps ( as illustrated in figure 3 ) : ( 1 ) 
 most of the work is done in parallel as a large number of independent 
 short md trajectories . ( 2 ) the parallel barrier : wait for all md steps 
 to complete then collect walker assignments and final positions , split 
 and merge walkers , and update weights .
each walker is then , in 
 parallel , executed and assigned to a cell .
the join waits for all 
 walkers to finish , before the resample procedure is applied .
if the 
 system has converged then the program halts , others walkers are rerun . in order 
 to use awe - wq one must define the cells 
 from which to run the walkers .
the general protocol , illustrated in 
 figure 4 , is as follows : ( 1 ) sample the search 
 space , ( 2 ) determine cell definitions , ( 3 ) prepare input files , ( 4 ) 
 run awe - wq , ( 5 ) monitor progress . first 
 some sampling must be done to explore 
 state space .
this can be accomplished with programs such as gromacs , 
 charmm , amber , namd , or using infrastructure such as folding@home . 

for this instance , sampling of alanine dipeptide has been run previously , 
 which can be extracted : to start , 
 we have collected 
 md data and stored them in the xtc directory .
the directory structure 
 is xtc / run#/frame0.xtc , where # is the trajectory number .
the sampling data is imported into msmbuilder 
 using the convertdatatohdf.py command , paramterized with a pdb file 
 containing the system coordinates and the path to the md data : diagram of the workflow 
 to run awe - wq with associated programs 
 and descriptions .
first a sampling must be run with some md software 
 such as gromacs , charmm , amber , etc .
prepare the input to awe - wq by converting the 
 files and defining metastable states .
run awe - wq by starting the master 
 process and allocating workers with various resources .
conformations are clustered 
 with 
 a hybrid k - centers / k - medoids algorithm 
 using the rmsd between atoms as the distance metric .
using a subset ( such as all non - hydrogens ) prevents 
 too fine a granularity from overfitting the data .
finally , we will 
 cluster the data into 100 groups . by inspecting the implied time scales ( not shown ) 

extract 
 a given number of 
 conformations from each of the cells defined above ( savestructures.py ) 
 and import the cell definitions from msmbuilder format to awe - wq format 
 ( awe - import - gens ) and setup the initial files for the awe - wq run ( awe - prepare ) . 

regions of metastable states need to then be determined , which are 
 ultimately given as a list of the cells belonging to each region ( e.g. , 
  folded  and  unfolded  ) and scripts are 
 given to do so using rmsd to a reference structure ( awe - calc - gens - rmsd , 
 awe - define - regions ) .
we plan to maintain 10 simulations in the cells , 
 so we need to extract conformations from the states using msmbuilder .
in order to run awe - wq 
 we must then import the cell 
 definitions which were written by msmbuilder to gens.lh5 . in order 
 to compute the fluxes we need specify regions of metastable 
 states .
awe - wq provides some commands to assist with this process : awe - calc - gens - rmsd and awe - define - regions .
we use awe - calc - gens - rmsd to compute the 
 rmsd of each cell to some reference conformation such as the native 
 state . by 
 plotting the distribution of values we can classify conformations 
 with rmsd  2.3  as folded and those 
 with rmsd  2.5  as unfolded .
the two 
 output files folded.dat and unfolded.dat now contain the integer indices 
 of the states belonging to these regions .
we can now prepare 
 for awe - wq by checking dependencies and populating 
 the directory with other necessary files by running awe - prepare .
there are two components to consider 
 when running awe - wq : the master process and the resources .
the master 
 is the driver of the algorithm , managing task definitions , scheduling , 
 processing , and the resampling procedure . in order to run the walkers
this 
 process loads the initial conformations ( walkers ) for each cell , and 
 is in charge of scheduling each task , processing the result , and the 
 resampling procedure .
this runs awe - wq maintaining 10 walkers in 100 
 cells , whose definition is provided in cells.dat with initial weights in data / populations.dat .
the metastable regions are provided in folded.dat and unfolded.dat as a list of cell i d numbers 
 belonging to each region .
finally , we give a name to the master (  awe - wq  ) 
 to that workers can easily locate the host and port .
resources 
 can be allocated either directly 
 using work_queue_worker or managed automatically 
 using work_queue_pool . using work_queue_worker also allows the worker to operate as a  foreman  , 
 enabling the hierarchical distribution of tasks .
the work_queue_pool program maintains workers for the master based on need and is in 
 charge of submission to various backends such as the local machine , 
 condor , sge , pbs , etc .
use work_queue_status to get the current resources runtime status ( number of workers , 
 number of tasks waiting / completed , etc . ) . by using awe - plot - wq - stats the plot of the resource usage over the runtime of the program can 
 be obtained .
additionally , the current status of the master , 
 such as workers busy and tasks complete can be viewed using the work_queue_status command .
the awe - flux command allows the convergence 
 of the calculated flux to be monitored .
for instance , the 
 energy surface for alanine dipeptide can be visualized as a function 
 of its dihedral angles .
as such , we can plot , as shown in figure 5 , the cell coordinates and the initial estimation 
 of the weights as well as computed weights after several iterations 
 of awe - wq .
since the awe method follows 
 a pattern of fan - out followed by fan - in we implement using a master / worker 
 paradigm .
awe - wq is designed to ( 1 ) support off - the - shelf md backend 
 software , ( 2 ) dynamically scale according to resource availability , 
 ( 3 ) support heterogeneous runtime environments , and ( 4 ) be robust 
 in the face of system failure .
application of awe to a nontrivial 
 biological system may have on the order of 10  000 simulations 
 per iteration .
on one cpu , if each walker requires 30 min to simulate , 
 7 months would be required to run a single iteration of awe . since 
 walkers are short , assuming 30 min per walker , 7 months would be required 
 to run 10  000 walkers for each awe iteration using a single 
 processor . with 1000 cpus , the expected time drops to 5 h , assuming 
 uniform performance .
additionally , the synchronization barrier allows 
 straggling workers to slow down the entire application .
since we expect 
 to need hundreds of iterations to converge , the challenge is to implement 
 awe such that it scales to 1000 + processors and can be sustained for 
 months with no major slowdown due to stragglers . rather than 
 reimplement the core md algorithms which are present in multiple software 
 packages ( i.e. , gromacs , amber , charmm , etc . )
this flexibility presented several challenges 
 by imposing a requirement of specifying the steps needed to run the 
 backend , the transfer of files , the execution environment , and the 
 interaction with the other steps of the procedure .
the approach taken 
 in the following : a walker is described as a task which consists of 
 a program to run and the input and output files . once specified , this 
 task is scheduled to run on an available worker , which is a process 
 that accepts any input files , executes the task s command , 
 and returns the results upon request .
the awe - prepare script sets 
 up files and scripts in the local directory necessary to run awe - wq .
one of these scripts , execute-task.sh , executes the commands on the worker . in order to use a different 

md backend execute-task.sh needs to define 
 how to run the md , and awe - wq needs to specify 
 the files to send to the worker .
our solution 
 to the worker fault - tolerance problem implies the elasticity  the 
 ability to dynamically scale to available resources  of the 
 program .
if a task fails during execution , it is rescheduled to run 
 on the next available worker .
as resources are added or removed , tasks 
 are started or reschedule , respectively .
this allows a project to 
 be started with the campus cluster , then ec2 machines run until budget 
 exhaustion , then continued with only the cluster .
managing the resource 
 pools are the interactions between the master , the catalog server 
 ( cs ) , and workqueue pool ( wqp ) processes .
the cs allows workers to 
 dynamically determine the location of the master by mapping a project 
 name ( e.g. , awe - wq ) to a hostname and port 
 ( e.g. , fah.crc.nd.edu:1024 ) .
the master registers 
 this information when started , which allows it to be restarted in 
 the event of machine failure without requiring the resources to be 
 reallocated .
in addition , the cs store information about the current 
 resource requirements of a master , such as the number of tasks waiting 
 to be scheduled . a wqp can thus use this information to automatically 
 start workers based on the runtime needs . this is further enabled 
 by the asynchronous model , which treats tasks as independent , arbitrarily 
 executable entities .
additionally , awe - wq 
 can use heterogeneous resources , allowing systems with different properties 
 and ( potentially ) operational semantics to be used .
for instance , 
 this supports using gpu - based as well as cpu - based md backend codes . 

this is accomplished by extending the task abstraction to dynamically 
 resolve file location upon worker connection rather than statically . 

for example , an input file can be specified for a task as binaries/$os/$arch / mdrun .
when a worker connects to the 
 master it reports the operating system and architecture of the remote 
 node in the $ os and $ arch variables , which may then be expanded to send binaries / linux / i686/mdrun and binaries / darwin / x86_64/mdrun to to workers 
 running on a 32-bit linux and 64-bit mac os machine , respectively . 
 finally , these variables can be overridden by the user to provide 
 additional granularity .
awe - wq is implemented 
 with the assumption 
 that machine failures are not only possible but in fact common , which 
 is the opposite of the paradigm underlying mpi infrastructures . as 
 such two broad categories of failures need to be accounted for : failure 
 of the worker or master processes . in the first case ,
worker failure 
 are handled transparently and tasks rescheduled to the next available 
 worker .
failure of the master node would normally result in 
 the termination with an inability to resume .
this is overcome by using 
 a transactional logging mechanism : a combination of persistent and 
 volatile logs .
each result from a walker is added to a log and after 
 each awe iteration the state of the system is checkpointed to disk 
 and the log cleared . if the process fails , then the state is initialized 
 from the checkpoint and the log used to replay and missed results 
 from that iteration .
the second involves task replication 
 to remove the effect of straggling workers . by monitoring the state 
 of the workers
 either busy or idle ( wi)and the number of tasks waiting 
 to be assigned a worker , task replication can be engaged when wi > 0 .
once a task returns , 
 all its replicates are canceled , results are accepted in a first - come - first - serve 
 fashion . by combining replication with preferential scheduling ,
workers 
 are continuously used and slow workers allowed to contribute results 
 until the fast machines are available , which then run the replicated 
 tasks .
the inner loop of awe - wq is an event - loop 
 with asynchronous task execution .
each iteration begins by translating 
 the internal data model of the walkers into task , which are submitted 
 to the wq scheduler . while tasks are executing , the program is idle 
 until a result is returned .
each result creates an event that must 
 be handled , either by translation or resubmission .
handlers are lightweight : 
 read the final coordinates and cell assignment into the data models . 

since task executions are asynchronous and results are processed serially , 
 the event model enforces load balancing by spreading out the calls 
 to the event handlers .
increasing 
 scalability is further 
 accomplished by minimizing the communication overhead and maximizing 
 the parallel work .
ideally , the master process will spend most of 
 its time waiting . by categorizing files as those common to all tasks 
 ( common ) or specific to the current one ( unique )
, the common files can be cached on the workers . 
 since we expect the number of walkers to be larger than resource availability 
 the caching mechanism shifts overhead to depend on resources rather 
 than problem size .
for instance , using a gromacs backend each task 
 requires 34 mb in input files , binaries , and miscellaneous files . in the case of the ww domain simulation
, each task required 34 mb 
 of gromacs - related files as overhead and 100 kb of task specific data . 
 assuming 10  000 tasks / iteration and 500 iterations , the total 
 amount of data transferred is 34 ( mb / task )  10  000 ( tasks / iteration ) 

 500 iterations  2 ( tb / mb ) = 162 
 tb . by caching input files , data transferred depends only on the number 
 of new workers seen , rather than total number of tasks , which may 
 be orders of magnitude greater . in this case , assuming 1000 workers 
 are connected the total data transferred becomes 34 mb  1000 
  2 ( gb / mb ) = 33 gb .
we have run awe - wq for several months where 3.5 million tasks have 
 executed over 600 compute years in a heterogeneous environment simulating 
 over 1.5 ms of time . in this section
available resources consisted 
 of the following : notre dame hpc grid with 6000 cores shared among 
 campus users ; notre dame condor pool with variable ( usually around 
 1200 ) cores and flocking capabilities with purdue university and the 
 university of wisconsin  madison ; stanford university s 
 institute for computational and mathematical engineering ( icme ) cluster 
 with 200 cpus and 100 nvidia gpus , and cloud virtual machines via 
 amazon ec2 and microsoft azure .
the gpu cluster executes tasks within minutes , while the nd hpc grid 
 and ec2 machines are the best performant cpu - based pools .
the range 
 of performance of condor machines is unsurprising and illustrative 
 of the cycle - scavenging nature of the platform .
the azure machines 
 were the least performant resources , partly due to the necessity of 
 running the tasks via cygwin .
distribution of task execution times indicate 
 several peaks associated 
 with the different underlying resources being used .
figure 7 demonstrates the fault - tolerance elasticity and 
 scalability 
 of the application as an expanded view of the number of busy resources 
 over the runtime  the inset shows the entire run .
unexpected 
 failures of the master ( months 4 , 6 ) initiated development and integration 
 of features , bug fixes , and monitoring capabilities . due to power 
 outages and infrastructure upgrades the master process restarted multiple 
 times . in spite of these failures
over 250 iterations were completed 
 without data loss due to the transactional logging mechanism .
the 
 periodically persisting state additionally allowed us to integrate 
 new features , such as task replication , and bug fixes over the course 
 of the project s execution .
focusing on the region between 
 days 865 and 900 shows a steady decrease in the number of busy resources 
 until a new pool of workers became available and awe - wq automatically 
 started using them .
it must be emphasized as well that minimal human 
 input was needed to take advantage of the additional resources : due 
 to the use of the catalog server all that was required was to start 
 a wq pool process for the additional resources .
the periodic 
 dips , every 5 days in the first half and every 20 
 h for the second , correspond to the global synchronization barrier , 
 each period requiring 12  000 tasks to complete . by scheduling 
 the tasks based on resource availability and progression through the 
 awe algorithm ,
resources are fully utilized and the long tail effect 
 from straggling workers is eliminated . due to the use of the caching 
 mechanism communication overhead is reduced : while the entire task 
 payload ( common and unique files ) is 34 mb it is only for the first 
 task a worker executes that his is transferred . each subsequent task 
 only requires 100 kb of data to be transferred ( both in- and outbound 
 traffic ) .
the ww domain is a 34 residue protein 
 domain with two conserved tryptophan ( w ) residues and comprised of 
 two antiparallel -sheets as shown in figure 8 .
simulation parameters were the amber96 force 
 field with generalized born implicit solvent , a viscosity parameter 
 (  ) of 1 ps , temperature of 395 k ( close 
 to the optimal folding temperature ) , 
 saving conformations every 1 ns .
the simulation was run for 200 s 
 during which multiple folding and unfolding events were observed . 

the conformations sampled by the simulation were then clustered using 
 the  and  carbons using a cutoff of 3  . in order 
 to reduce the time spent sampling the unfolded space we established 
 that conformations with an rmsd to the folded structure less than 
 3 and 6  , and greater as folded , intermediate , and unfolded , 
 respectively . using these definitions the number of unfolded states 
 was reduced , for a total of 601 cells .
the awe - wq parameters were 
 thus 601 cells , 20 walkers per cell , and walker length of 500 ps using 
 the same md parameters as the source data .
expanded view of the entire run ( inset ) demonstrating 
 elasticity , 
 scalability , and fault - tolerance . in the inset ,
long gaps ( such as 
 the sudden failures at month 4 and 6 ) in activity indicate times when 
 major development was done on the code ( features added , bugs fixed ) . 
 despite the troubles experienced by the master node , calcuation of 
 the results resumed once the code was deployed . plotting the number 
 of busy workers over many days shows the ability of the software to 
 adapt to a changing environment .
the periodic dips indicate the global 
 synchronization barrier , each requiring 12  000 tasks to complete . 
 a large increase after day 885 occurred when many nodes were suddenly 
 available and were automatically incorporated .
over the course of many months awe - wq completed 
 345 iterations , with an aggregate 2.3 ms of time simulated using 3.6 
 million tasks , and a peak performance of 1000 ns / h .
analysis indicates 
 that the awe - wq results converge within acceptable error tolerance . 

table 1 displays the forward and backward rate 
 estimations computed using awe - wq and the 200 s simulation . 

the awe estimation lies within the confidence interval  with 
 smaller statistical error  than the brute force simulation .
figure 9 illustrates that forward and backward 
 fluxes took 30 iterations to converge to the confidence interval built 
 from brute force simulation and are stabilized within the interval 
 afterward .
the resampling procedure of awe introduces no bias in the 
 estimation of transition rates while exhibiting smaller statistical 
 error , showing an improvement over the brute force method .
transitions among the states are examined with a time lag of  
 = t = 0.5 ns .
three significant pathways from 
 extended to folded state are extracted from the network and shown 
 in figure 10 .
starting from an unfolded conformation , 
 two pathways involve multiple helical intermediate structures before 
 -sheet formation occurs via rearrangement .
the third pathway 
 occurs as a rapid collapse before the formation of loop 2 .
loop 1 
 is then able to form , after which the molecule undergoes further refinement 
 to arrive at the native state .
average forward ( top ) and backward ( bottom ) 
 flux from awe ( blue 
 curve ) with error bar of one standard deviation estimated by the block 
 averaging method , compared with confidence interval built from the 
 brute - force simulation .
the red solid line represents the mean flux 
 from the brute - force simulation , and the dashed lines are mean flux 
 plus / minus one standard deviation .
first three 
 significant pathways from extended to folded conformations . 
 two pathways ( top , middle ) involve helix formation as an intermediary 
 while the third pathway ( bottom ) consists of a collapse then refinement .
awe - wq 
 meets our criteria : allowing the use of off - the - shelf components 
 for md backends , pooling of heterogeneous resources , scalability to 
 over 1000 workers , and ability to cope with failures of both worker 
 and master processes .
the ww domain was used to validate awe - wq by 
 applying the method to a long md simulation in which multiple folding 
 events are observed .
the forward and backward reaction rates computed 
 by awe - wq are within the range and have smaller error , of the rates 
 calculated from the reference simulation .
 a limitation of traditional molecular 
 dynamics ( md ) is that reaction 
 rates are difficult to compute . 
 this is due to the rarity of observing 
 transitions between metastable states since high energy barriers trap 
 the system in these states . 
 recently the weighted ensemble ( we ) family 
 of methods have emerged which can flexibly and efficiently sample 
 conformational space without being trapped and allow calculation of 
 unbiased rates . however , while we can sample correctly and efficiently , 
 a scalable implementation applicable to interesting biomolecular systems 
 is not available . 
 we provide here a gplv2 implementation called awe - wq 
 of a we algorithm using the master / worker distributed computing workqueue 
 ( wq ) framework . 
 awe - wq is scalable to thousands of nodes and supports 
 dynamic allocation of computer resources , heterogeneous resource usage 
 ( such as central processing units ( cpu ) and graphical processing units 
 ( gpus ) concurrently ) , seamless heterogeneous cluster usage ( i.e. , 
 campus grids and cloud providers ) , and support for arbitrary md codes 
 such as gromacs , while ensuring that all statistics are unbiased . 
 
 we applied awe - wq to a 34 residue protein which simulated 1.5 ms over 
 8 months with peak aggregate performance of 1000 ns / h . 
 comparison 
 was done with a 200 s simulation collected on a gpu over a 
 similar timespan . 
 the folding and unfolded rates were of comparable 
 accuracy . 
Introduction
Design and Implementation
Results
Conclusions