075a028761ae41b642f1b492fe090530c245c332
[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,2021, by the GROMACS development team, led by
6  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
7  * and including many others, as listed in the AUTHORS file in the
8  * top-level source directory and at http://www.gromacs.org.
9  *
10  * GROMACS is free software; you can redistribute it and/or
11  * modify it under the terms of the GNU Lesser General Public License
12  * as published by the Free Software Foundation; either version 2.1
13  * of the License, or (at your option) any later version.
14  *
15  * GROMACS is distributed in the hope that it will be useful,
16  * but WITHOUT ANY WARRANTY; without even the implied warranty of
17  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
18  * Lesser General Public License for more details.
19  *
20  * You should have received a copy of the GNU Lesser General Public
21  * License along with GROMACS; if not, see
22  * http://www.gnu.org/licenses, or write to the Free Software Foundation,
23  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
24  *
25  * If you want to redistribute modifications to GROMACS, please
26  * consider that scientific software is very special. Version
27  * control is crucial - bugs must be traceable. We will be happy to
28  * consider code for inclusion in the official distribution, but
29  * derived work must not be called official GROMACS. Details are found
30  * in the README & COPYING files - if they are missing, get the
31  * official version at http://www.gromacs.org.
32  *
33  * To help us fund GROMACS development, we humbly ask that you cite
34  * the research papers on the package. Check out http://www.gromacs.org.
35  */
36 /*! \internal \file
37  * \brief
38  * Implements classes 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) const;
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), bSelfTestMode_(bSelfTestMode), bInUse_(false)
292 {
293     const std::string dirname  = bSelfTestMode ? TestFileManager::getGlobalOutputTempDirectory()
294                                                : TestFileManager::getInputDataDirectory();
295     const std::string filename = TestFileManager::getTestSpecificFileName(".xml");
296     fullFilename_              = Path::join(dirname, "refdata", filename);
297
298     switch (mode)
299     {
300         case ReferenceDataMode::Compare:
301             if (File::exists(fullFilename_, File::throwOnError))
302             {
303                 compareRootEntry_ = readReferenceDataFile(fullFilename_);
304             }
305             break;
306         case ReferenceDataMode::CreateMissing:
307             if (File::exists(fullFilename_, File::throwOnError))
308             {
309                 compareRootEntry_ = readReferenceDataFile(fullFilename_);
310             }
311             else
312             {
313                 compareRootEntry_ = ReferenceDataEntry::createRoot();
314                 outputRootEntry_  = ReferenceDataEntry::createRoot();
315             }
316             break;
317         case ReferenceDataMode::UpdateChanged:
318             if (File::exists(fullFilename_, File::throwOnError))
319             {
320                 compareRootEntry_ = readReferenceDataFile(fullFilename_);
321             }
322             else
323             {
324                 compareRootEntry_ = ReferenceDataEntry::createRoot();
325             }
326             outputRootEntry_          = ReferenceDataEntry::createRoot();
327             updateMismatchingEntries_ = true;
328             break;
329         case ReferenceDataMode::UpdateAll:
330             compareRootEntry_ = ReferenceDataEntry::createRoot();
331             outputRootEntry_  = ReferenceDataEntry::createRoot();
332             break;
333         case ReferenceDataMode::Count: GMX_THROW(InternalError("Invalid reference data mode"));
334     }
335 }
336
337 void TestReferenceDataImpl::onTestEnd(bool testPassed) const
338 {
339     if (!bInUse_)
340     {
341         return;
342     }
343     // TODO: Only write the file with update-changed if there were actual changes.
344     if (outputRootEntry_)
345     {
346         if (testPassed)
347         {
348             std::string dirname = Path::getParentPath(fullFilename_);
349             if (!Directory::exists(dirname))
350             {
351                 if (Directory::create(dirname) != 0)
352                 {
353                     GMX_THROW(TestException("Creation of reference data directory failed: " + dirname));
354                 }
355             }
356             writeReferenceDataFile(fullFilename_, *outputRootEntry_);
357         }
358     }
359     else if (compareRootEntry_)
360     {
361         checkUnusedEntries(*compareRootEntry_, "");
362     }
363 }
364
365 } // namespace internal
366
367
368 /********************************************************************
369  * TestReferenceChecker::Impl
370  */
371
372 /*! \internal \brief
373  * Private implementation class for TestReferenceChecker.
374  *
375  * \ingroup module_testutils
376  */
377 class TestReferenceChecker::Impl
378 {
379 public:
380     //! String constant for naming XML elements for boolean values.
381     static const char* const cBooleanNodeName;
382     //! String constant for naming XML elements for string values.
383     static const char* const cStringNodeName;
384     //! String constant for naming XML elements for unsigned char values.
385     static const char* const cUCharNodeName;
386     //! String constant for naming XML elements for integer values.
387     static const char* const cIntegerNodeName;
388     //! String constant for naming XML elements for int32 values.
389     static const char* const cInt32NodeName;
390     //! String constant for naming XML elements for unsigned int32 values.
391     static const char* const cUInt32NodeName;
392     //! String constant for naming XML elements for int32 values.
393     static const char* const cInt64NodeName;
394     //! String constant for naming XML elements for unsigned int64 values.
395     static const char* const cUInt64NodeName;
396     //! String constant for naming XML elements for floating-point values.
397     static const char* const cRealNodeName;
398     //! String constant for naming XML attribute for value identifiers.
399     static const char* const cIdAttrName;
400     //! String constant for naming compounds for vectors.
401     static const char* const cVectorType;
402     //! String constant for naming compounds for key-value tree objects.
403     static const char* const cObjectType;
404     //! String constant for naming compounds for sequences.
405     static const char* const cSequenceType;
406     //! String constant for value identifier for sequence length.
407     static const char* const cSequenceLengthName;
408
409     //! Creates a checker that does nothing.
410     explicit Impl(bool initialized);
411     //! Creates a checker with a given root entry.
412     Impl(const std::string&            path,
413          ReferenceDataEntry*           compareRootEntry,
414          ReferenceDataEntry*           outputRootEntry,
415          bool                          updateMismatchingEntries,
416          bool                          bSelfTestMode,
417          const FloatingPointTolerance& defaultTolerance);
418
419     //! Returns the path of this checker with \p id appended.
420     std::string appendPath(const char* id) const;
421
422     //! Creates an entry with given parameters and fills it with \p checker.
423     static ReferenceDataEntry::EntryPointer createEntry(const char*                       type,
424                                                         const char*                       id,
425                                                         const IReferenceDataEntryChecker& checker)
426     {
427         ReferenceDataEntry::EntryPointer entry(new ReferenceDataEntry(type, id));
428         checker.fillEntry(entry.get());
429         return entry;
430     }
431     //! Checks an entry for correct type and using \p checker.
432     static ::testing::AssertionResult checkEntry(const ReferenceDataEntry&         entry,
433                                                  const std::string&                fullId,
434                                                  const char*                       type,
435                                                  const IReferenceDataEntryChecker& checker)
436     {
437         if (entry.type() != type)
438         {
439             return ::testing::AssertionFailure() << "Mismatching reference data item type" << std::endl
440                                                  << "  In item: " << fullId << std::endl
441                                                  << "   Actual: " << type << std::endl
442                                                  << "Reference: " << entry.type();
443         }
444         return checker.checkEntry(entry, fullId);
445     }
446     //! Finds an entry by id and updates the last found entry pointer.
447     ReferenceDataEntry* findEntry(const char* id);
448     /*! \brief
449      * Finds/creates a reference data entry to match against.
450      *
451      * \param[in]  type   Type of entry to create.
452      * \param[in]  id     Unique identifier of the entry (can be NULL, in
453      *      which case the next entry without an id is matched).
454      * \param[out] checker  Checker to use for filling out created entries.
455      * \returns    Matching entry, or NULL if no matching entry found
456      *      (NULL is never returned in write mode; new entries are created
457      *      instead).
458      */
459     ReferenceDataEntry* findOrCreateEntry(const char*                       type,
460                                           const char*                       id,
461                                           const IReferenceDataEntryChecker& checker);
462     /*! \brief
463      * Helper method for checking a reference data value.
464      *
465      * \param[in]  name   Type of entry to find.
466      * \param[in]  id     Unique identifier of the entry (can be NULL, in
467      *     which case the next entry without an id is matched).
468      * \param[in]  checker  Checker that provides logic specific to the
469      *     type of the entry.
470      * \returns    Whether the reference data matched, including details
471      *     of the mismatch if the comparison failed.
472      * \throws     TestException if there is a problem parsing the
473      *     reference data.
474      *
475      * Performs common tasks in checking a reference value, such as
476      * finding or creating the correct entry.
477      * Caller needs to provide a checker object that provides the string
478      * value for a newly created entry and performs the actual comparison
479      * against a found entry.
480      */
481     ::testing::AssertionResult processItem(const char*                       name,
482                                            const char*                       id,
483                                            const IReferenceDataEntryChecker& checker);
484     /*! \brief
485      * Whether the checker is initialized.
486      */
487     bool initialized() const { return initialized_; }
488     /*! \brief
489      * Whether the checker should ignore all validation calls.
490      *
491      * This is used to ignore any calls within compounds for which
492      * reference data could not be found, such that only one error is
493      * issued for the missing compound, instead of every individual value.
494      */
495     bool shouldIgnore() const
496     {
497         GMX_RELEASE_ASSERT(initialized(), "Accessing uninitialized reference data checker.");
498         return compareRootEntry_ == nullptr;
499     }
500
501     //! Whether initialized with other means than the default constructor.
502     bool initialized_;
503     //! Default floating-point comparison tolerance.
504     FloatingPointTolerance defaultTolerance_;
505     /*! \brief
506      * Human-readable path to the root node of this checker.
507      *
508      * For the root checker, this will be "/", and for each compound, the
509      * id of the compound is added.  Used for reporting comparison
510      * mismatches.
511      */
512     std::string path_;
513     /*! \brief
514      * Current entry under which reference data is searched for comparison.
515      *
516      * Points to either the TestReferenceDataImpl::compareRootEntry_, or to
517      * a compound entry in the tree rooted at that entry.
518      *
519      * Can be NULL, in which case this checker does nothing (doesn't even
520      * report errors, see shouldIgnore()).
521      */
522     ReferenceDataEntry* compareRootEntry_;
523     /*! \brief
524      * Current entry under which entries for writing are created.
525      *
526      * Points to either the TestReferenceDataImpl::outputRootEntry_, or to
527      * a compound entry in the tree rooted at that entry.  NULL if only
528      * comparing, or if shouldIgnore() returns `false`.
529      */
530     ReferenceDataEntry* outputRootEntry_;
531     /*! \brief
532      * Iterator to a child of \a compareRootEntry_ that was last found.
533      *
534      * If `compareRootEntry_->isValidChild()` returns false, no entry has
535      * been found yet.
536      * After every check, is updated to point to the entry that was used
537      * for the check.
538      * Subsequent checks start the search for the matching node on this
539      * node.
540      */
541     ReferenceDataEntry::ChildIterator lastFoundEntry_;
542     /*! \brief
543      * Whether the reference data is being written (true) or compared
544      * (false).
545      */
546     bool updateMismatchingEntries_;
547     //! `true` if self-testing (enables extra failure messages).
548     bool bSelfTestMode_;
549     /*! \brief
550      * Current number of unnamed elements in a sequence.
551      *
552      * It is the index of the current unnamed element.
553      */
554     int seqIndex_;
555 };
556
557 const char* const TestReferenceChecker::Impl::cBooleanNodeName    = "Bool";
558 const char* const TestReferenceChecker::Impl::cStringNodeName     = "String";
559 const char* const TestReferenceChecker::Impl::cUCharNodeName      = "UChar";
560 const char* const TestReferenceChecker::Impl::cIntegerNodeName    = "Int";
561 const char* const TestReferenceChecker::Impl::cInt32NodeName      = "Int32";
562 const char* const TestReferenceChecker::Impl::cUInt32NodeName     = "UInt32";
563 const char* const TestReferenceChecker::Impl::cInt64NodeName      = "Int64";
564 const char* const TestReferenceChecker::Impl::cUInt64NodeName     = "UInt64";
565 const char* const TestReferenceChecker::Impl::cRealNodeName       = "Real";
566 const char* const TestReferenceChecker::Impl::cIdAttrName         = "Name";
567 const char* const TestReferenceChecker::Impl::cVectorType         = "Vector";
568 const char* const TestReferenceChecker::Impl::cObjectType         = "Object";
569 const char* const TestReferenceChecker::Impl::cSequenceType       = "Sequence";
570 const char* const TestReferenceChecker::Impl::cSequenceLengthName = "Length";
571
572
573 TestReferenceChecker::Impl::Impl(bool initialized) :
574     initialized_(initialized),
575     defaultTolerance_(defaultRealTolerance()),
576     compareRootEntry_(nullptr),
577     outputRootEntry_(nullptr),
578     updateMismatchingEntries_(false),
579     bSelfTestMode_(false),
580     seqIndex_(-1)
581 {
582 }
583
584
585 TestReferenceChecker::Impl::Impl(const std::string&            path,
586                                  ReferenceDataEntry*           compareRootEntry,
587                                  ReferenceDataEntry*           outputRootEntry,
588                                  bool                          updateMismatchingEntries,
589                                  bool                          bSelfTestMode,
590                                  const FloatingPointTolerance& defaultTolerance) :
591     initialized_(true),
592     defaultTolerance_(defaultTolerance),
593     path_(path),
594     compareRootEntry_(compareRootEntry),
595     outputRootEntry_(outputRootEntry),
596     lastFoundEntry_(compareRootEntry->children().end()),
597     updateMismatchingEntries_(updateMismatchingEntries),
598     bSelfTestMode_(bSelfTestMode),
599     seqIndex_(-1)
600 {
601 }
602
603
604 std::string TestReferenceChecker::Impl::appendPath(const char* id) const
605 {
606     return id != nullptr ? formatEntryPath(path_, id) : formatSequenceEntryPath(path_, seqIndex_);
607 }
608
609
610 ReferenceDataEntry* TestReferenceChecker::Impl::findEntry(const char* id)
611 {
612     ReferenceDataEntry::ChildIterator entry = compareRootEntry_->findChild(id, lastFoundEntry_);
613     seqIndex_                               = (id == nullptr) ? seqIndex_ + 1 : -1;
614     if (compareRootEntry_->isValidChild(entry))
615     {
616         lastFoundEntry_ = entry;
617         return entry->get();
618     }
619     return nullptr;
620 }
621
622 ReferenceDataEntry* TestReferenceChecker::Impl::findOrCreateEntry(const char* type,
623                                                                   const char* id,
624                                                                   const IReferenceDataEntryChecker& checker)
625 {
626     ReferenceDataEntry* entry = findEntry(id);
627     if (entry == nullptr && outputRootEntry_ != nullptr)
628     {
629         lastFoundEntry_ = compareRootEntry_->addChild(createEntry(type, id, checker));
630         entry           = lastFoundEntry_->get();
631     }
632     return entry;
633 }
634
635 ::testing::AssertionResult TestReferenceChecker::Impl::processItem(const char* type,
636                                                                    const char* id,
637                                                                    const IReferenceDataEntryChecker& checker)
638 {
639     if (shouldIgnore())
640     {
641         return ::testing::AssertionSuccess();
642     }
643     std::string         fullId = appendPath(id);
644     ReferenceDataEntry* entry  = findOrCreateEntry(type, id, checker);
645     if (entry == nullptr)
646     {
647         return ::testing::AssertionFailure() << "Reference data item " << fullId << " not found";
648     }
649     entry->setChecked();
650     ::testing::AssertionResult result(checkEntry(*entry, fullId, type, checker));
651     if (outputRootEntry_ != nullptr && entry->correspondingOutputEntry() == nullptr)
652     {
653         if (!updateMismatchingEntries_ || result)
654         {
655             outputRootEntry_->addChild(entry->cloneToOutputEntry());
656         }
657         else
658         {
659             ReferenceDataEntry::EntryPointer outputEntry(createEntry(type, id, checker));
660             entry->setCorrespondingOutputEntry(outputEntry.get());
661             outputRootEntry_->addChild(move(outputEntry));
662             return ::testing::AssertionSuccess();
663         }
664     }
665     if (bSelfTestMode_ && !result)
666     {
667         ReferenceDataEntry expected(type, id);
668         checker.fillEntry(&expected);
669         result << std::endl
670                << "String value: '" << expected.value() << "'" << std::endl
671                << " Ref. string: '" << entry->value() << "'";
672     }
673     return result;
674 }
675
676
677 /********************************************************************
678  * TestReferenceData
679  */
680
681 TestReferenceData::TestReferenceData() : impl_(initReferenceDataInstance()) {}
682
683
684 TestReferenceData::TestReferenceData(ReferenceDataMode mode) :
685     impl_(initReferenceDataInstanceForSelfTest(mode))
686 {
687 }
688
689
690 TestReferenceData::~TestReferenceData() {}
691
692
693 TestReferenceChecker TestReferenceData::rootChecker()
694 {
695     if (!impl_->bInUse_ && !impl_->compareRootEntry_)
696     {
697         ADD_FAILURE() << "Reference data file not found: " << impl_->fullFilename_;
698     }
699     impl_->bInUse_ = true;
700     if (!impl_->compareRootEntry_)
701     {
702         return TestReferenceChecker(new TestReferenceChecker::Impl(true));
703     }
704     impl_->compareRootEntry_->setChecked();
705     return TestReferenceChecker(new TestReferenceChecker::Impl("",
706                                                                impl_->compareRootEntry_.get(),
707                                                                impl_->outputRootEntry_.get(),
708                                                                impl_->updateMismatchingEntries_,
709                                                                impl_->bSelfTestMode_,
710                                                                defaultRealTolerance()));
711 }
712
713
714 /********************************************************************
715  * TestReferenceChecker
716  */
717
718 TestReferenceChecker::TestReferenceChecker() : impl_(new Impl(false)) {}
719
720 TestReferenceChecker::TestReferenceChecker(Impl* impl) : impl_(impl) {}
721
722 TestReferenceChecker::TestReferenceChecker(const TestReferenceChecker& other) :
723     impl_(new Impl(*other.impl_))
724 {
725 }
726
727 TestReferenceChecker::TestReferenceChecker(TestReferenceChecker&& other) noexcept :
728     impl_(std::move(other.impl_))
729 {
730 }
731
732 TestReferenceChecker& TestReferenceChecker::operator=(TestReferenceChecker&& other) noexcept
733 {
734     impl_ = std::move(other.impl_);
735     return *this;
736 }
737
738 TestReferenceChecker::~TestReferenceChecker() {}
739
740 bool TestReferenceChecker::isValid() const
741 {
742     return impl_->initialized();
743 }
744
745
746 void TestReferenceChecker::setDefaultTolerance(const FloatingPointTolerance& tolerance)
747 {
748     impl_->defaultTolerance_ = tolerance;
749 }
750
751
752 void TestReferenceChecker::checkUnusedEntries()
753 {
754     if (impl_->compareRootEntry_)
755     {
756         gmx::test::checkUnusedEntries(*impl_->compareRootEntry_, impl_->path_);
757         // Mark them checked so that they are reported only once.
758         impl_->compareRootEntry_->setCheckedIncludingChildren();
759     }
760 }
761
762 void TestReferenceChecker::disableUnusedEntriesCheck()
763 {
764     if (impl_->compareRootEntry_)
765     {
766         impl_->compareRootEntry_->setCheckedIncludingChildren();
767     }
768 }
769
770
771 bool TestReferenceChecker::checkPresent(bool bPresent, const char* id)
772 {
773     if (impl_->shouldIgnore() || impl_->outputRootEntry_ != nullptr)
774     {
775         return bPresent;
776     }
777     ReferenceDataEntry::ChildIterator entry =
778             impl_->compareRootEntry_->findChild(id, impl_->lastFoundEntry_);
779     const bool bFound = impl_->compareRootEntry_->isValidChild(entry);
780     if (bFound != bPresent)
781     {
782         ADD_FAILURE() << "Mismatch while checking reference data item '" << impl_->appendPath(id) << "'\n"
783                       << "Expected: " << (bPresent ? "it is present.\n" : "it is absent.\n")
784                       << "  Actual: " << (bFound ? "it is present." : "it is absent.");
785     }
786     if (bFound && bPresent)
787     {
788         impl_->lastFoundEntry_ = entry;
789         return true;
790     }
791     return false;
792 }
793
794
795 TestReferenceChecker TestReferenceChecker::checkCompound(const char* type, const char* id)
796 {
797     if (impl_->shouldIgnore())
798     {
799         return TestReferenceChecker(new Impl(true));
800     }
801     std::string         fullId = impl_->appendPath(id);
802     NullChecker         checker;
803     ReferenceDataEntry* entry = impl_->findOrCreateEntry(type, id, checker);
804     if (entry == nullptr)
805     {
806         ADD_FAILURE() << "Reference data item " << fullId << " not found";
807         return TestReferenceChecker(new Impl(true));
808     }
809     entry->setChecked();
810     if (impl_->updateMismatchingEntries_)
811     {
812         entry->makeCompound(type);
813     }
814     else
815     {
816         ::testing::AssertionResult result(impl_->checkEntry(*entry, fullId, type, checker));
817         EXPECT_PLAIN(result);
818         if (!result)
819         {
820             return TestReferenceChecker(new Impl(true));
821         }
822     }
823     if (impl_->outputRootEntry_ != nullptr && entry->correspondingOutputEntry() == nullptr)
824     {
825         impl_->outputRootEntry_->addChild(entry->cloneToOutputEntry());
826     }
827     return TestReferenceChecker(new Impl(fullId,
828                                          entry,
829                                          entry->correspondingOutputEntry(),
830                                          impl_->updateMismatchingEntries_,
831                                          impl_->bSelfTestMode_,
832                                          impl_->defaultTolerance_));
833 }
834
835 TestReferenceChecker TestReferenceChecker::checkCompound(const char* type, const std::string& id)
836 {
837     return checkCompound(type, id.c_str());
838 }
839
840 /*! \brief Throw a TestException if the caller tries to write particular refdata that can't work.
841  *
842  * If the string to write is non-empty and has only whitespace,
843  * TinyXML2 can't read it correctly, so throw an exception for this
844  * case, so that we can't accidentally use it and run into mysterious
845  * problems.
846  *
847  * \todo Eliminate this limitation of TinyXML2. See
848  * e.g. https://github.com/leethomason/tinyxml2/issues/432
849  */
850 static void throwIfNonEmptyAndOnlyWhitespace(const std::string& s, const char* id)
851 {
852     if (!s.empty() && std::all_of(s.cbegin(), s.cend(), [](const char& c) { return std::isspace(c); }))
853     {
854         std::string message("String '" + s + "' with ");
855         message += (id != nullptr) ? "null " : "";
856         message += "ID ";
857         message += (id != nullptr) ? "" : id;
858         message +=
859                 " cannot be handled. We must refuse to write a refdata String"
860                 "field for a non-empty string that contains only whitespace, "
861                 "because it will not be read correctly by TinyXML2.";
862         GMX_THROW(TestException(message));
863     }
864 }
865
866 void TestReferenceChecker::checkBoolean(bool value, const char* id)
867 {
868     EXPECT_PLAIN(impl_->processItem(
869             Impl::cBooleanNodeName, id, ExactStringChecker(value ? "true" : "false")));
870 }
871
872
873 void TestReferenceChecker::checkString(const char* value, const char* id)
874 {
875     throwIfNonEmptyAndOnlyWhitespace(value, id);
876     EXPECT_PLAIN(impl_->processItem(Impl::cStringNodeName, id, ExactStringChecker(value)));
877 }
878
879
880 void TestReferenceChecker::checkString(const std::string& value, const char* id)
881 {
882     throwIfNonEmptyAndOnlyWhitespace(value, id);
883     EXPECT_PLAIN(impl_->processItem(Impl::cStringNodeName, id, ExactStringChecker(value)));
884 }
885
886
887 void TestReferenceChecker::checkTextBlock(const std::string& value, const char* id)
888 {
889     EXPECT_PLAIN(impl_->processItem(Impl::cStringNodeName, id, ExactStringBlockChecker(value)));
890 }
891
892
893 void TestReferenceChecker::checkUChar(unsigned char value, const char* id)
894 {
895     EXPECT_PLAIN(impl_->processItem(
896             Impl::cUCharNodeName, id, ExactStringChecker(formatString("%d", value))));
897 }
898
899 void TestReferenceChecker::checkInteger(int value, const char* id)
900 {
901     EXPECT_PLAIN(impl_->processItem(
902             Impl::cIntegerNodeName, id, ExactStringChecker(formatString("%d", value))));
903 }
904
905 void TestReferenceChecker::checkInt32(int32_t value, const char* id)
906 {
907     EXPECT_PLAIN(impl_->processItem(
908             Impl::cInt32NodeName, id, ExactStringChecker(formatString("%" PRId32, value))));
909 }
910
911 void TestReferenceChecker::checkUInt32(uint32_t value, const char* id)
912 {
913     EXPECT_PLAIN(impl_->processItem(
914             Impl::cUInt32NodeName, id, ExactStringChecker(formatString("%" PRIu32, value))));
915 }
916
917 void TestReferenceChecker::checkInt64(int64_t value, const char* id)
918 {
919     EXPECT_PLAIN(impl_->processItem(
920             Impl::cInt64NodeName, id, ExactStringChecker(formatString("%" PRId64, value))));
921 }
922
923 void TestReferenceChecker::checkUInt64(uint64_t value, const char* id)
924 {
925     EXPECT_PLAIN(impl_->processItem(
926             Impl::cUInt64NodeName, id, ExactStringChecker(formatString("%" PRIu64, value))));
927 }
928
929 void TestReferenceChecker::checkDouble(double value, const char* id)
930 {
931     FloatingPointChecker<double> checker(value, impl_->defaultTolerance_);
932     EXPECT_PLAIN(impl_->processItem(Impl::cRealNodeName, id, checker));
933 }
934
935
936 void TestReferenceChecker::checkFloat(float value, const char* id)
937 {
938     FloatingPointChecker<float> checker(value, impl_->defaultTolerance_);
939     EXPECT_PLAIN(impl_->processItem(Impl::cRealNodeName, id, checker));
940 }
941
942
943 void TestReferenceChecker::checkReal(float value, const char* id)
944 {
945     checkFloat(value, id);
946 }
947
948
949 void TestReferenceChecker::checkReal(double value, const char* id)
950 {
951     checkDouble(value, id);
952 }
953
954
955 void TestReferenceChecker::checkRealFromString(const std::string& value, const char* id)
956 {
957     FloatingPointFromStringChecker<real> checker(value, impl_->defaultTolerance_);
958     EXPECT_PLAIN(impl_->processItem(Impl::cRealNodeName, id, checker));
959 }
960
961
962 void TestReferenceChecker::checkVector(const BasicVector<int>& value, const char* id)
963 {
964     TestReferenceChecker compound(checkCompound(Impl::cVectorType, id));
965     compound.checkInteger(value[0], "X");
966     compound.checkInteger(value[1], "Y");
967     compound.checkInteger(value[2], "Z");
968 }
969
970
971 void TestReferenceChecker::checkVector(const BasicVector<float>& value, const char* id)
972 {
973     TestReferenceChecker compound(checkCompound(Impl::cVectorType, id));
974     compound.checkReal(value[0], "X");
975     compound.checkReal(value[1], "Y");
976     compound.checkReal(value[2], "Z");
977 }
978
979
980 void TestReferenceChecker::checkVector(const BasicVector<double>& value, const char* id)
981 {
982     TestReferenceChecker compound(checkCompound(Impl::cVectorType, id));
983     compound.checkReal(value[0], "X");
984     compound.checkReal(value[1], "Y");
985     compound.checkReal(value[2], "Z");
986 }
987
988
989 void TestReferenceChecker::checkVector(const int value[3], const char* id)
990 {
991     checkVector(BasicVector<int>(value), id);
992 }
993
994
995 void TestReferenceChecker::checkVector(const float value[3], const char* id)
996 {
997     checkVector(BasicVector<float>(value), id);
998 }
999
1000
1001 void TestReferenceChecker::checkVector(const double value[3], const char* id)
1002 {
1003     checkVector(BasicVector<double>(value), id);
1004 }
1005
1006
1007 void TestReferenceChecker::checkAny(const Any& any, const char* id)
1008 {
1009     if (any.isType<bool>())
1010     {
1011         checkBoolean(any.cast<bool>(), id);
1012     }
1013     else if (any.isType<int>())
1014     {
1015         checkInteger(any.cast<int>(), id);
1016     }
1017     else if (any.isType<int32_t>())
1018     {
1019         checkInt32(any.cast<int32_t>(), id);
1020     }
1021     else if (any.isType<uint32_t>())
1022     {
1023         checkInt32(any.cast<uint32_t>(), id);
1024     }
1025     else if (any.isType<int64_t>())
1026     {
1027         checkInt64(any.cast<int64_t>(), id);
1028     }
1029     else if (any.isType<uint64_t>())
1030     {
1031         checkInt64(any.cast<uint64_t>(), id);
1032     }
1033     else if (any.isType<float>())
1034     {
1035         checkFloat(any.cast<float>(), id);
1036     }
1037     else if (any.isType<double>())
1038     {
1039         checkDouble(any.cast<double>(), id);
1040     }
1041     else if (any.isType<std::string>())
1042     {
1043         checkString(any.cast<std::string>(), id);
1044     }
1045     else
1046     {
1047         GMX_THROW(TestException("Unsupported any type"));
1048     }
1049 }
1050
1051
1052 void TestReferenceChecker::checkKeyValueTreeObject(const KeyValueTreeObject& tree, const char* id)
1053 {
1054     TestReferenceChecker compound(checkCompound(Impl::cObjectType, id));
1055     for (const auto& prop : tree.properties())
1056     {
1057         compound.checkKeyValueTreeValue(prop.value(), prop.key().c_str());
1058     }
1059     compound.checkUnusedEntries();
1060 }
1061
1062
1063 void TestReferenceChecker::checkKeyValueTreeValue(const KeyValueTreeValue& value, const char* id)
1064 {
1065     if (value.isObject())
1066     {
1067         checkKeyValueTreeObject(value.asObject(), id);
1068     }
1069     else if (value.isArray())
1070     {
1071         const auto& values = value.asArray().values();
1072         checkSequence(values.begin(), values.end(), id);
1073     }
1074     else
1075     {
1076         checkAny(value.asAny(), id);
1077     }
1078 }
1079
1080
1081 TestReferenceChecker TestReferenceChecker::checkSequenceCompound(const char* id, size_t length)
1082 {
1083     TestReferenceChecker compound(checkCompound(Impl::cSequenceType, id));
1084     compound.checkInteger(static_cast<int>(length), Impl::cSequenceLengthName);
1085     return compound;
1086 }
1087
1088
1089 unsigned char TestReferenceChecker::readUChar(const char* id)
1090 {
1091     if (impl_->shouldIgnore())
1092     {
1093         GMX_THROW(TestException("Trying to read from non-existent reference data value"));
1094     }
1095     int value = 0;
1096     EXPECT_PLAIN(impl_->processItem(Impl::cUCharNodeName, id, ValueExtractor<int>(&value)));
1097     return value;
1098 }
1099
1100
1101 int TestReferenceChecker::readInteger(const char* id)
1102 {
1103     if (impl_->shouldIgnore())
1104     {
1105         GMX_THROW(TestException("Trying to read from non-existent reference data value"));
1106     }
1107     int value = 0;
1108     EXPECT_PLAIN(impl_->processItem(Impl::cIntegerNodeName, id, ValueExtractor<int>(&value)));
1109     return value;
1110 }
1111
1112
1113 int32_t TestReferenceChecker::readInt32(const char* id)
1114 {
1115     if (impl_->shouldIgnore())
1116     {
1117         GMX_THROW(TestException("Trying to read from non-existent reference data value"));
1118     }
1119     int32_t value = 0;
1120     EXPECT_PLAIN(impl_->processItem(Impl::cInt32NodeName, id, ValueExtractor<int32_t>(&value)));
1121     return value;
1122 }
1123
1124
1125 int64_t TestReferenceChecker::readInt64(const char* id)
1126 {
1127     if (impl_->shouldIgnore())
1128     {
1129         GMX_THROW(TestException("Trying to read from non-existent reference data value"));
1130     }
1131     int64_t value = 0;
1132     EXPECT_PLAIN(impl_->processItem(Impl::cInt64NodeName, id, ValueExtractor<int64_t>(&value)));
1133     return value;
1134 }
1135
1136
1137 float TestReferenceChecker::readFloat(const char* id)
1138 {
1139     if (impl_->shouldIgnore())
1140     {
1141         GMX_THROW(TestException("Trying to read from non-existent reference data value"));
1142     }
1143     float value = 0;
1144     EXPECT_PLAIN(impl_->processItem(Impl::cRealNodeName, id, ValueExtractor<float>(&value)));
1145     return value;
1146 }
1147
1148
1149 double TestReferenceChecker::readDouble(const char* id)
1150 {
1151     if (impl_->shouldIgnore())
1152     {
1153         GMX_THROW(TestException("Trying to read from non-existent reference data value"));
1154     }
1155     double value = 0;
1156     EXPECT_PLAIN(impl_->processItem(Impl::cRealNodeName, id, ValueExtractor<double>(&value)));
1157     return value;
1158 }
1159
1160
1161 std::string TestReferenceChecker::readString(const char* id)
1162 {
1163     if (impl_->shouldIgnore())
1164     {
1165         GMX_THROW(TestException("Trying to read from non-existent reference data value"));
1166     }
1167     std::string value;
1168     EXPECT_PLAIN(impl_->processItem(Impl::cStringNodeName, id, ValueExtractor<std::string>(&value)));
1169     return value;
1170 }
1171
1172 } // namespace test
1173 } // namespace gmx