diff --git a/bin/mrdna b/bin/mrdna
index ed63432009bccfbc580eab3d2c2fe4bbb0d4ae4c..4364474450fe9488b9f9f623f42d558e216d2b45 100755
--- a/bin/mrdna
+++ b/bin/mrdna
@@ -24,6 +24,9 @@ parser.add_argument('--fine-steps', type=float, default=1e7,
 parser.add_argument('--backbone-scale', type=float, default=1.0,
                     help='Factor to scale DNA backbone in atomic model; try 0.25 to avoid clashes for atomistic simulations')
 
+parser.add_argument('--draw-cylinders', type=bool, default=False,
+                    help='whether or not to draw the cylinders')
+
 parser.add_argument('input_file', type=str,
                     help="""Any of the following:
                     (1) a cadnano JSON file;
@@ -54,19 +57,29 @@ if __name__ == '__main__':
 
     if args.output_prefix is not None:
         prefix = args.output_prefix
-
-    run_args = dict(
-        model = model,
-        output_name = prefix,
-        job_id = "job-" + prefix,
-        directory = args.directory,
-        gpu = args.gpu,
-        coarse_output_period = int(args.output_period),
-        fine_output_period = int(args.output_period),
-        coarse_steps = int(args.coarse_steps),
-        fine_steps = int(args.fine_steps),
-        backbone_scale = args.backbone_scale
-    )
+    
+    if args.draw_cylinders is not False:
+        model = model
+        model.vmd_cylinder_tcl()
+        output_name = prefix+ ".cylinders.tcl"
+        model._clear_beads()
+        model._generate_atomic_model(scale=args.backbone_scale)
+        model.write_atomic_ENM( prefix )
+        model.atomic_simulate( output_name = prefix )
+        
+    else:
+        run_args = dict(
+            model = model,
+            output_name = prefix,
+            job_id = "job-" + prefix,
+            directory = args.directory,
+            gpu = args.gpu,
+            coarse_output_period = int(args.output_period),
+            fine_output_period = int(args.output_period),
+            coarse_steps = int(args.coarse_steps),
+            fine_steps = int(args.fine_steps),
+            backbone_scale = args.backbone_scale
+        )
 
 
-    simulate( **run_args )
+        simulate( **run_args )
diff --git a/mrdna/segmentmodel.py b/mrdna/segmentmodel.py
index c8efc6a95105ffafdec639ceec9fee512aeac35d..fd78e543c6778d06a440c9949a1af7abff7d0eee 100644
--- a/mrdna/segmentmodel.py
+++ b/mrdna/segmentmodel.py
@@ -2813,3 +2813,25 @@ proc calcforces {} {
                 raise(Exception)
 
             file_handle.write("")
+
+    def vmd_cylinder_tcl(self, file_name="drawCylinders.tcl"):
+        #raise NotImplementedError
+        with open(file_name, 'w') as tclFile:
+            tclFile.write("## beginning TCL script \n")
+            def draw_cylinder(segment,radius_value=10,color="cyan"):
+                tclFile.write("## cylinder being drawn... \n")
+                r0 = segment.contour_to_position(0)
+                r1 = segment.contour_to_position(1)
+                
+                x0,y0,z0 = r0
+                x1,y1,z1 = r1
+                
+                radius_value = str(radius_value)
+                color = str(color)
+                
+                tclFile.write("graphics top color {} \n".format(color))
+                tclFile.write("graphics top cylinder {{ {} {} {} }} {{ {} {} {} }} radius {} resolution 30 filled yes \n".format(r0[0], r0[1], r0[2], r1[0], r1[1], r1[2], radius_value))
+
+            ## material
+            tclFile.write("graphics top materials on \n")
+            tclFile.write("graphics top material AOEdgy \n")