Updated simulate.py function; renamed examples

dnarbd/examples/run.py → dnarbd/examples/from_cadnano.py

@@ -4,7 +4,7 @@ from glob import glob
 import re
 
 from dnarbd.readers import cadnano_reader()
-from dnarbd.simulate import run_simulation_protocol
+from dnarbd.simulate import run_multiresolution_simulation as simulate
 
 if __name__ == '__main__':
     if len(sys.argv) > 1:
@@ -22,13 +22,11 @@ if __name__ == '__main__':
             print("WARNING: skipping unreadable json file {}".format(f))
             continue
             
-        gpu = 0
         directory = "run.out.7bd92d9.twist-angles"
-        run_simulation_protocol( out, "job-"+out+"-", model, gpu=gpu,
-                                 directory = directory)
-        # try:
-        #     run_simulation_protocol( out, "job-"+out+"-", model, gpu=gpu,
-        #                              directory = directory)
-        # except:
-        #     print("WARNING: failed to simulate json file {}".format(f))
-        #     continue
+        simulate( model,
+                  output_name = out,
+                  coarse_steps = 1e7,
+                  fine_steps = 1e7,
+                  gpu=1,
+                  directory = directory
+        )

dnarbd/examples/segmentmodel_example.py → dnarbd/examples/from_code.py

 import numpy as np
 from dnarbd import SegmentModel, SingleStrandedSegment, DoubleStrandedSegment
+from dnarbd.simulate import run_multiresolution_simulation as simulate
 
 """
-Example of using segmentmodel to construct a CG
-
-May need to update your shell PATH variable before running this script, e.g.:
-PATH=/home/cmaffeo2/miniconda3/bin:$PATH
+Example using segmentmodel API to construct a model
 
 """
 
-
-
 if __name__ == "__main__":
-
+    
+    """ First create the individual helices and strands """
     seg1 = DoubleStrandedSegment("strand", num_bp = 46)
 
     seg2 = SingleStrandedSegment("strand", 
@@ -24,20 +21,37 @@ if __name__ == "__main__":
                                  end_position = seg1.end_position + np.array((-5,0,5)),
                                  num_nt = 128)
 
+    """ Each "Segment" will has a "start" and "end"
+    - start3 and end5 are an attribute for all Segments
+    - start5 and end3 are an attribute for only DoubleStrandedSegment
+    The 3/5 refers to 3' or 5' end of a DNA strand
+    """
+
     seg1.start3
     seg1.start5
     seg1.end3
     seg1.end5
 
+    """ Add intrahelical connections """
     seg1.connect_end3(seg2)     # equivalent for ssDNA to: seg1.connect_end3(seg2.end5)
     seg1.connect_end5(seg3)
     seg1.connect_start3(seg3)
 
+    """ Add crossover connections """
+    ...
+
+
+    """ Create a model """
     model = SegmentModel( [seg1, seg2, seg3], 
                           local_twist = True,
                           max_basepairs_per_bead = 5,
                           max_nucleotides_per_bead = 5,
                           dimensions=(5000,5000,5000),
                       )
-
-    model.simulate( outputPrefix = 'strand-test', outputPeriod=1e4, numSteps=1e6, gpu=1 )
+    
+    """ Simulate """
+    simulate( model,
+              output_name = "strand-test", 
+              coarse_steps = 1e7,
+              fine_steps = 1e7,
+              gpu=1 )

dnarbd/simulate.py

@@ -5,19 +5,28 @@ from .coords import readArbdCoords, readAvgArbdCoords
 arbd="/home/cmaffeo2/development/cuda/arbd.dbg/src/arbd" # reduced the mem footprint cause vmd
 namd="/home/cmaffeo2/development/namd-bin/NAMD_Git-2017-07-06_Linux-x86_64-multicore-CUDA/namd2"
 
-def multiresolution_simulation( output_name, job_id, model,
-                                remove_long_bonds=False,
+def multiresolution_simulation( model, output_name, 
+                                job_id=None,
                                 gpu = 0,
                                 directory=None,
-                                coarse_steps = 5e7+1,
-                                fine_steps = 5e7+1
+                                coarse_steps = 5e7,
+                                fine_steps = 5e7,
+                                coarse_output_period = 1e5,
+                                fine_output_period = 1e5,
+                                remove_long_bonds=False,
+                                backbone_scale=1.0,
                             ):
 
+    ## Round steps up to nearest multiple of output_period, plus 1
+    coarse_steps = ((coarse_steps//coarse_output_period)+1)*coarse_output_period+1
+    fine_steps = ((fine_steps//fine_output_period)+1)*fine_output_period+1
+
     ret = None
     d_orig = os.getcwd()
     try:
         if directory is None:
-            directory = tempfile.mkdtemp(prefix='origami-%s-' % job_id, dir='/dev/shm/')
+            tmp = job_id if job_id is not None else output_name
+            directory = tempfile.mkdtemp(prefix='origami-%s-' % tmp, dir='/dev/shm/')
         elif not os.path.exists(directory):
             os.makedirs(directory)
         os.chdir(directory)
@@ -27,8 +36,7 @@ def multiresolution_simulation( output_name, job_id, model,
         # TODO: add the output directory in to the readArbdCoords functions or make it an attribute of the model object 
 
         """ Coarse simulation """
-        # simargs = dict(timestep=200e-6, outputPeriod=1e5, gpu=gpu )
-        simargs = dict(timestep=200e-6, outputPeriod=1e5, gpu=gpu, arbd=arbd )
+        simargs = dict(timestep=200e-6, outputPeriod=coarse_output_period, gpu=gpu)
 
         model._clear_beads()
         model._generate_bead_model( 5, 5, local_twist=False )
@@ -53,6 +61,7 @@ def multiresolution_simulation( output_name, job_id, model,
         output_prefix = "%s-2" % output_name
         full_output_prefix = "%s/%s" % (output_directory,output_prefix)
         simargs['timestep'] = 50e-6
+        simargs['outputPeriod'] = fine_output_period
         model._update_segment_positions(coordinates)
         model._clear_beads()
         model._generate_bead_model( 1, 1, local_twist=True, escapable_twist=True )
@@ -75,8 +84,7 @@ def multiresolution_simulation( output_name, job_id, model,
         full_output_prefix = "%s/%s" % (output_directory,output_prefix)
         model._update_segment_positions(coordinates)
         model._clear_beads()
-        # model._generate_atomic_model(scale=0.25)
-        model._generate_atomic_model(scale=1.0)
+        model._generate_atomic_model(scale=backbone_scale)
         model.atomic_simulate( outputPrefix = output_prefix )
 
         ret = directory