Merge pull request #1428 from glotzerlab/remove-inaccurate-validation…

…-tests

Remove inaccurate validation tests.

README.md

@@ -30,6 +30,8 @@ types of particle simulations.
   Additional information and publications.
 - [HOOMD-blue benchmark scripts](https://github.com/glotzerlab/hoomd-benchmarks):
   Scripts to evaluate the performance of HOOMD-blue simulations.
+- [HOOMD-blue validation tests](https://github.com/glotzerlab/hoomd-validation):
+  Scripts to validate that HOOMD-blue performs accurate simulations.
 
 ## Related tools
 

hoomd/conftest.py

@@ -14,8 +14,6 @@
 import atexit
 import os
 import numpy
-import math
-import warnings
 from hoomd.logging import LoggerCategories
 from hoomd.snapshot import Snapshot
 from hoomd import Simulation
@@ -548,108 +546,6 @@ def autotuned_kernel_parameter_check(instance, activate, all_optional=False):
         assert instance.kernel_parameters == initial_kernel_parameters
 
 
-class BlockAverage:
-    """Block average method for estimating standard deviation of the mean.
-
-    Args:
-        data: List of values
-    """
-
-    def __init__(self, data):
-        # round down to the nearest power of 2
-        N = 2**int(math.log(len(data)) / math.log(2))
-        if N != len(data):
-            warnings.warn(
-                "Ignoring some data. Data array should be a power of 2.")
-
-        block_sizes = []
-        block_mean = []
-        block_variance = []
-
-        # take means of blocks and the mean/variance of all blocks, growing
-        # blocks by factors of 2
-        block_size = 1
-        while block_size <= N // 8:
-            num_blocks = N // block_size
-            block_data = numpy.zeros(num_blocks)
-
-            for i in range(0, num_blocks):
-                start = i * block_size
-                end = start + block_size
-                block_data[i] = numpy.mean(data[start:end])
-
-            block_mean.append(numpy.mean(block_data))
-            block_variance.append(numpy.var(block_data) / (num_blocks - 1))
-
-            block_sizes.append(block_size)
-            block_size *= 2
-
-        self._block_mean = numpy.array(block_mean)
-        self._block_variance = numpy.array(block_variance)
-        self._block_sizes = numpy.array(block_sizes)
-        self.data = numpy.array(data)
-
-        # check for a plateau in the relative error before the last data point
-        block_relative_error = numpy.sqrt(self._block_variance) / numpy.fabs(
-            self._block_mean)
-        relative_error_derivative = (numpy.diff(block_relative_error)
-                                     / numpy.diff(self._block_sizes))
-        if numpy.all(relative_error_derivative > 0):
-            warnings.warn("Block averaging failed to plateau, run longer")
-
-    def get_hierarchical_errors(self):
-        """Get details on the hierarchical errors."""
-        return (self._block_sizes, self._block_mean, self._block_variance)
-
-    @property
-    def standard_deviation(self):
-        """float: The error estimate on the mean."""
-        if numpy.all(self.data == self.data[0]):
-            return 0
-
-        return numpy.sqrt(numpy.max(self._block_variance))
-
-    @property
-    def mean(self):
-        """float: The mean."""
-        return self._block_mean[-1]
-
-    @property
-    def relative_error(self):
-        """float: The relative error."""
-        return self.standard_deviation / numpy.fabs(self.mean)
-
-    def assert_close(self,
-                     reference_mean,
-                     reference_deviation,
-                     z=6,
-                     max_relative_error=0.02):
-        """Assert that the distribution is constent with a given reference.
-
-        Also assert that the relative error of the distribution is small.
-        Otherwise, test runs with massive fluctuations would likely lead to
-        passing tests.
-
-        Args:
-            reference_mean: Known good mean value
-            reference_deviation: Standard deviation of the known good value
-            z: Number of standard deviations
-            max_relative_error: Maximum relative error to allow
-        """
-        sample_mean = self.mean
-        sample_deviation = self.standard_deviation
-
-        assert sample_deviation / sample_mean <= max_relative_error
-
-        # compare if 0 is within the confidence interval around the difference
-        # of the means
-        deviation_diff = ((sample_deviation**2
-                           + reference_deviation**2)**(1 / 2.))
-        mean_diff = math.fabs(sample_mean - reference_mean)
-        deviation_allowed = z * deviation_diff
-        assert mean_diff <= deviation_allowed
-
-
 class ListWriter(hoomd.custom.Action):
     """Log a single quantity to a list.
 

hoomd/hpmc/pytest/CMakeLists.txt

@@ -20,7 +20,6 @@ set(files __init__.py
           test_small_box_2d.py
           test_small_box_3d.py
           conftest.py
-          test_nec.py
           depletants_sphere.py
     )
 

hoomd/md/pytest/CMakeLists.txt

@@ -19,7 +19,6 @@ set(files __init__.py
     test_dihedral.py
     test_flags.py
     test_kernel_parameters.py
-    test_lj_equation_of_state.py
     test_potential.py
     test_pppm_coulomb.py
     test_manifolds.py

hoomd/md/pytest/test_alj.py

@@ -2,133 +2,13 @@
 # Part of HOOMD-blue, released under the BSD 3-Clause License.
 
 import numpy as np
-import numpy.testing as npt
 import pytest
 
 import hoomd
 import hoomd.conftest
 from hoomd import md
 
 
-@pytest.mark.validate
-def test_conservation(simulation_factory, lattice_snapshot_factory):
-    # For test, use a unit area hexagon.
-    particle_vertices = np.array([[6.20403239e-01, 0.00000000e+00, 0],
-                                  [3.10201620e-01, 5.37284966e-01, 0],
-                                  [-3.10201620e-01, 5.37284966e-01, 0],
-                                  [-6.20403239e-01, 7.59774841e-17, 0],
-                                  [-3.10201620e-01, -5.37284966e-01, 0],
-                                  [3.10201620e-01, -5.37284966e-01, 0]])
-    area = 1.0
-    circumcircle_radius = 0.6204032392788702
-    incircle_radius = 0.5372849659264116
-    num_vertices = len(particle_vertices)
-
-    circumcircle_diameter = 2 * circumcircle_radius
-
-    # Just initialize in a simple cubic lattice.
-    sim = simulation_factory(
-        lattice_snapshot_factory(a=4 * circumcircle_diameter,
-                                 n=10,
-                                 dimensions=2))
-    sim.seed = 175
-
-    # Initialize moments of inertia since original simulation was HPMC.
-    mass = area
-    # https://math.stackexchange.com/questions/2004798/moment-of-inertia-for-a-n-sided-regular-polygon # noqa
-    moment_inertia = (mass * circumcircle_diameter**2 / 6)
-    moment_inertia *= (1 + 2 * np.cos(np.pi / num_vertices)**2)
-
-    with sim.state.cpu_local_snapshot as snapshot:
-        snapshot.particles.mass[:] = mass
-        snapshot.particles.moment_inertia[:] = np.array([0, 0, moment_inertia])
-        # Not sure if this should be incircle or circumcircle;
-        # probably doesn't matter based on current usage, but may
-        # matter in the future for the potential if it's modified
-        # to actually use diameter.
-        snapshot.particles.diameter[:] = circumcircle_diameter
-    kT = 0.3
-    sim.state.thermalize_particle_momenta(hoomd.filter.All(), kT)
-
-    # Create box resize updater
-    packing_fraction = 0.4
-    final_area = area * sim.state.N_particles / packing_fraction
-    L_final = np.sqrt(final_area)
-    final_box = hoomd.Box.square(L_final)
-
-    n_compression_start = int(1e4)
-    n_compression_end = int(1e5)
-    n_compression_total = n_compression_end - n_compression_start
-
-    box_resize = hoomd.update.BoxResize(
-        box1=sim.state.box,
-        box2=final_box,
-        trigger=int(n_compression_total / 10000),
-        variant=hoomd.variant.Ramp(0, 1, n_compression_start,
-                                   n_compression_total),
-        filter=hoomd.filter.All())
-    sim.operations += box_resize
-
-    # Define forces and methods
-    r_cut_scale = 1.3
-    kernel_scale = (1 / np.cos(np.pi / num_vertices))
-    incircle_diameter = 2 * incircle_radius
-    r_cut_set = incircle_diameter * kernel_scale * r_cut_scale
-
-    alj = md.pair.aniso.ALJ(default_r_cut=r_cut_set, nlist=md.nlist.Cell(0.4))
-
-    alj.shape["A"] = {
-        "vertices": particle_vertices,
-        "faces": [],
-        "rounding_radii": 0
-    }
-
-    alpha = 0  # Make it WCA-only (no attraction)
-    eps_att = 1.0
-    alj.params[("A", "A")] = {
-        "epsilon": eps_att,
-        "sigma_i": incircle_diameter,
-        "sigma_j": incircle_diameter,
-        "alpha": alpha
-    }
-
-    nve = md.methods.NVE(filter=hoomd.filter.All())
-    integrator = md.Integrator(dt=1e-4,
-                               forces=[alj],
-                               methods=[nve],
-                               integrate_rotational_dof=True)
-    sim.operations.integrator = integrator
-
-    # Compress box
-    sim.run(n_compression_end)
-
-    thermo = md.compute.ThermodynamicQuantities(hoomd.filter.All())
-    sim.operations += thermo
-
-    # Reset velocities after the compression, and equilibriate
-    sim.state.thermalize_particle_momenta(hoomd.filter.All(), kT)
-    sim.run(1000)
-
-    # run sim and get values back
-    w = hoomd.conftest.ManyListWriter([(thermo, 'potential_energy'),
-                                       (thermo, 'kinetic_energy'),
-                                       (integrator, 'linear_momentum')])
-    writer = hoomd.write.CustomWriter(action=w, trigger=1)
-    sim.operations.writers.append(writer)
-    sim.run(1000)
-    pe, ke, momentum = w.data
-    total_energies = np.array(pe) + np.array(ke)
-
-    # Ensure energy conservation up to the 3 digit per-particle.
-    npt.assert_allclose(total_energies,
-                        total_energies[0],
-                        atol=0.03 * sim.state.N_particles)
-
-    # Test momentum conservation.
-    p_magnitude = np.linalg.norm(momentum, axis=-1)
-    npt.assert_allclose(p_magnitude, p_magnitude[0], atol=1e-13)
-
-
 def test_type_shapes(simulation_factory, two_particle_snapshot_factory):
     alj = md.pair.aniso.ALJ(md.nlist.Cell(buffer=0.1))
     sim = simulation_factory(two_particle_snapshot_factory(d=2.0))