Merge branch release-2021 into master
[alexxy/gromacs.git] / src / testutils / refdata.cpp
1 /*
2  * This file is part of the GROMACS molecular simulation package.
3  *
4  * Copyright (c) 2011-2018, The GROMACS development team.
5  * Copyright (c) 2019,2020, 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 and functions from refdata.h.
39  *
40  * \author Teemu Murtola <teemu.murtola@gmail.com>
41  * \ingroup module_testutils
42  */
43 #include "gmxpre.h"
44
45 #include "testutils/refdata.h"
46
47 #include <cctype>
48 #include <cstdlib>
49
50 #include <algorithm>
51 #include <limits>
52 #include <string>
53
54 #include <gtest/gtest.h>
55
56 #include "gromacs/math/vectypes.h"
57 #include "gromacs/options/basicoptions.h"
58 #include "gromacs/options/ioptionscontainer.h"
59 #include "gromacs/utility/any.h"
60 #include "gromacs/utility/exceptions.h"
61 #include "gromacs/utility/gmxassert.h"
62 #include "gromacs/utility/keyvaluetree.h"
63 #include "gromacs/utility/path.h"
64 #include "gromacs/utility/real.h"
65 #include "gromacs/utility/stringutil.h"
66
67 #include "testutils/testasserts.h"
68 #include "testutils/testexceptions.h"
69 #include "testutils/testfilemanager.h"
70
71 #include "refdata_checkers.h"
72 #include "refdata_impl.h"
73 #include "refdata_xml.h"
74
75 namespace gmx
76 {
77 namespace test
78 {
79
80 /********************************************************************
81  * TestReferenceData::Impl declaration
82  */
83
84 namespace internal
85 {
86
87 /*! \internal \brief
88  * Private implementation class for TestReferenceData.
89  *
90  * \ingroup module_testutils
91  */
92 class TestReferenceDataImpl
93 {
94 public:
95     //! Initializes a checker in the given mode.
96     TestReferenceDataImpl(ReferenceDataMode mode, bool bSelfTestMode);
97
98     //! Performs final reference data processing when test ends.
99     void onTestEnd(bool testPassed);
100
101     //! Full path of the reference data file.
102     std::string fullFilename_;
103     /*! \brief
104      * Root entry for comparing the reference data.
105      *
106      * Null after construction iff in compare mode and reference data was
107      * not loaded successfully.
108      * In all write modes, copies are present for nodes added to
109      * \a outputRootEntry_, and ReferenceDataEntry::correspondingOutputEntry()
110      * points to the copy in the output tree.
111      */
112     ReferenceDataEntry::EntryPointer compareRootEntry_;
113     /*! \brief
114      * Root entry for writing new reference data.
115      *
116      * Null if only comparing against existing data.  Otherwise, starts
117      * always as empty.
118      * When creating new reference data, this is maintained as a copy of
119      * \a compareRootEntry_.
120      * When updating existing data, entries are added either by copying
121      * from \a compareRootEntry_ (if they exist and comparison passes), or
122      * by creating new ones.
123      */
124     ReferenceDataEntry::EntryPointer outputRootEntry_;
125     /*! \brief
126      * Whether updating existing reference data.
127      */
128     bool updateMismatchingEntries_;
129     //! `true` if self-testing (enables extra failure messages).
130     bool bSelfTestMode_;
131     /*! \brief
132      * Whether any reference checkers have been created for this data.
133      */
134     bool bInUse_;
135 };
136
137 } // namespace internal
138
139 /********************************************************************
140  * Internal helpers
141  */
142
143 namespace
144 {
145
146 //! Convenience typedef for a smart pointer to TestReferenceDataImpl.
147 typedef std::shared_ptr<internal::TestReferenceDataImpl> TestReferenceDataImplPointer;
148
149 /*! \brief
150  * Global reference data instance.
151  *
152  * The object is created when the test creates a TestReferenceData, and the
153  * object is destructed (and other post-processing is done) at the end of each
154  * test by ReferenceDataTestEventListener (which is installed as a Google Test
155  * test listener).
156  */
157 TestReferenceDataImplPointer g_referenceData;
158 //! Global reference data mode set with setReferenceDataMode().
159 ReferenceDataMode g_referenceDataMode = ReferenceDataMode::Compare;
160
161 //! Returns the global reference data mode.
162 ReferenceDataMode getReferenceDataMode()
163 {
164     return g_referenceDataMode;
165 }
166
167 //! Returns a reference to the global reference data object.
168 TestReferenceDataImplPointer initReferenceDataInstance()
169 {
170     GMX_RELEASE_ASSERT(!g_referenceData, "Test cannot create multiple TestReferenceData instances");
171     g_referenceData.reset(new internal::TestReferenceDataImpl(getReferenceDataMode(), false));
172     return g_referenceData;
173 }
174
175 //! Handles reference data creation for self-tests.
176 TestReferenceDataImplPointer initReferenceDataInstanceForSelfTest(ReferenceDataMode mode)
177 {
178     if (g_referenceData)
179     {
180         GMX_RELEASE_ASSERT(g_referenceData.unique(),
181                            "Test cannot create multiple TestReferenceData instances");
182         g_referenceData->onTestEnd(true);
183         g_referenceData.reset();
184     }
185     g_referenceData.reset(new internal::TestReferenceDataImpl(mode, true));
186     return g_referenceData;
187 }
188
189 class ReferenceDataTestEventListener : public ::testing::EmptyTestEventListener
190 {
191 public:
192     void OnTestEnd(const ::testing::TestInfo& test_info) override
193     {
194         if (g_referenceData)
195         {
196             GMX_RELEASE_ASSERT(g_referenceData.unique(), "Test leaked TestRefeferenceData objects");
197             g_referenceData->onTestEnd(test_info.result()->Passed());
198             g_referenceData.reset();
199         }
200     }
201
202     void OnTestProgramEnd(const ::testing::UnitTest& /*unused*/) override
203     {
204         // Could be used e.g. to free internal buffers allocated by an XML parsing library
205     }
206 };
207
208 //! Formats a path to a reference data entry with a non-null id.
209 std::string formatEntryPath(const std::string& prefix, const std::string& id)
210 {
211     return prefix + "/" + id;
212 }
213
214 //! Formats a path to a reference data entry with a null id.
215 std::string formatSequenceEntryPath(const std::string& prefix, int seqIndex)
216 {
217     return formatString("%s/[%d]", prefix.c_str(), seqIndex + 1);
218 }
219
220 //! Finds all entries that have not been checked under a given root.
221 void gatherUnusedEntries(const ReferenceDataEntry& root,
222                          const std::string&        rootPath,
223                          std::vector<std::string>* unusedPaths)
224 {
225     if (!root.hasBeenChecked())
226     {
227         unusedPaths->push_back(rootPath);
228         return;
229     }
230     int seqIndex = 0;
231     for (const auto& child : root.children())
232     {
233         std::string path;
234         if (child->id().empty())
235         {
236             path = formatSequenceEntryPath(rootPath, seqIndex);
237             ++seqIndex;
238         }
239         else
240         {
241             path = formatEntryPath(rootPath, child->id());
242         }
243         gatherUnusedEntries(*child, path, unusedPaths);
244     }
245 }
246
247 //! Produces a GTest assertion of any entries under given root have not been checked.
248 void checkUnusedEntries(const ReferenceDataEntry& root, const std::string& rootPath)
249 {
250     std::vector<std::string> unusedPaths;
251     gatherUnusedEntries(root, rootPath, &unusedPaths);
252     if (!unusedPaths.empty())
253     {
254         std::string paths;
255         if (unusedPaths.size() > 5)
256         {
257             paths = joinStrings(unusedPaths.begin(), unusedPaths.begin() + 5, "\n  ");
258             paths = "  " + paths + "\n  ...";
259         }
260         else
261         {
262             paths = joinStrings(unusedPaths.begin(), unusedPaths.end(), "\n  ");
263             paths = "  " + paths;
264         }
265         ADD_FAILURE() << "Reference data items not used in test:" << std::endl << paths;
266     }
267 }
268
269 } // namespace
270
271 void initReferenceData(IOptionsContainer* options)
272 {
273     static const gmx::EnumerationArray<ReferenceDataMode, const char*> s_refDataNames = {
274         { "check", "create", "update-changed", "update-all" }
275     };
276     options->addOption(EnumOption<ReferenceDataMode>("ref-data")
277                                .enumValue(s_refDataNames)
278                                .store(&g_referenceDataMode)
279                                .description("Operation mode for tests that use reference data"));
280     ::testing::UnitTest::GetInstance()->listeners().Append(new ReferenceDataTestEventListener);
281 }
282
283 /********************************************************************
284  * TestReferenceDataImpl definition
285  */
286
287 namespace internal
288 {
289
290 TestReferenceDataImpl::TestReferenceDataImpl(ReferenceDataMode mode, bool bSelfTestMode) :
291     updateMismatchingEntries_(false),
292     bSelfTestMode_(bSelfTestMode),
293     bInUse_(false)
294 {
295     const std::string dirname = bSelfTestMode ? TestFileManager::getGlobalOutputTempDirectory()
296                                               : TestFileManager::getInputDataDirectory();
297     const std::string filename = TestFileManager::getTestSpecificFileName(".xml");
298     fullFilename_              = Path::join(dirname, "refdata", filename);
299
300     switch (mode)
301     {
302         case ReferenceDataMode::Compare:
303             if (File::exists(fullFilename_, File::throwOnError))
304             {
305                 compareRootEntry_ = readReferenceDataFile(fullFilename_);
306             }
307             break;
308         case ReferenceDataMode::CreateMissing:
309             if (File::exists(fullFilename_, File::throwOnError))
310             {
311                 compareRootEntry_ = readReferenceDataFile(fullFilename_);
312             }
313             else
314             {
315                 compareRootEntry_ = ReferenceDataEntry::createRoot();
316                 outputRootEntry_  = ReferenceDataEntry::createRoot();
317             }
318             break;
319         case ReferenceDataMode::UpdateChanged:
320             if (File::exists(fullFilename_, File::throwOnError))
321             {
322                 compareRootEntry_ = readReferenceDataFile(fullFilename_);
323             }
324             else
325             {
326                 compareRootEntry_ = ReferenceDataEntry::createRoot();
327             }
328             outputRootEntry_          = ReferenceDataEntry::createRoot();
329             updateMismatchingEntries_ = true;
330             break;
331         case ReferenceDataMode::UpdateAll:
332             compareRootEntry_ = ReferenceDataEntry::createRoot();
333             outputRootEntry_  = ReferenceDataEntry::createRoot();
334             break;
335         case ReferenceDataMode::Count: GMX_THROW(InternalError("Invalid reference data mode"));
336     }
337 }
338
339 void TestReferenceDataImpl::onTestEnd(bool testPassed)
340 {
341     if (!bInUse_)
342     {
343         return;
344     }
345     // TODO: Only write the file with update-changed if there were actual changes.
346     if (outputRootEntry_)
347     {
348         if (testPassed)
349         {
350             std::string dirname = Path::getParentPath(fullFilename_);
351             if (!Directory::exists(dirname))
352             {
353                 if (Directory::create(dirname) != 0)
354                 {
355                     GMX_THROW(TestException("Creation of reference data directory failed: " + dirname));
356                 }
357             }
358             writeReferenceDataFile(fullFilename_, *outputRootEntry_);
359         }
360     }
361     else if (compareRootEntry_)
362     {
363         checkUnusedEntries(*compareRootEntry_, "");
364     }
365 }
366
367 } // namespace internal
368
369
370 /********************************************************************
371  * TestReferenceChecker::Impl
372  */
373
374 /*! \internal \brief
375  * Private implementation class for TestReferenceChecker.
376  *
377  * \ingroup module_testutils
378  */
379 class TestReferenceChecker::Impl
380 {
381 public:
382     //! String constant for naming XML elements for boolean values.
383     static const char* const cBooleanNodeName;
384     //! String constant for naming XML elements for string values.
385     static const char* const cStringNodeName;
386     //! String constant for naming XML elements for unsigned char values.
387     static const char* const cUCharNodeName;
388     //! String constant for naming XML elements for integer values.
389     static const char* const cIntegerNodeName;
390     //! String constant for naming XML elements for int32 values.
391     static const char* const cInt32NodeName;
392     //! String constant for naming XML elements for unsigned int32 values.
393     static const char* const cUInt32NodeName;
394     //! String constant for naming XML elements for int32 values.
395     static const char* const cInt64NodeName;
396     //! String constant for naming XML elements for unsigned int64 values.
397     static const char* const cUInt64NodeName;
398     //! String constant for naming XML elements for floating-point values.
399     static const char* const cRealNodeName;
400     //! String constant for naming XML attribute for value identifiers.
401     static const char* const cIdAttrName;
402     //! String constant for naming compounds for vectors.
403     static const char* const cVectorType;
404     //! String constant for naming compounds for key-value tree objects.
405     static const char* const cObjectType;
406     //! String constant for naming compounds for sequences.
407     static const char* const cSequenceType;
408     //! String constant for value identifier for sequence length.
409     static const char* const cSequenceLengthName;
410
411     //! Creates a checker that does nothing.
412     explicit Impl(bool initialized);
413     //! Creates a checker with a given root entry.
414     Impl(const std::string&            path,
415          ReferenceDataEntry*           compareRootEntry,
416          ReferenceDataEntry*           outputRootEntry,
417          bool                          updateMismatchingEntries,
418          bool                          bSelfTestMode,
419          const FloatingPointTolerance& defaultTolerance);
420
421     //! Returns the path of this checker with \p id appended.
422     std::string appendPath(const char* id) const;
423
424     //! Creates an entry with given parameters and fills it with \p checker.
425     static ReferenceDataEntry::EntryPointer createEntry(const char*                       type,
426                                                         const char*                       id,
427                                                         const IReferenceDataEntryChecker& checker)
428     {
429         ReferenceDataEntry::EntryPointer entry(new ReferenceDataEntry(type, id));
430         checker.fillEntry(entry.get());
431         return entry;
432     }
433     //! Checks an entry for correct type and using \p checker.
434     static ::testing::AssertionResult checkEntry(const ReferenceDataEntry&         entry,
435                                                  const std::string&                fullId,
436                                                  const char*                       type,
437                                                  const IReferenceDataEntryChecker& checker)
438     {
439         if (entry.type() != type)
440         {
441             return ::testing::AssertionFailure() << "Mismatching reference data item type" << std::endl
442                                                  << "  In item: " << fullId << std::endl
443                                                  << "   Actual: " << type << std::endl
444                                                  << "Reference: " << entry.type();
445         }
446         return checker.checkEntry(entry, fullId);
447     }
448     //! Finds an entry by id and updates the last found entry pointer.
449     ReferenceDataEntry* findEntry(const char* id);
450     /*! \brief
451      * Finds/creates a reference data entry to match against.
452      *
453      * \param[in]  type   Type of entry to create.
454      * \param[in]  id     Unique identifier of the entry (can be NULL, in
455      *      which case the next entry without an id is matched).
456      * \param[out] checker  Checker to use for filling out created entries.
457      * \returns    Matching entry, or NULL if no matching entry found
458      *      (NULL is never returned in write mode; new entries are created
459      *      instead).
460      */
461     ReferenceDataEntry* findOrCreateEntry(const char*                       type,
462                                           const char*                       id,
463                                           const IReferenceDataEntryChecker& checker);
464     /*! \brief
465      * Helper method for checking a reference data value.
466      *
467      * \param[in]  name   Type of entry to find.
468      * \param[in]  id     Unique identifier of the entry (can be NULL, in
469      *     which case the next entry without an id is matched).
470      * \param[in]  checker  Checker that provides logic specific to the
471      *     type of the entry.
472      * \returns    Whether the reference data matched, including details
473      *     of the mismatch if the comparison failed.
474      * \throws     TestException if there is a problem parsing the
475      *     reference data.
476      *
477      * Performs common tasks in checking a reference value, such as
478      * finding or creating the correct entry.
479      * Caller needs to provide a checker object that provides the string
480      * value for a newly created entry and performs the actual comparison
481      * against a found entry.
482      */
483     ::testing::AssertionResult processItem(const char*                       name,
484                                            const char*                       id,
485                                            const IReferenceDataEntryChecker& checker);
486     /*! \brief
487      * Whether the checker is initialized.
488      */
489     bool initialized() const { return initialized_; }
490     /*! \brief
491      * Whether the checker should ignore all validation calls.
492      *
493      * This is used to ignore any calls within compounds for which
494      * reference data could not be found, such that only one error is
495      * issued for the missing compound, instead of every individual value.
496      */
497     bool shouldIgnore() const
498     {
499         GMX_RELEASE_ASSERT(initialized(), "Accessing uninitialized reference data checker.");
500         return compareRootEntry_ == nullptr;
501     }
502
503     //! Whether initialized with other means than the default constructor.
504     bool initialized_;
505     //! Default floating-point comparison tolerance.
506     FloatingPointTolerance defaultTolerance_;
507     /*! \brief
508      * Human-readable path to the root node of this checker.
509      *
510      * For the root checker, this will be "/", and for each compound, the
511      * id of the compound is added.  Used for reporting comparison
512      * mismatches.
513      */
514     std::string path_;
515     /*! \brief
516      * Current entry under which reference data is searched for comparison.
517      *
518      * Points to either the TestReferenceDataImpl::compareRootEntry_, or to
519      * a compound entry in the tree rooted at that entry.
520      *
521      * Can be NULL, in which case this checker does nothing (doesn't even
522      * report errors, see shouldIgnore()).
523      */
524     ReferenceDataEntry* compareRootEntry_;
525     /*! \brief
526      * Current entry under which entries for writing are created.
527      *
528      * Points to either the TestReferenceDataImpl::outputRootEntry_, or to
529      * a compound entry in the tree rooted at that entry.  NULL if only
530      * comparing, or if shouldIgnore() returns `false`.
531      */
532     ReferenceDataEntry* outputRootEntry_;
533     /*! \brief
534      * Iterator to a child of \a compareRootEntry_ that was last found.
535      *
536      * If `compareRootEntry_->isValidChild()` returns false, no entry has
537      * been found yet.
538      * After every check, is updated to point to the entry that was used
539      * for the check.
540      * Subsequent checks start the search for the matching node on this
541      * node.
542      */
543     ReferenceDataEntry::ChildIterator lastFoundEntry_;
544     /*! \brief
545      * Whether the reference data is being written (true) or compared
546      * (false).
547      */
548     bool updateMismatchingEntries_;
549     //! `true` if self-testing (enables extra failure messages).
550     bool bSelfTestMode_;
551     /*! \brief
552      * Current number of unnamed elements in a sequence.
553      *
554      * It is the index of the current unnamed element.
555      */
556     int seqIndex_;
557 };
558
559 const char* const TestReferenceChecker::Impl::cBooleanNodeName    = "Bool";
560 const char* const TestReferenceChecker::Impl::cStringNodeName     = "String";
561 const char* const TestReferenceChecker::Impl::cUCharNodeName      = "UChar";
562 const char* const TestReferenceChecker::Impl::cIntegerNodeName    = "Int";
563 const char* const TestReferenceChecker::Impl::cInt32NodeName      = "Int32";
564 const char* const TestReferenceChecker::Impl::cUInt32NodeName     = "UInt32";
565 const char* const TestReferenceChecker::Impl::cInt64NodeName      = "Int64";
566 const char* const TestReferenceChecker::Impl::cUInt64NodeName     = "UInt64";
567 const char* const TestReferenceChecker::Impl::cRealNodeName       = "Real";
568 const char* const TestReferenceChecker::Impl::cIdAttrName         = "Name";
569 const char* const TestReferenceChecker::Impl::cVectorType         = "Vector";
570 const char* const TestReferenceChecker::Impl::cObjectType         = "Object";
571 const char* const TestReferenceChecker::Impl::cSequenceType       = "Sequence";
572 const char* const TestReferenceChecker::Impl::cSequenceLengthName = "Length";
573
574
575 TestReferenceChecker::Impl::Impl(bool initialized) :
576     initialized_(initialized),
577     defaultTolerance_(defaultRealTolerance()),
578     compareRootEntry_(nullptr),
579     outputRootEntry_(nullptr),
580     updateMismatchingEntries_(false),
581     bSelfTestMode_(false),
582     seqIndex_(-1)
583 {
584 }
585
586
587 TestReferenceChecker::Impl::Impl(const std::string&            path,
588                                  ReferenceDataEntry*           compareRootEntry,
589                                  ReferenceDataEntry*           outputRootEntry,
590                                  bool                          updateMismatchingEntries,
591                                  bool                          bSelfTestMode,
592                                  const FloatingPointTolerance& defaultTolerance) :
593     initialized_(true),
594     defaultTolerance_(defaultTolerance),
595     path_(path),
596     compareRootEntry_(compareRootEntry),
597     outputRootEntry_(outputRootEntry),
598     lastFoundEntry_(compareRootEntry->children().end()),
599     updateMismatchingEntries_(updateMismatchingEntries),
600     bSelfTestMode_(bSelfTestMode),
601     seqIndex_(-1)
602 {
603 }
604
605
606 std::string TestReferenceChecker::Impl::appendPath(const char* id) const
607 {
608     return id != nullptr ? formatEntryPath(path_, id) : formatSequenceEntryPath(path_, seqIndex_);
609 }
610
611
612 ReferenceDataEntry* TestReferenceChecker::Impl::findEntry(const char* id)
613 {
614     ReferenceDataEntry::ChildIterator entry = compareRootEntry_->findChild(id, lastFoundEntry_);
615     seqIndex_                               = (id == nullptr) ? seqIndex_ + 1 : -1;
616     if (compareRootEntry_->isValidChild(entry))
617     {
618         lastFoundEntry_ = entry;
619         return entry->get();
620     }
621     return nullptr;
622 }
623
624 ReferenceDataEntry* TestReferenceChecker::Impl::findOrCreateEntry(const char* type,
625                                                                   const char* id,
626                                                                   const IReferenceDataEntryChecker& checker)
627 {
628     ReferenceDataEntry* entry = findEntry(id);
629     if (entry == nullptr && outputRootEntry_ != nullptr)
630     {
631         lastFoundEntry_ = compareRootEntry_->addChild(createEntry(type, id, checker));
632         entry           = lastFoundEntry_->get();
633     }
634     return entry;
635 }
636
637 ::testing::AssertionResult TestReferenceChecker::Impl::processItem(const char* type,
638                                                                    const char* id,
639                                                                    const IReferenceDataEntryChecker& checker)
640 {
641     if (shouldIgnore())
642     {
643         return ::testing::AssertionSuccess();
644     }
645     std::string         fullId = appendPath(id);
646     ReferenceDataEntry* entry  = findOrCreateEntry(type, id, checker);
647     if (entry == nullptr)
648     {
649         return ::testing::AssertionFailure() << "Reference data item " << fullId << " not found";
650     }
651     entry->setChecked();
652     ::testing::AssertionResult result(checkEntry(*entry, fullId, type, checker));
653     if (outputRootEntry_ != nullptr && entry->correspondingOutputEntry() == nullptr)
654     {
655         if (!updateMismatchingEntries_ || result)
656         {
657             outputRootEntry_->addChild(entry->cloneToOutputEntry());
658         }
659         else
660         {
661             ReferenceDataEntry::EntryPointer outputEntry(createEntry(type, id, checker));
662             entry->setCorrespondingOutputEntry(outputEntry.get());
663             outputRootEntry_->addChild(move(outputEntry));
664             return ::testing::AssertionSuccess();
665         }
666     }
667     if (bSelfTestMode_ && !result)
668     {
669         ReferenceDataEntry expected(type, id);
670         checker.fillEntry(&expected);
671         result << std::endl
672                << "String value: '" << expected.value() << "'" << std::endl
673                << " Ref. string: '" << entry->value() << "'";
674     }
675     return result;
676 }
677
678
679 /********************************************************************
680  * TestReferenceData
681  */
682
683 TestReferenceData::TestReferenceData() : impl_(initReferenceDataInstance()) {}
684
685
686 TestReferenceData::TestReferenceData(ReferenceDataMode mode) :
687     impl_(initReferenceDataInstanceForSelfTest(mode))
688 {
689 }
690
691
692 TestReferenceData::~TestReferenceData() {}
693
694
695 TestReferenceChecker TestReferenceData::rootChecker()
696 {
697     if (!impl_->bInUse_ && !impl_->compareRootEntry_)
698     {
699         ADD_FAILURE() << "Reference data file not found: " << impl_->fullFilename_;
700     }
701     impl_->bInUse_ = true;
702     if (!impl_->compareRootEntry_)
703     {
704         return TestReferenceChecker(new TestReferenceChecker::Impl(true));
705     }
706     impl_->compareRootEntry_->setChecked();
707     return TestReferenceChecker(new TestReferenceChecker::Impl("",
708                                                                impl_->compareRootEntry_.get(),
709                                                                impl_->outputRootEntry_.get(),
710                                                                impl_->updateMismatchingEntries_,
711                                                                impl_->bSelfTestMode_,
712                                                                defaultRealTolerance()));
713 }
714
715
716 /********************************************************************
717  * TestReferenceChecker
718  */
719
720 TestReferenceChecker::TestReferenceChecker() : impl_(new Impl(false)) {}
721
722 TestReferenceChecker::TestReferenceChecker(Impl* impl) : impl_(impl) {}
723
724 TestReferenceChecker::TestReferenceChecker(const TestReferenceChecker& other) :
725     impl_(new Impl(*other.impl_))
726 {
727 }
728
729 TestReferenceChecker::TestReferenceChecker(TestReferenceChecker&& other) noexcept :
730     impl_(std::move(other.impl_))
731 {
732 }
733
734 TestReferenceChecker& TestReferenceChecker::operator=(TestReferenceChecker&& other) noexcept
735 {
736     impl_ = std::move(other.impl_);
737     return *this;
738 }
739
740 TestReferenceChecker::~TestReferenceChecker() {}
741
742 bool TestReferenceChecker::isValid() const
743 {
744     return impl_->initialized();
745 }
746
747
748 void TestReferenceChecker::setDefaultTolerance(const FloatingPointTolerance& tolerance)
749 {
750     impl_->defaultTolerance_ = tolerance;
751 }
752
753
754 void TestReferenceChecker::checkUnusedEntries()
755 {
756     if (impl_->compareRootEntry_)
757     {
758         gmx::test::checkUnusedEntries(*impl_->compareRootEntry_, impl_->path_);
759         // Mark them checked so that they are reported only once.
760         impl_->compareRootEntry_->setCheckedIncludingChildren();
761     }
762 }
763
764 void TestReferenceChecker::disableUnusedEntriesCheck()
765 {
766     if (impl_->compareRootEntry_)
767     {
768         impl_->compareRootEntry_->setCheckedIncludingChildren();
769     }
770 }
771
772
773 bool TestReferenceChecker::checkPresent(bool bPresent, const char* id)
774 {
775     if (impl_->shouldIgnore() || impl_->outputRootEntry_ != nullptr)
776     {
777         return bPresent;
778     }
779     ReferenceDataEntry::ChildIterator entry =
780             impl_->compareRootEntry_->findChild(id, impl_->lastFoundEntry_);
781     const bool bFound = impl_->compareRootEntry_->isValidChild(entry);
782     if (bFound != bPresent)
783     {
784         ADD_FAILURE() << "Mismatch while checking reference data item '" << impl_->appendPath(id) << "'\n"
785                       << "Expected: " << (bPresent ? "it is present.\n" : "it is absent.\n")
786                       << "  Actual: " << (bFound ? "it is present." : "it is absent.");
787     }
788     if (bFound && bPresent)
789     {
790         impl_->lastFoundEntry_ = entry;
791         return true;
792     }
793     return false;
794 }
795
796
797 TestReferenceChecker TestReferenceChecker::checkCompound(const char* type, const char* id)
798 {
799     if (impl_->shouldIgnore())
800     {
801         return TestReferenceChecker(new Impl(true));
802     }
803     std::string         fullId = impl_->appendPath(id);
804     NullChecker         checker;
805     ReferenceDataEntry* entry = impl_->findOrCreateEntry(type, id, checker);
806     if (entry == nullptr)
807     {
808         ADD_FAILURE() << "Reference data item " << fullId << " not found";
809         return TestReferenceChecker(new Impl(true));
810     }
811     entry->setChecked();
812     if (impl_->updateMismatchingEntries_)
813     {
814         entry->makeCompound(type);
815     }
816     else
817     {
818         ::testing::AssertionResult result(impl_->checkEntry(*entry, fullId, type, checker));
819         EXPECT_PLAIN(result);
820         if (!result)
821         {
822             return TestReferenceChecker(new Impl(true));
823         }
824     }
825     if (impl_->outputRootEntry_ != nullptr && entry->correspondingOutputEntry() == nullptr)
826     {
827         impl_->outputRootEntry_->addChild(entry->cloneToOutputEntry());
828     }
829     return TestReferenceChecker(new Impl(fullId,
830                                          entry,
831                                          entry->correspondingOutputEntry(),
832                                          impl_->updateMismatchingEntries_,
833                                          impl_->bSelfTestMode_,
834                                          impl_->defaultTolerance_));
835 }
836
837 TestReferenceChecker TestReferenceChecker::checkCompound(const char* type, const std::string& id)
838 {
839     return checkCompound(type, id.c_str());
840 }
841
842 /*! \brief Throw a TestException if the caller tries to write particular refdata that can't work.
843  *
844  * If the string to write is non-empty and has only whitespace,
845  * TinyXML2 can't read it correctly, so throw an exception for this
846  * case, so that we can't accidentally use it and run into mysterious
847  * problems.
848  *
849  * \todo Eliminate this limitation of TinyXML2. See
850  * e.g. https://github.com/leethomason/tinyxml2/issues/432
851  */
852 static void throwIfNonEmptyAndOnlyWhitespace(const std::string& s, const char* id)
853 {
854     if (!s.empty() && std::all_of(s.cbegin(), s.cend(), [](const char& c) { return std::isspace(c); }))
855     {
856         std::string message("String '" + s + "' with ");
857         message += (id != nullptr) ? "null " : "";
858         message += "ID ";
859         message += (id != nullptr) ? "" : id;
860         message +=
861                 " cannot be handled. We must refuse to write a refdata String"
862                 "field for a non-empty string that contains only whitespace, "
863                 "because it will not be read correctly by TinyXML2.";
864         GMX_THROW(TestException(message));
865     }
866 }
867
868 void TestReferenceChecker::checkBoolean(bool value, const char* id)
869 {
870     EXPECT_PLAIN(impl_->processItem(
871             Impl::cBooleanNodeName, id, ExactStringChecker(value ? "true" : "false")));
872 }
873
874
875 void TestReferenceChecker::checkString(const char* value, const char* id)
876 {
877     throwIfNonEmptyAndOnlyWhitespace(value, id);
878     EXPECT_PLAIN(impl_->processItem(Impl::cStringNodeName, id, ExactStringChecker(value)));
879 }
880
881
882 void TestReferenceChecker::checkString(const std::string& value, const char* id)
883 {
884     throwIfNonEmptyAndOnlyWhitespace(value, id);
885     EXPECT_PLAIN(impl_->processItem(Impl::cStringNodeName, id, ExactStringChecker(value)));
886 }
887
888
889 void TestReferenceChecker::checkTextBlock(const std::string& value, const char* id)
890 {
891     EXPECT_PLAIN(impl_->processItem(Impl::cStringNodeName, id, ExactStringBlockChecker(value)));
892 }
893
894
895 void TestReferenceChecker::checkUChar(unsigned char value, const char* id)
896 {
897     EXPECT_PLAIN(impl_->processItem(
898             Impl::cUCharNodeName, id, ExactStringChecker(formatString("%d", value))));
899 }
900
901 void TestReferenceChecker::checkInteger(int value, const char* id)
902 {
903     EXPECT_PLAIN(impl_->processItem(
904             Impl::cIntegerNodeName, id, ExactStringChecker(formatString("%d", value))));
905 }
906
907 void TestReferenceChecker::checkInt32(int32_t value, const char* id)
908 {
909     EXPECT_PLAIN(impl_->processItem(
910             Impl::cInt32NodeName, id, ExactStringChecker(formatString("%" PRId32, value))));
911 }
912
913 void TestReferenceChecker::checkUInt32(uint32_t value, const char* id)
914 {
915     EXPECT_PLAIN(impl_->processItem(
916             Impl::cUInt32NodeName, id, ExactStringChecker(formatString("%" PRIu32, value))));
917 }
918
919 void TestReferenceChecker::checkInt64(int64_t value, const char* id)
920 {
921     EXPECT_PLAIN(impl_->processItem(
922             Impl::cInt64NodeName, id, ExactStringChecker(formatString("%" PRId64, value))));
923 }
924
925 void TestReferenceChecker::checkUInt64(uint64_t value, const char* id)
926 {
927     EXPECT_PLAIN(impl_->processItem(
928             Impl::cUInt64NodeName, id, ExactStringChecker(formatString("%" PRIu64, value))));
929 }
930
931 void TestReferenceChecker::checkDouble(double value, const char* id)
932 {
933     FloatingPointChecker<double> checker(value, impl_->defaultTolerance_);
934     EXPECT_PLAIN(impl_->processItem(Impl::cRealNodeName, id, checker));
935 }
936
937
938 void TestReferenceChecker::checkFloat(float value, const char* id)
939 {
940     FloatingPointChecker<float> checker(value, impl_->defaultTolerance_);
941     EXPECT_PLAIN(impl_->processItem(Impl::cRealNodeName, id, checker));
942 }
943
944
945 void TestReferenceChecker::checkReal(float value, const char* id)
946 {
947     checkFloat(value, id);
948 }
949
950
951 void TestReferenceChecker::checkReal(double value, const char* id)
952 {
953     checkDouble(value, id);
954 }
955
956
957 void TestReferenceChecker::checkRealFromString(const std::string& value, const char* id)
958 {
959     FloatingPointFromStringChecker<real> checker(value, impl_->defaultTolerance_);
960     EXPECT_PLAIN(impl_->processItem(Impl::cRealNodeName, id, checker));
961 }
962
963
964 void TestReferenceChecker::checkVector(const BasicVector<int>& value, const char* id)
965 {
966     TestReferenceChecker compound(checkCompound(Impl::cVectorType, id));
967     compound.checkInteger(value[0], "X");
968     compound.checkInteger(value[1], "Y");
969     compound.checkInteger(value[2], "Z");
970 }
971
972
973 void TestReferenceChecker::checkVector(const BasicVector<float>& value, const char* id)
974 {
975     TestReferenceChecker compound(checkCompound(Impl::cVectorType, id));
976     compound.checkReal(value[0], "X");
977     compound.checkReal(value[1], "Y");
978     compound.checkReal(value[2], "Z");
979 }
980
981
982 void TestReferenceChecker::checkVector(const BasicVector<double>& value, const char* id)
983 {
984     TestReferenceChecker compound(checkCompound(Impl::cVectorType, id));
985     compound.checkReal(value[0], "X");
986     compound.checkReal(value[1], "Y");
987     compound.checkReal(value[2], "Z");
988 }
989
990
991 void TestReferenceChecker::checkVector(const int value[3], const char* id)
992 {
993     checkVector(BasicVector<int>(value), id);
994 }
995
996
997 void TestReferenceChecker::checkVector(const float value[3], const char* id)
998 {
999     checkVector(BasicVector<float>(value), id);
1000 }
1001
1002
1003 void TestReferenceChecker::checkVector(const double value[3], const char* id)
1004 {
1005     checkVector(BasicVector<double>(value), id);
1006 }
1007
1008
1009 void TestReferenceChecker::checkAny(const Any& any, const char* id)
1010 {
1011     if (any.isType<bool>())
1012     {
1013         checkBoolean(any.cast<bool>(), id);
1014     }
1015     else if (any.isType<int>())
1016     {
1017         checkInteger(any.cast<int>(), id);
1018     }
1019     else if (any.isType<int32_t>())
1020     {
1021         checkInt32(any.cast<int32_t>(), id);
1022     }
1023     else if (any.isType<uint32_t>())
1024     {
1025         checkInt32(any.cast<uint32_t>(), id);
1026     }
1027     else if (any.isType<int64_t>())
1028     {
1029         checkInt64(any.cast<int64_t>(), id);
1030     }
1031     else if (any.isType<uint64_t>())
1032     {
1033         checkInt64(any.cast<uint64_t>(), id);
1034     }
1035     else if (any.isType<float>())
1036     {
1037         checkFloat(any.cast<float>(), id);
1038     }
1039     else if (any.isType<double>())
1040     {
1041         checkDouble(any.cast<double>(), id);
1042     }
1043     else if (any.isType<std::string>())
1044     {
1045         checkString(any.cast<std::string>(), id);
1046     }
1047     else
1048     {
1049         GMX_THROW(TestException("Unsupported any type"));
1050     }
1051 }
1052
1053
1054 void TestReferenceChecker::checkKeyValueTreeObject(const KeyValueTreeObject& tree, const char* id)
1055 {
1056     TestReferenceChecker compound(checkCompound(Impl::cObjectType, id));
1057     for (const auto& prop : tree.properties())
1058     {
1059         compound.checkKeyValueTreeValue(prop.value(), prop.key().c_str());
1060     }
1061     compound.checkUnusedEntries();
1062 }
1063
1064
1065 void TestReferenceChecker::checkKeyValueTreeValue(const KeyValueTreeValue& value, const char* id)
1066 {
1067     if (value.isObject())
1068     {
1069         checkKeyValueTreeObject(value.asObject(), id);
1070     }
1071     else if (value.isArray())
1072     {
1073         const auto& values = value.asArray().values();
1074         checkSequence(values.begin(), values.end(), id);
1075     }
1076     else
1077     {
1078         checkAny(value.asAny(), id);
1079     }
1080 }
1081
1082
1083 TestReferenceChecker TestReferenceChecker::checkSequenceCompound(const char* id, size_t length)
1084 {
1085     TestReferenceChecker compound(checkCompound(Impl::cSequenceType, id));
1086     compound.checkInteger(static_cast<int>(length), Impl::cSequenceLengthName);
1087     return compound;
1088 }
1089
1090
1091 unsigned char TestReferenceChecker::readUChar(const char* id)
1092 {
1093     if (impl_->shouldIgnore())
1094     {
1095         GMX_THROW(TestException("Trying to read from non-existent reference data value"));
1096     }
1097     int value = 0;
1098     EXPECT_PLAIN(impl_->processItem(Impl::cUCharNodeName, id, ValueExtractor<int>(&value)));
1099     return value;
1100 }
1101
1102
1103 int TestReferenceChecker::readInteger(const char* id)
1104 {
1105     if (impl_->shouldIgnore())
1106     {
1107         GMX_THROW(TestException("Trying to read from non-existent reference data value"));
1108     }
1109     int value = 0;
1110     EXPECT_PLAIN(impl_->processItem(Impl::cIntegerNodeName, id, ValueExtractor<int>(&value)));
1111     return value;
1112 }
1113
1114
1115 int32_t TestReferenceChecker::readInt32(const char* id)
1116 {
1117     if (impl_->shouldIgnore())
1118     {
1119         GMX_THROW(TestException("Trying to read from non-existent reference data value"));
1120     }
1121     int32_t value = 0;
1122     EXPECT_PLAIN(impl_->processItem(Impl::cInt32NodeName, id, ValueExtractor<int32_t>(&value)));
1123     return value;
1124 }
1125
1126
1127 int64_t TestReferenceChecker::readInt64(const char* id)
1128 {
1129     if (impl_->shouldIgnore())
1130     {
1131         GMX_THROW(TestException("Trying to read from non-existent reference data value"));
1132     }
1133     int64_t value = 0;
1134     EXPECT_PLAIN(impl_->processItem(Impl::cInt64NodeName, id, ValueExtractor<int64_t>(&value)));
1135     return value;
1136 }
1137
1138
1139 float TestReferenceChecker::readFloat(const char* id)
1140 {
1141     if (impl_->shouldIgnore())
1142     {
1143         GMX_THROW(TestException("Trying to read from non-existent reference data value"));
1144     }
1145     float value = 0;
1146     EXPECT_PLAIN(impl_->processItem(Impl::cRealNodeName, id, ValueExtractor<float>(&value)));
1147     return value;
1148 }
1149
1150
1151 double TestReferenceChecker::readDouble(const char* id)
1152 {
1153     if (impl_->shouldIgnore())
1154     {
1155         GMX_THROW(TestException("Trying to read from non-existent reference data value"));
1156     }
1157     double value = 0;
1158     EXPECT_PLAIN(impl_->processItem(Impl::cRealNodeName, id, ValueExtractor<double>(&value)));
1159     return value;
1160 }
1161
1162
1163 std::string TestReferenceChecker::readString(const char* id)
1164 {
1165     if (impl_->shouldIgnore())
1166     {
1167         GMX_THROW(TestException("Trying to read from non-existent reference data value"));
1168     }
1169     std::string value;
1170     EXPECT_PLAIN(impl_->processItem(Impl::cStringNodeName, id, ValueExtractor<std::string>(&value)));
1171     return value;
1172 }
1173
1174 } // namespace test
1175 } // namespace gmx