diff --git a/mrdna/model/arbdmodel.py b/mrdna/model/arbdmodel.py
index cb85a16dabbdf20a054dc38e55a4788fdbdace88..f10babeb5b2fa71a37ed1145dd67d5858b718ff5 100644
--- a/mrdna/model/arbdmodel.py
+++ b/mrdna/model/arbdmodel.py
@@ -684,7 +684,7 @@ class ArbdModel(PdbModel):
         if typeA != typeB:
             self.nbSchemes.append( (nbScheme, typeB, typeA) )
 
-    def simulate(self, output_name, output_directory='output', num_steps=100000000, timestep=None, gpu=0, output_period=1e4, arbd=None, directory='.'):
+    def simulate(self, output_name, output_directory='output', num_steps=100000000, timestep=None, gpu=0, output_period=1e4, arbd=None, directory='.', replicas=1):
         assert(type(gpu) is int)
         num_steps = int(num_steps)
 
@@ -732,7 +732,11 @@ class ArbdModel(PdbModel):
             os.sync()
 
             ## http://stackoverflow.com/questions/18421757/live-output-from-subprocess-command
-            cmd = (arbd, '-g', "%d" % gpu, "%s.bd" % output_name, "%s/%s" % (output_directory, output_name))
+
+            cmd = [arbd, '-g', "%d" % gpu]
+            if replicas > 1:
+                cmd = cmd + ['-r',replicas]
+            cmd = cmd + ["%s.bd" % output_name, "%s/%s" % (output_directory, output_name)]
             cmd = tuple(str(x) for x in cmd)
 
             print("Running ARBD with: %s" % " ".join(cmd))