From fb248cd11c751b4a5bc7791e8b608ad5d5de3fef Mon Sep 17 00:00:00 2001 From: Christian Blau Date: Mon, 20 Aug 2018 11:07:26 +0200 Subject: [PATCH] Add test for structural similarity measures rmsd and rho. Tests rmsd and size-independent rho-factor calculation routines. Change-Id: I1c74cb4e3185e38aafc5cb04352af2789b40eca9 --- src/gromacs/math/tests/CMakeLists.txt | 3 +- src/gromacs/math/tests/dofit.cpp | 101 ++++++++++++++++++++++++++++++++++ 2 files changed, 103 insertions(+), 1 deletion(-) create mode 100644 src/gromacs/math/tests/dofit.cpp diff --git a/src/gromacs/math/tests/CMakeLists.txt b/src/gromacs/math/tests/CMakeLists.txt index 781248c023..537c2b60c9 100644 --- a/src/gromacs/math/tests/CMakeLists.txt +++ b/src/gromacs/math/tests/CMakeLists.txt @@ -1,7 +1,7 @@ # # This file is part of the GROMACS molecular simulation package. # -# Copyright (c) 2014,2015, by the GROMACS development team, led by +# Copyright (c) 2014,2015,2018, by the GROMACS development team, led by # Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, # and including many others, as listed in the AUTHORS file in the # top-level source directory and at http://www.gromacs.org. @@ -33,6 +33,7 @@ # the research papers on the package. Check out http://www.gromacs.org. gmx_add_unit_test(MathUnitTests math-test + dofit.cpp functions.cpp invertmatrix.cpp vectypes.cpp diff --git a/src/gromacs/math/tests/dofit.cpp b/src/gromacs/math/tests/dofit.cpp new file mode 100644 index 0000000000..c971c5bb00 --- /dev/null +++ b/src/gromacs/math/tests/dofit.cpp @@ -0,0 +1,101 @@ +/* + * This file is part of the GROMACS molecular simulation package. + * + * Copyright (c) 2018, by the GROMACS development team, led by + * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, + * and including many others, as listed in the AUTHORS file in the + * top-level source directory and at http://www.gromacs.org. + * + * GROMACS is free software; you can redistribute it and/or + * modify it under the terms of the GNU Lesser General Public License + * as published by the Free Software Foundation; either version 2.1 + * of the License, or (at your option) any later version. + * + * GROMACS is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + * Lesser General Public License for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with GROMACS; if not, see + * http://www.gnu.org/licenses, or write to the Free Software Foundation, + * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. + * + * If you want to redistribute modifications to GROMACS, please + * consider that scientific software is very special. Version + * control is crucial - bugs must be traceable. We will be happy to + * consider code for inclusion in the official distribution, but + * derived work must not be called official GROMACS. Details are found + * in the README & COPYING files - if they are missing, get the + * official version at http://www.gromacs.org. + * + * To help us fund GROMACS development, we humbly ask that you cite + * the research papers on the package. Check out http://www.gromacs.org. + */ +/*! \internal \file + * \brief + * Tests structure similarity measures rmsd and size-independent rho factor. + * + * \author Christian Blau + * \ingroup module_math + */ +#include "gmxpre.h" + +#include + +#include + +#include "gromacs/math/do_fit.h" +#include "gromacs/math/vec.h" + +#include "testutils/testasserts.h" + +namespace +{ + +using gmx::RVec; +using gmx::test::defaultRealTolerance; +class StructureSimilarityTest : public ::testing::Test +{ + protected: + static constexpr int c_nAtoms = 4; + std::array structureA_ {{{1, 0, 0}, {0, 1, 0}, {0, 0, 1}, {0, 0, 0}}}; + std::array structureB_ {{{0, 1, 0}, {0, 0, 1}, {1, 0, 0}, {0, 0, 0}}}; + std::array masses_ {{1, 1, 1, 0}}; + std::array index_ {{0, 1, 2}}; + rvec * x1_ = gmx::as_rvec_array(structureA_.data()); + rvec * x2_ = gmx::as_rvec_array(structureB_.data()); + real * m_ = masses_.data(); +}; + +TEST_F(StructureSimilarityTest, StructureComparedToSelfHasZeroRMSD) +{ + EXPECT_REAL_EQ_TOL(0., rmsdev(c_nAtoms, m_, x1_, x1_), defaultRealTolerance()); +} + +TEST_F(StructureSimilarityTest, StructureComparedToSelfHasZeroRho) +{ + EXPECT_REAL_EQ_TOL(0., rhodev(c_nAtoms, m_, x1_, x1_), defaultRealTolerance()); +} + +TEST_F(StructureSimilarityTest, YieldsCorrectRMSD) +{ + EXPECT_REAL_EQ_TOL(sqrt(2.0), rmsdev(c_nAtoms, m_, x1_, x2_), defaultRealTolerance()); +} + +TEST_F(StructureSimilarityTest, YieldsCorrectRho) +{ + EXPECT_REAL_EQ_TOL(2., rhodev(c_nAtoms, m_, x1_, x2_), defaultRealTolerance()); +} + +TEST_F(StructureSimilarityTest, YieldsCorrectRMSDWithIndex) +{ + EXPECT_REAL_EQ_TOL(sqrt(2.0), rmsdev_ind(index_.size(), index_.data(), m_, x1_, x2_), defaultRealTolerance()); +} + +TEST_F(StructureSimilarityTest, YieldsCorrectRhoWidthIndex) +{ + EXPECT_REAL_EQ_TOL(2., rhodev_ind(index_.size(), index_.data(), m_, x1_, x2_), defaultRealTolerance()); +} + +} // namespace -- 2.11.4.GIT