Repository: choderalab/perses
Issue: Specifying timestep in FAH YAML leads to error

Manually specifying a timestep in the FAH YAML file:
```YAML
protein_pdb: protein.pdb
ligand_file: DOCKING_4MDS_POSES.sdf
old_ligand_index: 0
new_ligand_index: 1
forcefield_files:
    - amber/ff14SB.xml # ff14SB protein force field
    - amber/tip3p_standard.xml # TIP3P and recommended monovalent ion parameters
    - amber/tip3p_HFE_multivalent.xml # for divalent ions
    - amber/phosaa10.xml # HANDLES THE TPO
# don't constrain bonds to hydrogen to work around pathological slowdown with CCMA
remove_constraints: not water 
atom_expression:
    - IntType
bond_expression:
    - DefaultBonds
small_molecule_forcefield:
    - openff-1.2.0
phase_project_ids:
    complex: 13424
    solvent: 13425
    vacuum: vacuum-sprint3
ionic_strength: 0.070 # molar, to match experimental conditions
timestep: 2.0 # 2 fs
```
leads to this error:
```
2020-08-26 23:29:31,741:(0.00s):root:Setting up phase complex
Traceback (most recent call last):
  File "/home/chodera/miniconda/envs/perses/lib/python3.7/site-packages/simtk/unit/quantity.py", line 202, in __getattr__
    ret_val = getattr(self._value, attribute)
AttributeError: 'float' object has no attribute 'exp'

The above exception was the direct cause of the following exception:

Traceback (most recent call last):
  File "/home/chodera/miniconda/envs/perses/bin/perses-fah", line 33, in <module>
    sys.exit(load_entry_point('perses==0.8.1.dev0', 'console_scripts', 'perses-fah')())
  File "/home/chodera/miniconda/envs/perses/lib/python3.7/site-packages/perses-0.8.1.dev0-py3.7.egg/perses/app/fah_generator.py", line 469, in run
    run_neq_fah_setup(**setup_options)
  File "/home/chodera/miniconda/envs/perses/lib/python3.7/site-packages/perses-0.8.1.dev0-py3.7.egg/perses/app/fah_generator.py", line 374, in run_neq_fah_setup
    integrator = make_neq_integrator(**setup_options)
  File "/home/chodera/miniconda/envs/perses/lib/python3.7/site-packages/perses-0.8.1.dev0-py3.7.egg/perses/app/fah_generator.py", line 81, in make_neq_integrator
    integrator = PeriodicNonequilibriumIntegrator(alchemical_functions, nsteps_eq, nsteps_neq, neq_splitting, timestep=timestep)
  File "/home/chodera/miniconda/envs/perses/lib/python3.7/site-packages/openmmtools/integrators.py", line 1980, in __init__
    super().__init__(alchemical_functions, splitting, nsteps_neq=nsteps_neq, **kwargs)
  File "/home/chodera/miniconda/envs/perses/lib/python3.7/site-packages/openmmtools/integrators.py", line 1783, in __init__
    super(AlchemicalNonequilibriumLangevinIntegrator, self).__init__(*args, **kwargs)
  File "/home/chodera/miniconda/envs/perses/lib/python3.7/site-packages/openmmtools/integrators.py", line 1578, in __init__
    super(NonequilibriumLangevinIntegrator, self).__init__(*args, **kwargs)
  File "/home/chodera/miniconda/envs/perses/lib/python3.7/site-packages/openmmtools/integrators.py", line 1138, in __init__
    self.addGlobalVariable("a", np.exp(-gamma * h))
TypeError: loop of ufunc does not support argument 0 of type Quantity which has no callable exp method
```

Changed files:
perses/app/setup_relative_calculation.py, perses/app/fah_generator.py

diff --git a/perses/app/fah_generator.py b/perses/app/fah_generator.py
index ef4c246d..292c1625 100644
--- a/perses/app/fah_generator.py
+++ b/perses/app/fah_generator.py
@@ -189,7 +189,7 @@ def run_neq_fah_setup(ligand_file,
                       trajectory_directory,
                       complex_box_dimensions=(9.8, 9.8, 9.8),
                       solvent_box_dimensions=(3.5, 3.5, 3.5),
-                      timestep=4.0 * unit.femtosecond,
+                      timestep=4.0,
                       eq_splitting='V R O R V',
                       neq_splitting='V R H O R V',
                       measure_shadow_work=False,
@@ -257,7 +257,7 @@ def run_neq_fah_setup(ligand_file,
             define box dimensions of complex phase (in nm)
         solvent_box_dimensions : Vec3, default=(3.5, 3.5, 3.5)
             define box dimensions of solvent phase (in nm)
-        timestep : simtk.unit.Quantity, default=4.*unit.femtosecond
+        timestep : float, default=4.
             step size of nonequilibrium integration
         eq_splitting : str, default = 'V R O R V'
             splitting string of relaxation dynamics
diff --git a/perses/app/setup_relative_calculation.py b/perses/app/setup_relative_calculation.py
index a9460041..8aa1b81f 100644
--- a/perses/app/setup_relative_calculation.py
+++ b/perses/app/setup_relative_calculation.py
@@ -316,7 +316,7 @@ def run_setup(setup_options, serialize_systems=True, build_samplers=True):
         - 'topology_proposals':
     """
     phases = setup_options['phases']
-    known_phases = ['complex','solvent','vacuum']
+    known_phases = ['complex', 'solvent', 'vacuum']
     for phase in phases:
         assert (phase in known_phases), f"Unknown phase, {phase} provided. run_setup() can be used with {known_phases}"
 
@@ -352,7 +352,10 @@ def run_setup(setup_options, serialize_systems=True, build_samplers=True):
     forcefield_files = setup_options['forcefield_files']
 
     if "timestep" in setup_options:
-        timestep = setup_options['timestep'] * unit.femtoseconds
+        if isinstance(setup_options['timestep'], float):
+            timestep = setup_options['timestep'] * unit.femtoseconds
+        else:
+            timestep = setup_options['timestep']
         _logger.info(f"\\ttimestep: {timestep}.")
     else:
         timestep = 1.0 * unit.femtoseconds
@@ -380,10 +383,18 @@ def run_setup(setup_options, serialize_systems=True, build_samplers=True):
     else:
         measure_shadow_work = False
         _logger.info(f"\\tno measure_shadow_work specified: defaulting to False.")
-
-    pressure = setup_options['pressure'] * unit.atmosphere
-    temperature = setup_options['temperature'] * unit.kelvin
-    solvent_padding_angstroms = setup_options['solvent_padding'] * unit.angstrom
+    if isinstance(setup_options['pressure'],float):
+        pressure = setup_options['pressure'] * unit.atmosphere
+    else:
+        pressure = setup_options['pressure']
+    if isinstance(setup_options['temperature'], float):
+        temperature = setup_options['temperature'] * unit.kelvin
+    else:
+        temperature = setup_options['temperature']
+    if isinstance(setup_options['solvent_padding'], float):
+        solvent_padding_angstroms = setup_options['solvent_padding'] * unit.angstrom
+    else:
+        solvent_padding_angstroms = setup_options['solvent_padding']
     if isinstance(setup_options['ionic_strength'], float):
         ionic_strength = setup_options['ionic_strength'] * unit.molar
     else: