src/gromacs/analysisdata/modules/frameaverager.h

   1 /*
   2  * This file is part of the GROMACS molecular simulation package.
   3  *
   4  * Copyright (c) 2013,2014,2015,2019, 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  * Declares gmx::AnalysisDataFrameAverager.
  38  *
  39  * \author Teemu Murtola <teemu.murtola@gmail.com>
  40  * \ingroup module_analysisdata
  41  */
  42 #ifndef GMX_ANALYSISDATA_MODULES_FRAMEAVERAGER_H
  43 #define GMX_ANALYSISDATA_MODULES_FRAMEAVERAGER_H
  44
  45 #include <vector>
  46
  47 #include "gromacs/utility/gmxassert.h"
  48 #include "gromacs/utility/real.h"
  49
  50 namespace gmx
  51 {
  52
  53 class AnalysisDataPointSetRef;
  54
  55 /*! \internal
  56  * \brief
  57  * Helper class for modules that average values over frames.
  58  *
  59  * This class implements common functionality for analysis data modules that
  60  * need to average a set of values over frames.  Currently, it is designed for
  61  * computing averages for each input column independently, but should be
  62  * relatively easy to make more general if required.
  63  *
  64  * This class takes care of accumulating the values and computing their
  65  * variance.  It allows different number of samples for each input column.
  66  * Accumulation is always in double precision and uses a formula that is
  67  * relatively stable numerically.  For now, does nothing fancy,
  68  * but provides ground for other implementation (e.g., related to
  69  * parallelization) that would benefit all such modules.
  70  *
  71  * Methods in this class do not throw unless otherwise indicated.
  72  *
  73  * \ingroup module_analysisdata
  74  */
  75 class AnalysisDataFrameAverager
  76 {
  77 public:
  78     AnalysisDataFrameAverager() : bFinished_(false) {}
  79
  80     /*! \brief
  81      * Returns the number of columns in this averager.
  82      */
  83     int columnCount() const { return values_.size(); }
  84
  85     /*! \brief
  86      * Sets the number of columns in the input data.
  87      *
  88      * \throws std::bad_alloc if out of memory.
  89      *
  90      * Typically called from IAnalysisDataModule::dataStarted().
  91      *
  92      * Must be called exactly once, before setting calling any other method
  93      * in the class.
  94      */
  95     void setColumnCount(int columnCount);
  96     /*! \brief
  97      * Adds a single value to the average for a given column.
  98      *
  99      * \param[in] index  Index of the column to add the value to.
 100      * \param[in] value  Value to add to the sample.
 101      */
 102     void addValue(int index, real value);
 103     /*! \brief
 104      * Accumulates data from a given point set into the average.
 105      *
 106      * Typically called from IAnalysisDataModule::pointsAdded().
 107      *
 108      * Each call accumulates the values for those columns that are present
 109      * in the point set.  Can be called multiple times for a frame, and
 110      * does not need to be called for every frame.
 111      */
 112     void addPoints(const AnalysisDataPointSetRef& points);
 113     /*! \brief
 114      * Finalizes the calculation of the averages and variances.
 115      *
 116      * Does any computation that is not done during the accumulation in
 117      * addPoints().  Currently, does nothing, but provided as a placeholder
 118      * for more complex implementation.
 119      *
 120      * Typically called from IAnalysisDataModule::dataFinished().
 121      */
 122     void finish();
 123
 124     /*! \brief
 125      * Returns the computed average for a given column.
 126      *
 127      * If called before finish(), the results are undefined.
 128      */
 129     real average(int index) const
 130     {
 131         GMX_ASSERT(index >= 0 && index < columnCount(), "Invalid column index");
 132         GMX_ASSERT(bFinished_, "Values available only after finished() has been called");
 133         return values_[index].average;
 134     }
 135     /*! \brief
 136      * Returns the computed (sample) variance for a given column.
 137      *
 138      * If called before finish(), the results are undefined.
 139      */
 140     real variance(int index) const
 141     {
 142         GMX_ASSERT(index >= 0 && index < columnCount(), "Invalid column index");
 143         GMX_ASSERT(bFinished_, "Values available only after finished() has been called");
 144         const AverageItem& item = values_[index];
 145         return item.samples > 1 ? item.squaredSum / (item.samples - 1) : 0.0;
 146     }
 147     /*! \brief
 148      * Returns the number of samples for a given column.
 149      *
 150      * If called before finish(), the results are undefined.
 151      */
 152     int sampleCount(int index) const
 153     {
 154         GMX_ASSERT(index >= 0 && index <= ssize(values_), "Invalid column index");
 155         GMX_ASSERT(bFinished_, "Values available only after finished() has been called");
 156         return values_[index].samples;
 157     }
 158
 159 private:
 160     struct AverageItem
 161     {
 162         AverageItem() : average(0.0), squaredSum(0.0), samples(0) {}
 163
 164         //! Average of the values so far.
 165         double average;
 166         //! Sum of squared deviations from the average for values so far.
 167         double squaredSum;
 168         //! Number of values so far.
 169         int samples;
 170     };
 171
 172     std::vector<AverageItem> values_;
 173     bool                     bFinished_;
 174 };
 175
 176 } // namespace gmx
 177
 178 #endif