Add characteristic polynomial for viscosity class

* Format `viscosity.py` with Black

transport_analysis/tests/test_viscosity.py

@@ -4,6 +4,9 @@
     ViscosityHelfand as VH,
 )
 import MDAnalysis as mda
+from MDAnalysis.units import constants
+from MDAnalysis.transformations import set_dimensions
+import numpy as np
 
 from MDAnalysis.exceptions import NoDataError
 from MDAnalysisTests.datafiles import PRM_NCBOX, TRJ_NCBOX, PSF, DCD
@@ -71,10 +74,42 @@ def step_vtraj_full(NSTEP):
         set_dimensions(dim)(u.trajectory.ts)
 
     # mass of 16.0
-    u.add_TopologyAttr('masses', [16.0])
+    u.add_TopologyAttr("masses", [16.0])
     return u
 
 
+def characteristic_poly_helfand(
+    total_frames, n_dim, temp_avg=300.0, vol_avg=8.0
+):
+    result = np.zeros(total_frames)
+
+    for lag in range(total_frames):
+        sum = 0
+        sum = np.dtype("float64").type(sum)
+
+        for curr in range((total_frames - lag)):
+            # mass * velocities * positions
+            helf_diff = 16.0 * (curr + lag) * 1 / 2 * (
+                (curr + lag) ** 2
+            ) - 16.0 * curr * 1 / 2 * (curr**2)
+            sum += helf_diff**2
+
+        vis_helf = (
+            sum
+            * n_dim
+            / (
+                (total_frames - lag)
+                * 2
+                * constants["Boltzmann_constant"]
+                * vol_avg
+                * temp_avg
+            )
+        )
+
+        result[lag] = vis_helf
+    return result
+
+
 class TestViscosityHelfand:
     # fixtures are helpful functions that set up a test
     # See more at https://docs.pytest.org/en/stable/how-to/fixtures.html

transport_analysis/viscosity.py

@@ -196,10 +196,7 @@ def _conclude(self):
 
         # divide by 2, boltzman constant, vol_avg, and temp_avg
         self.results.visc_by_particle = self.results.visc_by_particle / (
-            2
-            * constants["Boltzmann_constant"]
-            * self._vol_avg
-            * self.temp_avg
+            2 * constants["Boltzmann_constant"] * self._vol_avg * self.temp_avg
         )
         # average over # particles and update results array
         self.results.timeseries = self.results.visc_by_particle.mean(axis=1)