diff --git a/atomistics/calculators/lammps/calculator.py b/atomistics/calculators/lammps/calculator.py
index 481a5f18..dc968870 100644
--- a/atomistics/calculators/lammps/calculator.py
+++ b/atomistics/calculators/lammps/calculator.py
@@ -31,7 +31,7 @@
     from atomistics.calculators.interface import TaskName
 
 
-def calc_energy_with_lammps(structure: Atoms, potential_dataframe: DataFrame, lmp=None):
+def calc_energy_with_lammps(structure: Atoms, potential_dataframe: DataFrame, lmp=None, **kwargs):
     template_str = LAMMPS_THERMO_STYLE + "\n" + LAMMPS_THERMO + "\n" + LAMMPS_RUN
     lmp_instance = lammps_run(
         structure=structure,
@@ -42,13 +42,14 @@ def calc_energy_with_lammps(structure: Atoms, potential_dataframe: DataFrame, lm
             run=0,
         ),
         lmp=lmp,
+        **kwargs,
     )
     energy_pot = lmp_instance.interactive_energy_pot_getter()
     lammps_shutdown(lmp_instance=lmp_instance, close_instance=lmp is None)
     return energy_pot
 
 
-def calc_forces_with_lammps(structure: Atoms, potential_dataframe: DataFrame, lmp=None):
+def calc_forces_with_lammps(structure: Atoms, potential_dataframe: DataFrame, lmp=None, **kwargs):
     template_str = LAMMPS_THERMO_STYLE + "\n" + LAMMPS_THERMO + "\n" + LAMMPS_RUN
     lmp_instance = lammps_run(
         structure=structure,
@@ -59,6 +60,7 @@ def calc_forces_with_lammps(structure: Atoms, potential_dataframe: DataFrame, lm
             run=0,
         ),
         lmp=lmp,
+        **kwargs,
     )
     forces = lmp_instance.interactive_forces_getter()
     lammps_shutdown(lmp_instance=lmp_instance, close_instance=lmp is None)
