From 026e7efd9ea023e192c6f3d6c312bf4f121cadf6 Mon Sep 17 00:00:00 2001 From: podlesny <podlesny@zedat.fu-berlin.de> Date: Fri, 19 Feb 2021 16:31:29 +0100 Subject: [PATCH] build file and writer objects only once, report weightedNormalStress possible now --- dune/tectonic/io/hdf5-bodywriter.hh | 4 ++ dune/tectonic/io/hdf5-writer.hh | 37 +++++++++++----- .../io/hdf5/frictionalboundary-writer.hh | 19 ++++++++- dune/tectonic/io/io-handler.hh | 42 ++++++++++--------- 4 files changed, 71 insertions(+), 31 deletions(-) diff --git a/dune/tectonic/io/hdf5-bodywriter.hh b/dune/tectonic/io/hdf5-bodywriter.hh index cf47eaad..8201b784 100644 --- a/dune/tectonic/io/hdf5-bodywriter.hh +++ b/dune/tectonic/io/hdf5-bodywriter.hh @@ -52,6 +52,10 @@ class HDF5BodyWriter { #endif*/ } + void reportWeightedNormalStress(const ProgramState& programState) { + frictionBoundaryWriter_->writeWeightedNormalStress(programState.timeStep, programState.weightedNormalStress[id_], programState.weights[id_]); + } + auto id() const { return id_; } diff --git a/dune/tectonic/io/hdf5-writer.hh b/dune/tectonic/io/hdf5-writer.hh index 00cc66ba..948645e7 100644 --- a/dune/tectonic/io/hdf5-writer.hh +++ b/dune/tectonic/io/hdf5-writer.hh @@ -16,6 +16,7 @@ class HDF5Writer { using VertexCoordinates = std::vector<typename HDF5BodyWriter::VertexCoordinates>; using VertexBases = std::vector<const VertexBasis*>; using FrictionPatches = std::vector<const typename HDF5BodyWriter::Patch*>; + using ScalarVector = typename ProgramState::ScalarVector; //using WeakPatches = std::vector<typename HDF5BodyWriter::WeakPatches>; //friend class HDF5NetworkWriter<ProgramState, VertexBasis, GridView>; @@ -26,11 +27,11 @@ class HDF5Writer { const FrictionPatches& frictionPatches) //const WeakPatches& weakPatches) : file_(file), + frictionPatches_(frictionPatches), iterationWriter_(file_), - timeWriter_(file_), - frictionPatches_(frictionPatches) { + timeWriter_(file_) { - for (size_t i=0; i<frictionPatches_.size(); i++) { + for (size_t i=0; i<vertexCoordinates.size(); i++) { if (frictionPatches_[i]->numVertices() > 0) bodyWriters_.push_back(std::make_unique<HDF5BodyWriter>(file_, i, vertexCoordinates[i], *vertexBases[i], *frictionPatches_[i])); //, weakPatches[i])); } @@ -41,7 +42,7 @@ class HDF5Writer { timeWriter_.write(programState); - //friction.updateAlpha(programState.alpha); + friction.updateAlpha(programState.alpha); // extract relative velocities using Vector = typename ProgramState::Vector; @@ -54,17 +55,23 @@ class HDF5Writer { using ScalarVector = typename ProgramState::ScalarVector; const auto frictionCoeff = friction.coefficientOfFriction(mortarV); - double norm = 0; + /*double norm = 0; const auto& bodyNodes = *frictionPatches_[0]->getVertices(); for (size_t i=bodyNodes.size(); i<frictionCoeff.size(); i++) { norm += frictionCoeff[i].two_norm(); - } - std::cout << std::setprecision(10) << "friction coefficients norm: " << norm << std::endl; + } */ + //std::cout << std::setprecision(10) << "friction coefficients norm: " << norm << std::endl; std::vector<ScalarVector> splitCoeff; split(frictionCoeff, splitCoeff); - //print(v_rel, "v_rel: "); + /*if (std::isnan(norm) or std::isinf(norm)) { + print(programState.alpha, "alpha"); + print(frictionCoeff, "frictionCoeff"); + print(splitCoeff, "splitCoeff"); + print(v_rel, "v_rel: "); + DUNE_THROW(Dune::Exception, "invalid state"); + }*/ for (size_t i=0; i<bodyWriters_.size(); i++) { auto bodyID = bodyWriters_[i]->id(); @@ -76,13 +83,21 @@ class HDF5Writer { iterationWriter_.write(programState.timeStep, iterationCount); } + + void reportWeightedNormalStress(const ProgramState& programState) { + for (size_t i=0; i<bodyWriters_.size(); i++) { + auto bodyID = bodyWriters_[i]->id(); + bodyWriters_[i]->reportWeightedNormalStress(programState); + } + } + template <class VectorType> void split(const VectorType& v, std::vector<VectorType>& splitV) const { - - const size_t nBodies = frictionPatches_.size(); + const auto nBodies = frictionPatches_.size(); size_t globalIdx = 0; splitV.resize(nBodies); + for (size_t bodyIdx=0; bodyIdx<nBodies; bodyIdx++) { const auto& bodyNodes = *frictionPatches_[bodyIdx]->getVertices(); @@ -100,11 +115,11 @@ class HDF5Writer { private: HDF5::File& file_; + const FrictionPatches& frictionPatches_; IterationWriter iterationWriter_; TimeWriter<ProgramState> timeWriter_; - const FrictionPatches& frictionPatches_; std::vector<std::unique_ptr<HDF5BodyWriter>> bodyWriters_; }; #endif diff --git a/dune/tectonic/io/hdf5/frictionalboundary-writer.hh b/dune/tectonic/io/hdf5/frictionalboundary-writer.hh index bc8ca848..b833dac5 100644 --- a/dune/tectonic/io/hdf5/frictionalboundary-writer.hh +++ b/dune/tectonic/io/hdf5/frictionalboundary-writer.hh @@ -25,10 +25,15 @@ template <class GridView> class FrictionalBoundaryWriter { frictionalBoundaryStateWriter_(group_, "state", frictionalBoundary.numVertices()), frictionalBoundaryCoefficientWriter_(group_, "coefficient", - frictionalBoundary.numVertices()) { + frictionalBoundary.numVertices()), + frictionalBoundaryWeightsWriter_(group_, "weights", + frictionalBoundary.numVertices()), + frictionalBoundaryWeightedNormalStressWriter_(group_, "weightedNormalStress", + frictionalBoundary.numVertices()){ auto const frictionalBoundaryCoordinates = restrictToSurface(vertexCoordinates, frictionalBoundary); + HDF5::SingletonWriter<2> frictionalBoundaryCoordinateWriter( group_, "coordinates", frictionalBoundaryCoordinates.size(), Vector::block_type::dimension); @@ -51,6 +56,16 @@ template <class GridView> class FrictionalBoundaryWriter { addEntry(frictionalBoundaryCoefficientWriter_, timeStep, frictionalBoundaryCoefficient); } + template <class ScalarVector> + void writeWeightedNormalStress(const size_t timeStep, const ScalarVector& weightedNormalStress, const ScalarVector& weights) { + + auto const frictionalBoundaryWeights = restrictToSurface(weights, frictionalBoundary_); + addEntry(frictionalBoundaryWeightsWriter_, timeStep, frictionalBoundaryWeights); + + auto const frictionalBoundaryWeightedNormalStress = restrictToSurface(weightedNormalStress, frictionalBoundary_); + addEntry(frictionalBoundaryWeightedNormalStressWriter_, timeStep, frictionalBoundaryWeightedNormalStress); + } + private: HDF5::Group& group_; @@ -60,6 +75,8 @@ template <class GridView> class FrictionalBoundaryWriter { HDF5::SequenceIO<2> frictionalBoundaryVelocityWriter_; HDF5::SequenceIO<1> frictionalBoundaryStateWriter_; HDF5::SequenceIO<1> frictionalBoundaryCoefficientWriter_; + HDF5::SequenceIO<1> frictionalBoundaryWeightsWriter_; + HDF5::SequenceIO<1> frictionalBoundaryWeightedNormalStressWriter_; }; #endif diff --git a/dune/tectonic/io/io-handler.hh b/dune/tectonic/io/io-handler.hh index 0e7802e4..be4e491a 100644 --- a/dune/tectonic/io/io-handler.hh +++ b/dune/tectonic/io/io-handler.hh @@ -2,6 +2,7 @@ #define SRC_IO_HANDLER_HH #include <filesystem> +#include <memory> #include <dune/common/parametertree.hh> @@ -24,7 +25,7 @@ #include "hdf5-writer.hh" #include "vtk.hh" -template <class Assembler, class ContactNetwork> +template <class Assembler, class ContactNetwork, class ProgramState> class IOHandler { private: using Vector = typename Assembler::Vector; @@ -32,6 +33,7 @@ class IOHandler { using MyVertexBasis = typename Assembler::VertexBasis; using MyCellBasis = typename Assembler::CellBasis; + using HDF5Writer = HDF5Writer<ProgramState, MyVertexBasis, typename Assembler::GV>; const size_t bodyCount_; std::vector<size_t> nVertices_; @@ -39,8 +41,12 @@ class IOHandler { std::vector<Vector> vertexCoordinates_; std::vector<const MyVertexBasis* > vertexBases_; std::vector<const MyCellBasis* > cellBases_; + std::vector<const typename ContactNetwork::BoundaryPatch*> frictionPatches_; std::vector<std::vector<const ConvexPolyhedron<LocalVector>*>> weakPatches_; + std::unique_ptr<HDF5::File> dataFile_; + std::unique_ptr<HDF5Writer> dataWriter_; + bool writeVTK_; bool writeData_; @@ -88,9 +94,15 @@ class IOHandler { if (!std::filesystem::is_directory("iterates/")) std::filesystem::create_directory("iterates"); } + + + contactNetwork.frictionPatches(frictionPatches_); + dataFile_ = std::make_unique<HDF5::File>("output.h5"); + dataWriter_ = std::make_unique<HDF5Writer>(*dataFile_, vertexCoordinates_, vertexBases_, + frictionPatches_); } - template <class ProgramState, class GlobalFriction> + template <class GlobalFriction> void write(const ProgramState& programState, const ContactNetwork& contactNetwork, const GlobalFriction& friction, const IterationRegister& iterationCount, bool initial = false) { if (writeData_) { writeData(programState, contactNetwork, friction, iterationCount, initial); @@ -112,7 +124,6 @@ class IOHandler { } } - template <class ProgramState> bool read(ProgramState& programState) { using ADS = RestartIO<ProgramState>; @@ -128,7 +139,6 @@ class IOHandler { } private: - template <class ProgramState> void writeVTK(const ProgramState& programState, const ContactNetwork& contactNetwork) { using ScalarVector = typename Assembler::ScalarVector; std::vector<ScalarVector> stress(bodyCount_); @@ -146,7 +156,6 @@ class IOHandler { programState.alpha, stress); } - template <class ProgramState> void writeRestarts(const ProgramState& programState) { if (programState.timeStep % restartSpacing_ == 0) { auto restartFile = std::make_unique<HDF5::File>("restarts.h5", HDF5::Access::READWRITE); @@ -157,21 +166,16 @@ class IOHandler { } } - template <class ProgramState, class GlobalFriction> + template <class GlobalFriction> void writeData(const ProgramState& programState, const ContactNetwork& contactNetwork, const GlobalFriction& friction, const IterationRegister& iterationCount, bool initial) { - std::vector<const typename ContactNetwork::BoundaryPatch*> frictionPatches; - contactNetwork.frictionPatches(frictionPatches); - - auto dataFile = std::make_unique<HDF5::File>("output.h5"); - auto dataWriter = std::make_unique< - HDF5Writer<ProgramState, MyVertexBasis, typename Assembler::GV>>( - *dataFile, vertexCoordinates_, vertexBases_, - frictionPatches); - - dataWriter->reportSolution(programState, contactNetwork, friction); - if (!initial) - dataWriter->reportIterations(programState, iterationCount); - dataFile->flush(); + dataWriter_->reportSolution(programState, contactNetwork, friction); + if (!initial) { + dataWriter_->reportIterations(programState, iterationCount); + } else { + dataWriter_->reportWeightedNormalStress(programState); + } + + dataFile_->flush(); } }; -- GitLab