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