@@ -66,7 +68,7 @@ def calc_forces_with_lammps(structure: Atoms, potential_dataframe: DataFrame, lm
 
 
 def calc_energy_and_forces_with_lammps(
-    structure, potential_dataframe: DataFrame, lmp=None
+    structure, potential_dataframe: DataFrame, lmp=None, **kwargs
 ):
     template_str = LAMMPS_THERMO_STYLE + "\n" + LAMMPS_THERMO + "\n" + LAMMPS_RUN
     lmp_instance = lammps_run(
@@ -78,6 +80,7 @@ def calc_energy_and_forces_with_lammps(
             run=0,
         ),
         lmp=lmp,
+        **kwargs,
     )
     energy_pot = lmp_instance.interactive_energy_pot_getter()
     forces = lmp_instance.interactive_forces_getter()
@@ -95,6 +98,7 @@ def optimize_positions_and_volume_with_lammps(
     maxeval=10000000,
     thermo=10,
     lmp=None,
+    **kwargs
 ):
     template_str = (
         LAMMPS_MINIMIZE_VOLUME
@@ -118,6 +122,7 @@ def optimize_positions_and_volume_with_lammps(
             thermo=thermo,
         ),
         lmp=lmp,
+        **kwargs,
     )
     structure_copy = structure.copy()
     structure_copy.set_cell(lmp_instance.interactive_cells_getter(), scale_atoms=True)
@@ -136,6 +141,7 @@ def optimize_positions_with_lammps(
     maxeval=10000000,
     thermo=10,
     lmp=None,
+    **kwargs
 ):
     template_str = LAMMPS_THERMO_STYLE + "\n" + LAMMPS_THERMO + "\n" + LAMMPS_MINIMIZE
     lmp_instance = lammps_run(
@@ -151,6 +157,7 @@ def optimize_positions_with_lammps(
             thermo=thermo,
         ),
         lmp=lmp,
+        **kwargs,
     )
     structure_copy = structure.copy()
     structure_copy.positions = lmp_instance.interactive_positions_getter()
@@ -159,7 +166,22 @@ def optimize_positions_with_lammps(
 
 
 def calc_molecular_dynamics_thermal_expansion_with_lammps(
-    structure, potential_dataframe, Tstart=15, Tstop=1500, Tstep=5, lmp=None
+    structure,
+    potential_dataframe,
+    Tstart=15,
+    Tstop=1500,
+    Tstep=5,
+    Tdamp=0.1,
+    run=100,
+    thermo=100,
+    timestep=0.001,
+    Pstart=0.0,
+    Pstop=0.0,
+    Pdamp=1.0,
+    seed=4928459,
+    dist="gaussian",
+    lmp=None,
+    **kwargs
 ):
     init_str = (
         LAMMPS_THERMO_STYLE
@@ -179,7 +201,17 @@ def calc_molecular_dynamics_thermal_expansion_with_lammps(
         init_str=init_str,
         run_str=run_str,
         temperature_lst=temperature_lst,
+        run=run,
+        thermo=thermo,
+        timestep=timestep,
+        Tdamp=Tdamp,
+        Pstart=Pstart,
+        Pstop=Pstop,
+        Pdamp=Pdamp,
+        seed=seed,
+        dist=dist,
         lmp=lmp,
+        **kwargs,
     )
     return temperature_lst, volume_md_lst
 
diff --git a/atomistics/calculators/lammps/helpers.py b/atomistics/calculators/lammps/helpers.py
index dfa13174..e13c67eb 100644
--- a/atomistics/calculators/lammps/helpers.py
+++ b/atomistics/calculators/lammps/helpers.py
@@ -36,12 +36,14 @@ def template_render_run(
     )
 
 
-def lammps_run(structure, potential_dataframe, input_template, lmp=None):
+def lammps_run(structure, potential_dataframe, input_template, lmp=None, **kwargs):
     potential_dataframe = validate_potential_dataframe(
         potential_dataframe=potential_dataframe
     )
     if lmp is None:
-        lmp = LammpsASELibrary()
+        lmp = LammpsASELibrary(
+            **kwargs
+        )
 
     # write structure to LAMMPS
     lmp.interactive_structure_setter(
@@ -65,31 +67,47 @@ def lammps_run(structure, potential_dataframe, input_template, lmp=None):
 
 
 def lammps_thermal_expansion_loop(
-    structure, potential_dataframe, init_str, run_str, temperature_lst, lmp=None
+    structure,
+    potential_dataframe,
+    init_str,
+    run_str,
+    temperature_lst,
+    run=100,
+    thermo=100,
+    timestep=0.001,
+    Tdamp=0.1,
+    Pstart=0.0,
+    Pstop=0.0,
+    Pdamp=1.0,
+    seed=4928459,
+    dist="gaussian",
+    lmp=None,
+    **kwargs,
 ):
     lmp_instance = lammps_run(
         structure=structure,
         potential_dataframe=potential_dataframe,
         input_template=Template(init_str).render(
-            thermo=100,
+            thermo=thermo,
             temp=temperature_lst[0],
-            timestep=0.001,
-            seed=4928459,
-            dist="gaussian",
+            timestep=timestep,
+            seed=seed,
+            dist=dist,
         ),
         lmp=lmp,
+        **kwargs,
     )
 
     volume_md_lst = []
     for temp in temperature_lst:
         run_str_rendered = Template(run_str).render(
-            run=100,
+            run=run,
             Tstart=temp - 5,
             Tstop=temp,
-            Tdamp=0.1,
-            Pstart=0.0,
-            Pstop=0.0,
-            Pdamp=1.0,
+            Tdamp=Tdamp,
+            Pstart=Pstart,
+            Pstop=Pstop,
+            Pdamp=Pdamp,
         )
         for l in run_str_rendered.split("\n"):
             lmp_instance.interactive_lib_command(l)