api/nblib/listed_forces/tests/linear_chain_input.hpp

   1 /*
   2  * This file is part of the GROMACS molecular simulation package.
   3  *
   4  * Copyright (c) 2020,2021, by the GROMACS development team, led by
   5  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
   6  * and including many others, as listed in the AUTHORS file in the
   7  * top-level source directory and at http://www.gromacs.org.
   8  *
   9  * GROMACS is free software; you can redistribute it and/or
  10  * modify it under the terms of the GNU Lesser General Public License
  11  * as published by the Free Software Foundation; either version 2.1
  12  * of the License, or (at your option) any later version.
  13  *
  14  * GROMACS is distributed in the hope that it will be useful,
  15  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  16  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  17  * Lesser General Public License for more details.
  18  *
  19  * You should have received a copy of the GNU Lesser General Public
  20  * License along with GROMACS; if not, see
  21  * http://www.gnu.org/licenses, or write to the Free Software Foundation,
  22  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
  23  *
  24  * If you want to redistribute modifications to GROMACS, please
  25  * consider that scientific software is very special. Version
  26  * control is crucial - bugs must be traceable. We will be happy to
  27  * consider code for inclusion in the official distribution, but
  28  * derived work must not be called official GROMACS. Details are found
  29  * in the README & COPYING files - if they are missing, get the
  30  * official version at http://www.gromacs.org.
  31  *
  32  * To help us fund GROMACS development, we humbly ask that you cite
  33  * the research papers on the package. Check out http://www.gromacs.org.
  34  */
  35 /*! \internal \file
  36  * \brief
  37  * This implements basic nblib utility tests
  38  *
  39  * \author Victor Holanda <victor.holanda@cscs.ch>
  40  * \author Joe Jordan <ejjordan@kth.se>
  41  * \author Prashanth Kanduri <kanduri@cscs.ch>
  42  * \author Sebastian Keller <keller@cscs.ch>
  43  */
  44
  45 #ifndef NBLIB_LINEAR_CHAIN_DATA_HPP
  46 #define NBLIB_LINEAR_CHAIN_DATA_HPP
  47
  48 #include "nblib/box.h"
  49 #include "nblib/topologyhelpers.h"
  50 #include "nblib/listed_forces/traits.h"
  51
  52 namespace nblib
  53 {
  54
  55 //! \brief sets up an interaction tuples for a linear chain with nParticles
  56 class LinearChainData
  57 {
  58 public:
  59     LinearChainData(int nP, float outlierRatio = 0.02) : nParticles(nP)
  60     {
  61         HarmonicBondType              bond1{ 376560, real(0.1001 * std::sqrt(3)) };
  62         HarmonicBondType              bond2{ 313800, real(0.1001 * std::sqrt(3)) };
  63         std::vector<HarmonicBondType> bonds{ bond1, bond2 };
  64         pickType<HarmonicBondType>(interactions).parameters = bonds;
  65
  66         HarmonicAngle              angle(Degrees(179.9), 397.5);
  67         std::vector<HarmonicAngle> angles{ angle };
  68         pickType<HarmonicAngle>(interactions).parameters = angles;
  69
  70         std::vector<InteractionIndex<HarmonicBondType>> bondIndices;
  71         for (int i = 0; i < nParticles - 1; ++i)
  72         {
  73             // third index: alternate between the two bond parameters
  74             bondIndices.push_back(InteractionIndex<HarmonicBondType>{ i, i + 1, i % 2 });
  75         }
  76         pickType<HarmonicBondType>(interactions).indices = bondIndices;
  77
  78         std::vector<InteractionIndex<HarmonicAngle>> angleIndices;
  79         for (int i = 0; i < nParticles - 2; ++i)
  80         {
  81             angleIndices.push_back(InteractionIndex<HarmonicAngle>{ i, i + 1, i + 2, 0 });
  82         }
  83         pickType<HarmonicAngle>(interactions).indices = angleIndices;
  84
  85         // initialize coordinates
  86         x.resize(nParticles);
  87         for (int i = 0; i < nParticles; ++i)
  88         {
  89             x[i] = real(i) * gmx::RVec{ 0.1, 0.1, 0.1 };
  90         }
  91
  92         forces = std::vector<gmx::RVec>(nParticles, gmx::RVec{ 0, 0, 0 });
  93
  94         box.reset(new Box(0.1 * (nParticles + 1), 0.1 * (nParticles + 1), 0.1 * (nParticles + 1)));
  95
  96         addOutliers(outlierRatio);
  97     }
  98
  99     // void createAggregates() { aggregateTransformations(interactions); }
 100
 101     void addOutliers(float outlierRatio)
 102     {
 103         HarmonicBondType dummyBond{ 1e-6, real(0.1001 * std::sqrt(3)) };
 104         pickType<HarmonicBondType>(interactions).parameters.push_back(dummyBond);
 105         int bondIndex = pickType<HarmonicBondType>(interactions).parameters.size() - 1;
 106
 107         int nOutliers = nParticles * outlierRatio;
 108         srand(42);
 109         for (int s = 0; s < nOutliers; ++s)
 110         {
 111             int from = rand() % nParticles;
 112             int to   = rand() % nParticles;
 113             if (from != to)
 114             {
 115                 pickType<HarmonicBondType>(interactions)
 116                         .indices.push_back(InteractionIndex<HarmonicBondType>{ from, to, bondIndex });
 117             }
 118         }
 119     }
 120
 121     int nParticles;
 122
 123     std::vector<gmx::RVec> x;
 124     std::vector<gmx::RVec> forces;
 125
 126     ListedInteractionData interactions;
 127
 128     std::shared_ptr<Box> box;
 129 };
 130
 131 } // namespace nblib
 132
 133 #endif // NBLIB_LINEAR_CHAIN_DATA_HPP
 134