src/gromacs/trajectoryanalysis/analysismodule.cpp

   1 /*
   2  * This file is part of the GROMACS molecular simulation package.
   3  *
   4  * Copyright (c) 2010-2018, The GROMACS development team.
   5  * Copyright (c) 2019,2020,2021, by the GROMACS development team, led by
   6  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
   7  * and including many others, as listed in the AUTHORS file in the
   8  * top-level source directory and at http://www.gromacs.org.
   9  *
  10  * GROMACS is free software; you can redistribute it and/or
  11  * modify it under the terms of the GNU Lesser General Public License
  12  * as published by the Free Software Foundation; either version 2.1
  13  * of the License, or (at your option) any later version.
  14  *
  15  * GROMACS is distributed in the hope that it will be useful,
  16  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  17  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  18  * Lesser General Public License for more details.
  19  *
  20  * You should have received a copy of the GNU Lesser General Public
  21  * License along with GROMACS; if not, see
  22  * http://www.gnu.org/licenses, or write to the Free Software Foundation,
  23  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
  24  *
  25  * If you want to redistribute modifications to GROMACS, please
  26  * consider that scientific software is very special. Version
  27  * control is crucial - bugs must be traceable. We will be happy to
  28  * consider code for inclusion in the official distribution, but
  29  * derived work must not be called official GROMACS. Details are found
  30  * in the README & COPYING files - if they are missing, get the
  31  * official version at http://www.gromacs.org.
  32  *
  33  * To help us fund GROMACS development, we humbly ask that you cite
  34  * the research papers on the package. Check out http://www.gromacs.org.
  35  */
  36 /*! \internal \file
  37  * \brief
  38  * Implements classes in analysismodule.h.
  39  *
  40  * \author Teemu Murtola <teemu.murtola@gmail.com>
  41  * \ingroup module_trajectoryanalysis
  42  */
  43 #include "gmxpre.h"
  44
  45 #include "analysismodule.h"
  46
  47 #include <map>
  48 #include <utility>
  49
  50 #include "gromacs/analysisdata/analysisdata.h"
  51 #include "gromacs/selection/selection.h"
  52 #include "gromacs/selection/selectioncollection.h"
  53 #include "gromacs/utility/exceptions.h"
  54 #include "gromacs/utility/gmxassert.h"
  55 #include "gromacs/utility/stringutil.h"
  56
  57 namespace gmx
  58 {
  59
  60 /********************************************************************
  61  * TrajectoryAnalysisModule::Impl
  62  */
  63
  64 /*! \internal \brief
  65  * Private implementation class for TrajectoryAnalysisModule.
  66  *
  67  * \ingroup module_trajectoryanalysis
  68  */
  69 class TrajectoryAnalysisModule::Impl
  70 {
  71 public:
  72     //! Container that associates a data set with its name.
  73     typedef std::map<std::string, AbstractAnalysisData*> DatasetContainer;
  74     //! Container that associates a AnalysisData object with its name.
  75     typedef std::map<std::string, AnalysisData*> AnalysisDatasetContainer;
  76
  77     //! List of registered data set names.
  78     std::vector<std::string> datasetNames_;
  79     /*! \brief
  80      * Keeps all registered data sets.
  81      *
  82      * This container also includes datasets from \a analysisDatasets_.
  83      */
  84     DatasetContainer datasets_;
  85     //! Keeps registered AnalysisData objects.
  86     AnalysisDatasetContainer analysisDatasets_;
  87 };
  88
  89 /********************************************************************
  90  * TrajectoryAnalysisModuleData::Impl
  91  */
  92
  93 /*! \internal \brief
  94  * Private implementation class for TrajectoryAnalysisModuleData.
  95  *
  96  * \ingroup module_trajectoryanalysis
  97  */
  98 class TrajectoryAnalysisModuleData::Impl
  99 {
 100 public:
 101     //! Container that associates a data handle to its AnalysisData object.
 102     typedef std::map<const AnalysisData*, AnalysisDataHandle> HandleContainer;
 103
 104     //! \copydoc TrajectoryAnalysisModuleData::TrajectoryAnalysisModuleData()
 105     Impl(TrajectoryAnalysisModule*          module,
 106          const AnalysisDataParallelOptions& opt,
 107          const SelectionCollection&         selections);
 108
 109     //! Checks whether the given AnalysisData has been initialized.
 110     static bool isInitialized(const AnalysisData& data);
 111
 112     //! Keeps a data handle for each AnalysisData object.
 113     HandleContainer handles_;
 114     //! Stores thread-local selections.
 115     const SelectionCollection& selections_;
 116 };
 117
 118 TrajectoryAnalysisModuleData::Impl::Impl(TrajectoryAnalysisModule*          module,
 119                                          const AnalysisDataParallelOptions& opt,
 120                                          const SelectionCollection&         selections) :
 121     selections_(selections)
 122 {
 123     TrajectoryAnalysisModule::Impl::AnalysisDatasetContainer::const_iterator i;
 124     for (i = module->impl_->analysisDatasets_.begin(); i != module->impl_->analysisDatasets_.end(); ++i)
 125     {
 126         AnalysisDataHandle handle;
 127         if (isInitialized(*i->second))
 128         {
 129             handle = i->second->startData(opt);
 130         }
 131         handles_.insert(std::make_pair(i->second, handle));
 132     }
 133 }
 134
 135 bool TrajectoryAnalysisModuleData::Impl::isInitialized(const AnalysisData& data)
 136 {
 137     for (int i = 0; i < data.dataSetCount(); ++i)
 138     {
 139         if (data.columnCount(i) > 0)
 140         {
 141             // If not all of the column counts are set, startData() in the
 142             // constructor asserts, so that does not need to be checked here.
 143             return true;
 144         }
 145     }
 146     return false;
 147 }
 148
 149
 150 /********************************************************************
 151  * TrajectoryAnalysisModuleData
 152  */
 153
 154 TrajectoryAnalysisModuleData::TrajectoryAnalysisModuleData(TrajectoryAnalysisModule* module,
 155                                                            const AnalysisDataParallelOptions& opt,
 156                                                            const SelectionCollection& selections) :
 157     impl_(new Impl(module, opt, selections))
 158 {
 159 }
 160
 161
 162 TrajectoryAnalysisModuleData::~TrajectoryAnalysisModuleData() {}
 163
 164
 165 void TrajectoryAnalysisModuleData::finishDataHandles()
 166 {
 167     // FIXME: Call finishData() for all handles even if one throws
 168     Impl::HandleContainer::iterator i;
 169     for (i = impl_->handles_.begin(); i != impl_->handles_.end(); ++i)
 170     {
 171         if (i->second.isValid())
 172         {
 173             i->second.finishData();
 174         }
 175     }
 176     impl_->handles_.clear();
 177 }
 178
 179
 180 AnalysisDataHandle TrajectoryAnalysisModuleData::dataHandle(const AnalysisData& data)
 181 {
 182     Impl::HandleContainer::const_iterator i = impl_->handles_.find(&data);
 183     GMX_RELEASE_ASSERT(i != impl_->handles_.end(), "Data handle requested on unknown dataset");
 184     return i->second;
 185 }
 186
 187
 188 Selection TrajectoryAnalysisModuleData::parallelSelection(const Selection& selection)
 189 {
 190     std::optional<Selection> sel = impl_->selections_.selection(selection.name());
 191     // Selections should never be missing in an analysis module, so this is an internal consistency check.
 192     GMX_RELEASE_ASSERT(sel.has_value(),
 193                        gmx::formatString("invalid selection %s", selection.name()).c_str());
 194     return sel.value();
 195 }
 196
 197
 198 SelectionList TrajectoryAnalysisModuleData::parallelSelections(const SelectionList& selections)
 199 {
 200     // TODO: Consider an implementation that does not allocate memory every time.
 201     SelectionList newSelections;
 202     newSelections.reserve(selections.size());
 203     SelectionList::const_iterator i = selections.begin();
 204     for (; i != selections.end(); ++i)
 205     {
 206         newSelections.push_back(parallelSelection(*i));
 207     }
 208     return newSelections;
 209 }
 210
 211
 212 /********************************************************************
 213  * TrajectoryAnalysisModuleDataBasic
 214  */
 215
 216 namespace
 217 {
 218
 219 /*! \brief
 220  * Basic thread-local trajectory analysis data storage class.
 221  *
 222  * Most simple tools should only require data handles and selections to be
 223  * thread-local, so this class implements just that.
 224  *
 225  * \ingroup module_trajectoryanalysis
 226  */
 227 class TrajectoryAnalysisModuleDataBasic : public TrajectoryAnalysisModuleData
 228 {
 229 public:
 230     /*! \brief
 231      * Initializes thread-local storage for data handles and selections.
 232      *
 233      * \param[in] module     Analysis module to use for data objects.
 234      * \param[in] opt        Data parallelization options.
 235      * \param[in] selections Thread-local selection collection.
 236      */
 237     TrajectoryAnalysisModuleDataBasic(TrajectoryAnalysisModule*          module,
 238                                       const AnalysisDataParallelOptions& opt,
 239                                       const SelectionCollection&         selections);
 240
 241     void finish() override;
 242 };
 243
 244 TrajectoryAnalysisModuleDataBasic::TrajectoryAnalysisModuleDataBasic(TrajectoryAnalysisModule* module,
 245                                                                      const AnalysisDataParallelOptions& opt,
 246                                                                      const SelectionCollection& selections) :
 247     TrajectoryAnalysisModuleData(module, opt, selections)
 248 {
 249 }
 250
 251
 252 void TrajectoryAnalysisModuleDataBasic::finish()
 253 {
 254     finishDataHandles();
 255 }
 256
 257 } // namespace
 258
 259
 260 /********************************************************************
 261  * TrajectoryAnalysisModule
 262  */
 263
 264 TrajectoryAnalysisModule::TrajectoryAnalysisModule() : impl_(new Impl()) {}
 265
 266
 267 TrajectoryAnalysisModule::~TrajectoryAnalysisModule() {}
 268
 269
 270 void TrajectoryAnalysisModule::optionsFinished(TrajectoryAnalysisSettings* /*settings*/) {}
 271
 272
 273 void TrajectoryAnalysisModule::initAfterFirstFrame(const TrajectoryAnalysisSettings& /*settings*/,
 274                                                    const t_trxframe& /*fr*/)
 275 {
 276 }
 277
 278
 279 TrajectoryAnalysisModuleDataPointer TrajectoryAnalysisModule::startFrames(const AnalysisDataParallelOptions& opt,
 280                                                                           const SelectionCollection& selections)
 281 {
 282     return TrajectoryAnalysisModuleDataPointer(new TrajectoryAnalysisModuleDataBasic(this, opt, selections));
 283 }
 284
 285
 286 void TrajectoryAnalysisModule::finishFrames(TrajectoryAnalysisModuleData* /*pdata*/) {}
 287
 288
 289 int TrajectoryAnalysisModule::datasetCount() const
 290 {
 291     return impl_->datasetNames_.size();
 292 }
 293
 294
 295 const std::vector<std::string>& TrajectoryAnalysisModule::datasetNames() const
 296 {
 297     return impl_->datasetNames_;
 298 }
 299
 300
 301 AbstractAnalysisData& TrajectoryAnalysisModule::datasetFromIndex(int index) const
 302 {
 303     if (index < 0 || index >= datasetCount())
 304     {
 305         GMX_THROW(APIError("Out of range data set index"));
 306     }
 307     Impl::DatasetContainer::const_iterator item = impl_->datasets_.find(impl_->datasetNames_[index]);
 308     GMX_RELEASE_ASSERT(item != impl_->datasets_.end(), "Inconsistent data set names");
 309     return *item->second;
 310 }
 311
 312
 313 AbstractAnalysisData& TrajectoryAnalysisModule::datasetFromName(const char* name) const
 314 {
 315     Impl::DatasetContainer::const_iterator item = impl_->datasets_.find(name);
 316     if (item == impl_->datasets_.end())
 317     {
 318         GMX_THROW(APIError("Unknown data set name"));
 319     }
 320     return *item->second;
 321 }
 322
 323
 324 void TrajectoryAnalysisModule::registerBasicDataset(AbstractAnalysisData* data, const char* name)
 325 {
 326     GMX_RELEASE_ASSERT(data != nullptr, "Attempting to register NULL data");
 327     // TODO: Strong exception safety should be possible to implement.
 328     GMX_RELEASE_ASSERT(impl_->datasets_.find(name) == impl_->datasets_.end(),
 329                        "Duplicate data set name registered");
 330     impl_->datasets_[name] = data;
 331     impl_->datasetNames_.emplace_back(name);
 332 }
 333
 334
 335 void TrajectoryAnalysisModule::registerAnalysisDataset(AnalysisData* data, const char* name)
 336 {
 337     // TODO: Strong exception safety should be possible to implement.
 338     registerBasicDataset(data, name);
 339     impl_->analysisDatasets_[name] = data;
 340 }
 341
 342
 343 void TrajectoryAnalysisModule::finishFrameSerial(int frameIndex)
 344 {
 345     Impl::AnalysisDatasetContainer::const_iterator data;
 346     for (data = impl_->analysisDatasets_.begin(); data != impl_->analysisDatasets_.end(); ++data)
 347     {
 348         data->second->finishFrameSerial(frameIndex);
 349     }
 350 }
 351
 352 } // namespace gmx