Forked from
agnumpde / dune-fufem
362 commits behind the upstream repository.
-
Elias Pipping authoredElias Pipping authored
Code owners
Assign users and groups as approvers for specific file changes. Learn more.
test-hdf5.cc 5.85 KiB
#ifdef HAVE_CONFIG_H
#include "config.h"
#endif
#include <cassert>
#include <dune/common/fmatrix.hh>
#include <dune/common/fvector.hh>
#include <dune/common/parametertree.hh>
#include <dune/istl/bvector.hh>
#include <dune/fufem/hdf5/attributes.hh>
#include <dune/fufem/hdf5/file.hh>
#include <dune/fufem/hdf5/sequenceio.hh>
#include <dune/fufem/hdf5/singletonwriter.hh>
void deleteIfNecessary(std::string filename) {
struct stat buffer;
if (stat(filename.c_str(), &buffer) == 0) {
remove(filename.c_str());
assert(stat(filename.c_str(), &buffer) != 0);
}
}
template <class Vector>
bool checkMissingFile(Vector const &written0) {
bool passed = true;
using ft = typename Vector::block_type::field_type;
std::string const filename("no-such-file.h5");
try {
deleteIfNecessary(filename);
HDF5::File file(filename.c_str(), HDF5::Access::READONLY);
HDF5::SequenceIO<2, ft> h5writer(file, "data", 2, 3);
addEntry(h5writer, 0, written0);
} catch (Dune::Exception e) {
} // Expecting an exception here!
return passed;
}
template <class Vector>
bool checkMissingDataset(Vector const &written0) {
bool passed = true;
using ft = typename Vector::block_type::field_type;
std::string const filename("out-missing.h5");
{
deleteIfNecessary(filename);
HDF5::File file(filename.c_str());
HDF5::SequenceIO<2, ft> h5writer(file, "data", 2, 3);
addEntry(h5writer, 0, written0);
}
Vector read0;
try {
HDF5::File file(filename.c_str(), HDF5::Access::READONLY);
HDF5::SequenceIO<2, ft> h5reader(file, "dataX", 2, 3);
readEntry(h5reader, 0, read0);
passed = false;
} catch (Dune::Exception e) {
} // Expecting an exception here!
return passed;
}
template <class Vector>
bool checkAddSequence(Vector const &written0, Vector const &written1) {
bool passed = true;
using ft = typename Vector::block_type::field_type;
std::string const filename("out-sequence-add.h5");
{
deleteIfNecessary(filename);
HDF5::File file(filename.c_str());
HDF5::SequenceIO<2, ft> h5writer(file, "data", 2, 3);
addEntry(h5writer, 1, written1);
addEntry(h5writer, 0, written0);
}
std::cout << "wrote: " << std::endl;
std::cout << written0 << std::endl;
std::cout << written1 << std::endl;
Vector read0;
Vector read1;
{
HDF5::File file(filename.c_str(), HDF5::Access::READONLY);
HDF5::SequenceIO<2, ft> h5reader(file, "data", 2, 3);
readEntry(h5reader, 0, read0);
readEntry(h5reader, 1, read1);
}
std::cout << "got: " << std::endl;
std::cout << read0 << std::endl;
std::cout << read1 << std::endl;
for (size_t i = 0; i < 2; ++i)
for (size_t j = 0; j < 3; ++j) {
passed = passed and written0[i][j] == read0[i][j];
passed = passed and written1[i][j] == read1[i][j];
}
return passed;
}
template <class Vector>
bool checkAppendSequence(Vector const &written0, Vector const &written1) {
bool passed = true;
using ft = typename Vector::block_type::field_type;
std::string const filename("out-sequence-append.h5");
{
deleteIfNecessary(filename);
HDF5::File file(filename.c_str());
HDF5::SequenceIO<2, ft> h5writer(file, "data", 2, 3);
appendEntry(h5writer, written0);
appendEntry(h5writer, written1);
}
std::cout << "wrote: " << std::endl;
std::cout << written0 << std::endl;
std::cout << written1 << std::endl;
Vector read0;
Vector read1;
{
HDF5::File file(filename.c_str());
HDF5::SequenceIO<2, ft> h5reader(file, "data", 2, 3);
readEntry(h5reader, 0, read0);
readEntry(h5reader, 1, read1);
}
std::cout << "got: " << std::endl;
std::cout << read0 << std::endl;
std::cout << read1 << std::endl;
for (size_t i = 0; i < 2; ++i)
for (size_t j = 0; j < 3; ++j) {
passed = passed and written0[i][j] == read0[i][j];
passed = passed and written1[i][j] == read1[i][j];
}
return passed;
}
template <class Vector>
bool checkSetSingleton(Vector const &written0) {
using ft = typename Vector::block_type::field_type;
std::string const filename = "out-singleton.h5";
{
deleteIfNecessary(filename);
HDF5::File file(filename);
HDF5::SingletonWriter<2, ft> h5writer(file, "data", 2, 3);
setEntry(h5writer, written0);
}
return true;
}
template <class Action>
void traverseParset(Dune::ParameterTree const &parset, std::string prefix,
Action actor) {
for (auto const &key : parset.getValueKeys())
actor(key, parset[key], prefix);
for (auto const &sub : parset.getSubKeys())
traverseParset(parset.sub(sub), prefix + sub + ".", actor);
}
bool testParameters() {
std::string const filename = "out-attributes.h5";
deleteIfNecessary(filename);
HDF5::File file(filename);
std::cout << "Content of dataset:" << std::endl
<< "===============" << std::endl;
Dune::ParameterTree parset;
parset["a"] = "3";
parset["c.d"] = "something";
parset.report();
std::cout << "===============" << std::endl;
traverseParset(parset, "",
[&](std::string key, std::string value, std::string prefix) {
HDF5::Attribute(file, prefix + key).set(value);
});
bool passed = true;
auto const a = HDF5::Attribute(file, "a");
passed &= a.exists();
passed &= a.get() == "3";
auto const cd = HDF5::Attribute(file, "c.d");
passed &= cd.exists();
passed &= cd.get() == "something";
return passed;
}
int main(int argc, char *argv[]) {
using Vector = Dune::BlockVector<Dune::FieldVector<int, 3>>;
Vector b0 = {
{0, 1, 2}, {3, 4, 5},
};
Vector b1 = {{0, 10, 20}, {30, 40, 50}};
bool passed = true;
passed = passed and testParameters();
passed = passed and checkMissingFile(b0);
passed = passed and checkMissingDataset(b0);
passed = passed and checkAddSequence(b0, b1);
#if H5_VERSION_GE(1, 10, 0)
passed = passed and checkAppendSequence(b0, b1);
#endif
passed = passed and checkSetSingleton(b0);
return passed ? 0 : 1;
}