From 8ab7a52a2114ded2c8d4166c79cc732817968ad1 Mon Sep 17 00:00:00 2001 From: Elias Pipping <elias.pipping@fu-berlin.de> Date: Thu, 24 Jul 2014 14:10:56 +0200 Subject: [PATCH] [Cleanup] Print grid info --- src/sand-wedge.cc | 30 ++++++++++++++++++++++++++++++ 1 file changed, 30 insertions(+) diff --git a/src/sand-wedge.cc b/src/sand-wedge.cc index 0dfdc6f0..cf831898 100644 --- a/src/sand-wedge.cc +++ b/src/sand-wedge.cc @@ -44,6 +44,7 @@ #include <dune/fufem/sharedpointermap.hh> #include <dune/solvers/norms/energynorm.hh> #include <dune/solvers/norms/sumnorm.hh> +#include <dune/solvers/norms/twonorm.hh> #include <dune/solvers/solvers/loopsolver.hh> #include <dune/solvers/solvers/solver.hh> #include <dune/tnnmg/nonlinearities/zerononlinearity.hh> @@ -81,6 +82,21 @@ void initPython() { Python::run("sys.path.append('" datadir "')"); } +template <class Geometry> double diameter(Geometry const &geometry) { + auto const numCorners = geometry.corners(); + std::vector<typename Geometry::GlobalCoordinate> corners(numCorners); + for (int i = 0; i < numCorners; ++i) + corners[i] = geometry.corner(i); + + TwoNorm<typename Geometry::GlobalCoordinate> twoNorm; + + double diameter = 0.0; + for (int i = 0; i < numCorners; ++i) + for (int j = 0; j < i; ++j) + diameter = std::max(diameter, twoNorm.diff(corners[i], corners[j])); + return diameter; +} + int main(int argc, char *argv[]) { try { Dune::ParameterTree parset; @@ -97,10 +113,24 @@ int main(int argc, char *argv[]) { auto const refinements = parset.get<size_t>("grid.refinements"); grid->globalRefine(refinements); + double minDiameter = std::numeric_limits<double>::infinity(); + double maxDiameter = 0.0; + for (auto it = grid->template leafbegin<0>(); + it != grid->template leafend<0>(); ++it) { + auto const geometry = it->geometry(); + auto const diam = diameter(geometry); + minDiameter = std::min(minDiameter, diam); + maxDiameter = std::max(maxDiameter, diam); + } + std::cout << "min diameter: " << minDiameter << std::endl; + std::cout << "max diameter: " << maxDiameter << std::endl; + using GridView = Grid::LeafGridView; GridView const leafView = grid->leafGridView(); size_t const leafVertexCount = leafView.size(dims); + std::cout << "Number of DOFs: " << leafVertexCount << std::endl; + auto myFaces = gridConstructor.constructFaces(leafView); // Neumann boundary -- GitLab