2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2012, by the GROMACS development team, led by
5 * David van der Spoel, Berk Hess, Erik Lindahl, and including many
6 * others, as listed in the AUTHORS file in the top-level source
7 * directory and at http://www.gromacs.org.
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.
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.
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.
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.
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.
37 * Implements classes in datastorage.h and paralleloptions.h.
39 * \author Teemu Murtola <teemu.murtola@gmail.com>
40 * \ingroup module_analysisdata
42 #include "datastorage.h"
47 #include "gromacs/analysisdata/abstractdata.h"
48 #include "gromacs/analysisdata/dataframe.h"
49 #include "gromacs/analysisdata/paralleloptions.h"
50 #include "gromacs/utility/exceptions.h"
51 #include "gromacs/utility/gmxassert.h"
52 #include "gromacs/utility/uniqueptr.h"
57 /********************************************************************
58 * AnalysisDataParallelOptions
61 AnalysisDataParallelOptions::AnalysisDataParallelOptions()
62 : parallelizationFactor_(1)
67 AnalysisDataParallelOptions::AnalysisDataParallelOptions(int parallelizationFactor)
68 : parallelizationFactor_(parallelizationFactor)
70 GMX_RELEASE_ASSERT(parallelizationFactor >= 1,
71 "Invalid parallelization factor");
75 /********************************************************************
76 * AnalysisDataStorage::Impl
80 * Private implementation class for AnalysisDataStorage.
82 * \ingroup module_analysisdata
84 class AnalysisDataStorage::Impl
87 //! Smart pointer type for managing a stored frame.
88 typedef gmx_unique_ptr<AnalysisDataStorageFrame>::type FramePointer;
91 * Stored information about a single stored frame.
93 * Methods in this class do not throw.
97 //! Indicates what operations have been performed on a frame.
100 eMissing, //!< Frame has not yet been started.
101 eStarted, //!< startFrame() has been called.
102 eFinished, //!< finishFrame() has been called.
103 eNotified //!< Appropriate notifications have been sent.
106 //! Constructs an object that manages a given frame object.
107 explicit StoredFrame(AnalysisDataStorageFrame *frame)
108 : frame(frame), status(eMissing)
111 //! Whether the frame has been started with startFrame().
112 bool isStarted() const { return status >= eStarted; }
113 //! Whether the frame has been finished with finishFrame().
114 bool isFinished() const { return status >= eFinished; }
115 //! Whether all notifications have been sent.
116 bool isNotified() const { return status >= eNotified; }
117 //! Whether the frame is ready to be available outside the storage.
118 bool isAvailable() const { return status >= eFinished; }
126 //! In what state the frame currently is.
130 //! Shorthand for a list of data frames that are currently stored.
131 typedef std::vector<StoredFrame> FrameList;
135 //! Returns the number of columns in the attached data.
136 int columnCount() const;
137 //! Returns whether the storage is set to use multipoint data.
138 bool isMultipoint() const;
140 * Whether storage of all frames has been requested.
142 * Storage of all frames also works as expected if \a storageLimit_ is
143 * used in comparisons directly, but this method should be used to
144 * check how to manage \a frames_.
146 bool storeAll() const
148 return storageLimit_ == std::numeric_limits<int>::max();
150 //! Returns the index of the oldest frame that may be currently stored.
151 int firstStoredIndex() const;
153 * Computes index into \a frames_ for accessing frame \p index.
155 * \param[in] index Zero-based frame index.
156 * \retval -1 if \p index is not available in \a frames_.
160 int computeStorageLocation(int index) const;
163 * Computes an index into \a frames_ that is one past the last frame
168 size_t endStorageLocation() const;
171 * Extends \a frames_ to a new size.
173 * \throws std::bad_alloc if out of memory.
175 void extendBuffer(AnalysisDataStorage *storage, size_t newSize);
177 * Remove oldest frame from the storage to make space for a new one.
179 * Increments \a firstFrameLocation_ and reinitializes the frame that
180 * was made unavailable by this operation.
189 * Calls notification method in \a data_.
191 * \throws unspecified Any exception thrown by
192 * AbstractAnalysisData::notifyPointsAdd().
194 void notifyPointSet(const AnalysisDataPointSetRef &points);
196 * Calls notification methods for new frames.
198 * \param[in] firstLocation First frame to consider.
199 * \throws unspecified Any exception thrown by frame notification
200 * methods in AbstractAnalysisData.
202 * Notifies \a data_ of new frames (from \p firstLocation and after
203 * that) if all previous frames have already been notified.
204 * Also rotates the \a frames_ buffer as necessary.
206 void notifyNextFrames(size_t firstLocation);
208 //! Data object to use for notification calls.
209 AbstractAnalysisData *data_;
211 * Whether the storage has been set to allow multipoint.
213 * Should be possible to remove once full support for multipoint data
214 * has been implemented; isMultipoint() can simply return
215 * \c data_->isMultipoint() in that case.
219 * Number of past frames that need to be stored.
221 * Always non-negative. If storage of all frames has been requested,
222 * this is set to a large number.
226 * Number of future frames that may need to be started.
228 * Should always be at least one.
230 * \see AnalysisDataStorage::startFrame()
234 * Data frames that are currently stored.
236 * If storage of all frames has been requested, this is simply a vector
237 * of frames up to the latest frame that has been started.
238 * In this case, \a firstFrameLocation_ is always zero.
240 * If storage of all frames is not requested, this is a ring buffer of
241 * frames of size \c n=storageLimit_+pendingLimit_+1. If a frame with
242 * index \c index is currently stored, its location is
243 * \c index%frames_.size().
244 * When at most \a storageLimit_ first frames have been finished,
245 * this contains storage for the first \c n-1 frames.
246 * When more than \a storageLimit_ first frames have been finished,
247 * the oldest stored frame is stored in the location
248 * \a firstFrameLocation_, and \a storageLimit_ frames starting from
249 * this location are the last finished frames. \a pendingLimit_ frames
250 * follow, and some of these may be in progress or finished.
251 * There is always one unused frame in the buffer, which is initialized
252 * such that when \a firstFrameLocation_ is incremented, it becomes
253 * valid. This makes it easier to rotate the buffer in concurrent
254 * access scenarions (which are not yet otherwise implemented).
257 //! Location of oldest frame in \a frames_.
258 size_t firstFrameLocation_;
260 * Index of next frame that will be added to \a frames_.
262 * If all frames are not stored, this will be the index of the unused
263 * frame (see \a frames_).
268 AnalysisDataStorage::Impl::Impl()
269 : data_(NULL), bMultipoint_(false),
270 storageLimit_(0), pendingLimit_(1), firstFrameLocation_(0), nextIndex_(0)
276 AnalysisDataStorage::Impl::columnCount() const
278 GMX_ASSERT(data_ != NULL, "columnCount() called too early");
279 return data_->columnCount();
284 AnalysisDataStorage::Impl::isMultipoint() const
291 AnalysisDataStorage::Impl::firstStoredIndex() const
293 return frames_[firstFrameLocation_].frame->frameIndex();
298 AnalysisDataStorage::Impl::computeStorageLocation(int index) const
300 if (index < firstStoredIndex() || index >= nextIndex_)
304 return index % frames_.size();
309 AnalysisDataStorage::Impl::endStorageLocation() const
313 return frames_.size();
315 if (frames_[0].frame->frameIndex() == 0 || firstFrameLocation_ == 0)
317 return frames_.size() - 1;
319 return firstFrameLocation_ - 1;
324 AnalysisDataStorage::Impl::extendBuffer(AnalysisDataStorage *storage,
327 frames_.reserve(newSize);
328 while (frames_.size() < newSize)
330 frames_.push_back(StoredFrame(
331 new AnalysisDataStorageFrame(storage, columnCount(), nextIndex_)));
334 // The unused frame should not be included in the count.
343 AnalysisDataStorage::Impl::rotateBuffer()
345 GMX_ASSERT(!storeAll(),
346 "No need to rotate internal buffer if everything is stored");
347 size_t prevFirst = firstFrameLocation_;
348 size_t nextFirst = prevFirst + 1;
349 if (nextFirst == frames_.size())
353 firstFrameLocation_ = nextFirst;
354 StoredFrame &prevFrame = frames_[prevFirst];
355 prevFrame.status = StoredFrame::eMissing;
356 prevFrame.frame->header_ = AnalysisDataFrameHeader(nextIndex_ + 1, 0.0, 0.0);
357 prevFrame.frame->clearValues();
363 AnalysisDataStorage::Impl::notifyPointSet(const AnalysisDataPointSetRef &points)
365 data_->notifyPointsAdd(points);
370 AnalysisDataStorage::Impl::notifyNextFrames(size_t firstLocation)
372 if (firstLocation != firstFrameLocation_)
374 // firstLocation can only be zero here if !storeAll() because
375 // firstFrameLocation_ is always zero for storeAll()
377 (firstLocation == 0 ? frames_.size() - 1 : firstLocation - 1);
378 if (!frames_[prevIndex].isNotified())
383 size_t i = firstLocation;
384 size_t end = endStorageLocation();
387 Impl::StoredFrame &storedFrame = frames_[i];
388 if (!storedFrame.isFinished())
392 if (storedFrame.status == StoredFrame::eFinished)
394 data_->notifyFrameStart(storedFrame.frame->header());
395 data_->notifyPointsAdd(storedFrame.frame->currentPoints());
396 data_->notifyFrameFinish(storedFrame.frame->header());
397 storedFrame.status = StoredFrame::eNotified;
398 if (storedFrame.frame->frameIndex() >= storageLimit_)
404 if (!storeAll() && i >= frames_.size())
412 /********************************************************************
413 * AnalysisDataStorageFrame
416 AnalysisDataStorageFrame::AnalysisDataStorageFrame(AnalysisDataStorage *storage,
417 int columnCount, int index)
418 : storage_(*storage), header_(index, 0.0, 0.0), values_(columnCount)
423 AnalysisDataStorageFrame::~AnalysisDataStorageFrame()
428 AnalysisDataPointSetRef
429 AnalysisDataStorageFrame::currentPoints() const
431 std::vector<AnalysisDataValue>::const_iterator begin = values_.begin();
432 std::vector<AnalysisDataValue>::const_iterator end = values_.end();
433 while (begin != end && !begin->isSet())
437 while (end != begin && !(end-1)->isSet())
441 int firstColumn = (begin != end) ? begin - values_.begin() : 0;
442 return AnalysisDataPointSetRef(header_, firstColumn,
443 AnalysisDataValuesRef(begin, end));
448 AnalysisDataStorageFrame::clearValues()
450 std::vector<AnalysisDataValue>::iterator i;
451 for (i = values_.begin(); i != values_.end(); ++i)
459 AnalysisDataStorageFrame::finishPointSet()
461 storage_.impl_->notifyPointSet(currentPoints());
466 /********************************************************************
467 * AnalysisDataStorage
470 AnalysisDataStorage::AnalysisDataStorage()
476 AnalysisDataStorage::~AnalysisDataStorage()
482 AnalysisDataStorage::setMultipoint(bool bMultipoint)
484 if (bMultipoint && impl_->storageLimit_ > 0)
486 GMX_THROW(APIError("Storage of multipoint data not supported"));
488 impl_->bMultipoint_ = bMultipoint;
493 AnalysisDataStorage::setParallelOptions(const AnalysisDataParallelOptions &opt)
495 impl_->pendingLimit_ = 2 * opt.parallelizationFactor() - 1;
500 AnalysisDataStorage::tryGetDataFrame(int index) const
502 if (impl_->isMultipoint())
504 return AnalysisDataFrameRef();
506 int storageIndex = impl_->computeStorageLocation(index);
507 if (storageIndex == -1)
509 return AnalysisDataFrameRef();
511 const Impl::StoredFrame &storedFrame = impl_->frames_[storageIndex];
512 if (!storedFrame.isAvailable())
514 return AnalysisDataFrameRef();
516 const Impl::FramePointer &frame = storedFrame.frame;
517 return AnalysisDataFrameRef(frame->header(), frame->values_);
522 AnalysisDataStorage::requestStorage(int nframes)
524 if (impl_->isMultipoint())
529 // Handle the case when everything needs to be stored.
532 impl_->storageLimit_ = std::numeric_limits<int>::max();
535 // Check whether an earlier call has requested more storage.
536 if (nframes < impl_->storageLimit_)
540 impl_->storageLimit_ = nframes;
546 AnalysisDataStorage::startDataStorage(AbstractAnalysisData *data)
548 // Data needs to be set before calling extendBuffer()
550 setMultipoint(data->isMultipoint());
551 if (!impl_->storeAll())
553 impl_->extendBuffer(this, impl_->storageLimit_ + impl_->pendingLimit_ + 1);
558 AnalysisDataStorageFrame &
559 AnalysisDataStorage::startFrame(const AnalysisDataFrameHeader &header)
561 GMX_ASSERT(header.isValid(), "Invalid header");
562 Impl::StoredFrame *storedFrame;
563 if (impl_->storeAll())
565 size_t size = header.index() + 1;
566 if (impl_->frames_.size() < size)
568 impl_->extendBuffer(this, size);
570 storedFrame = &impl_->frames_[header.index()];
574 int storageIndex = impl_->computeStorageLocation(header.index());
575 if (storageIndex == -1)
577 GMX_THROW(APIError("Out of bounds frame index"));
579 storedFrame = &impl_->frames_[storageIndex];
581 GMX_RELEASE_ASSERT(!storedFrame->isStarted(),
582 "startFrame() called twice for the same frame");
583 GMX_RELEASE_ASSERT(storedFrame->frame->frameIndex() == header.index(),
584 "Inconsistent internal frame indexing");
585 storedFrame->status = Impl::StoredFrame::eStarted;
586 storedFrame->frame->header_ = header;
587 if (impl_->isMultipoint())
589 impl_->data_->notifyFrameStart(header);
591 return *storedFrame->frame;
595 AnalysisDataStorageFrame &
596 AnalysisDataStorage::startFrame(int index, real x, real dx)
598 return startFrame(AnalysisDataFrameHeader(index, x, dx));
602 AnalysisDataStorageFrame &
603 AnalysisDataStorage::currentFrame(int index)
605 int storageIndex = impl_->computeStorageLocation(index);
606 GMX_RELEASE_ASSERT(storageIndex >= 0, "Out of bounds frame index");
607 Impl::StoredFrame &storedFrame = impl_->frames_[storageIndex];
608 GMX_RELEASE_ASSERT(storedFrame.isStarted(),
609 "currentFrame() called for frame before startFrame()");
610 GMX_RELEASE_ASSERT(!storedFrame.isFinished(),
611 "currentFrame() called for frame after finishFrame()");
612 GMX_RELEASE_ASSERT(storedFrame.frame->frameIndex() == index,
613 "Inconsistent internal frame indexing");
614 return *storedFrame.frame;
619 AnalysisDataStorage::finishFrame(int index)
621 int storageIndex = impl_->computeStorageLocation(index);
622 GMX_RELEASE_ASSERT(storageIndex >= 0, "Out of bounds frame index");
623 Impl::StoredFrame &storedFrame = impl_->frames_[storageIndex];
624 GMX_RELEASE_ASSERT(storedFrame.isStarted(),
625 "finishFrame() called for frame before startFrame()");
626 GMX_RELEASE_ASSERT(!storedFrame.isFinished(),
627 "finishFrame() called twice for the same frame");
628 GMX_RELEASE_ASSERT(storedFrame.frame->frameIndex() == index,
629 "Inconsistent internal frame indexing");
630 storedFrame.status = Impl::StoredFrame::eFinished;
631 if (impl_->isMultipoint())
633 // TODO: Check that the last point set has been finished
634 impl_->data_->notifyFrameFinish(storedFrame.frame->header());
635 if (storedFrame.frame->frameIndex() >= impl_->storageLimit_)
637 impl_->rotateBuffer();
642 impl_->notifyNextFrames(storageIndex);
648 AnalysisDataStorage::finishFrame(const AnalysisDataStorageFrame &frame)
650 finishFrame(frame.frameIndex());