From e146075b126833e2adb2ef13cd30f8076a89827e Mon Sep 17 00:00:00 2001 From: Lei Qiao Date: Thu, 24 Sep 2015 18:20:26 -0500 Subject: [PATCH] add test distributed_grids/3d_refinement_12 --- tests/distributed_grids/3d_refinement_12.cc | 335 +++++++++++++++ .../distributed_grids/3d_refinement_12.output | 386 ++++++++++++++++++ 2 files changed, 721 insertions(+) create mode 100644 tests/distributed_grids/3d_refinement_12.cc create mode 100644 tests/distributed_grids/3d_refinement_12.output diff --git a/tests/distributed_grids/3d_refinement_12.cc b/tests/distributed_grids/3d_refinement_12.cc new file mode 100644 index 0000000000..5e0d8a9477 --- /dev/null +++ b/tests/distributed_grids/3d_refinement_12.cc @@ -0,0 +1,335 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2015 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + + +#include "../tests.h" +#include +#include +#include +#include +#include +#include + +using namespace dealii; + +// Purpose of this test: +// Demonstrate that mesh smoothing option limit_level_difference_at_vertices +// is always effective across different initial cells (different trees in p4est context) +// not matter it is set or not. +// +// Design of this test: +// The initial mesh is set to 8 cells in unite cube. Then the center most cell +// in each level is refined consecutively for 4 times. Finally, total active +// cell number in each refine level and cell center coordinates of the final +// active cells are compared for cases with and without flag mesh smoothing +// limit_level_difference_at_vertices. +// +// Expected result: +// In principle, the result of the two refinement runs should be different. +// However, under the current implementation, their results are the same. +// So, if some day this test failed, it might be a good news. + +// Now the output of mesh and refine indicators in vtu format is disabled to +// prevent unnecessary outputs during regression test. It you want to see the +// mesh, uncomment the following macro definition. +// #define __WRITE_VTU__ + +template class Location; + +template +class TriaTest +{ +private: + typedef typename parallel::distributed::Triangulation TypeTria; +public: + TriaTest(const typename dealii::Triangulation::MeshSmoothing smoothing_option = dealii::Triangulation::none); + void run(std::vector &n_cell, + std::set > &position_list); +private: + void write_vtu(const unsigned int counter) const; + + const MPI_Comm mpi_communicator; + TypeTria triangulation; + const unsigned int myid; + const bool I_am_host; + const std::string case_name; +}; + +// Equip class Point with operator < in order to sort the cell center coordinates. +template +class Location: public Point +{ +public: + Location(const Point &in) + : + Point(in) + {} + + inline + bool operator<(const Location &op) const + { + for (unsigned int d=0; d n_cell_smooth; + std::vector n_cell_no_smooth; + + std::set > final_cell_center_loactions_smooth; + std::set > final_cell_center_loactions_no_smooth; + if (I_am_host) + { + deallog << "Flag limit_level_difference_at_vertices set:" << std::endl; + } + { + TriaTest tria_test(dealii::Triangulation::limit_level_difference_at_vertices); + tria_test.run(n_cell_smooth, final_cell_center_loactions_smooth); + } + if (I_am_host) + { + deallog << "Flag limit_level_difference_at_vertices unset:" << std::endl; + } + { + TriaTest tria_test(dealii::Triangulation::none); + tria_test.run(n_cell_no_smooth, final_cell_center_loactions_no_smooth); + } + if (I_am_host) + { + deallog << "Compare result:" << std::endl; + } + { + bool n_cells_are_same = (n_cell_smooth.size() == n_cell_no_smooth.size()); + + for (unsigned int i = 0; n_cells_are_same && (i repetitions; + Point p1; + Point p2; + + for (unsigned int d=0; d +void TriaTest::run(std::vector &n_cell, + std::set > &position_list) +{ + n_cell.clear(); + position_list.clear(); + + unsigned int counter = 0; + if (I_am_host) + { + deallog << "n_loop n_cell" << std::endl; + deallog << counter << " " + << triangulation.n_global_active_cells() << std::endl; + } + ++counter; + + for (; counter < 5; ++counter) + { + { + Point p; + for (unsigned int d=0; dis_locally_owned() && ((cell->center()).distance(p) < 1e-4)) + { + cell->set_refine_flag(); + } + } + + triangulation.prepare_coarsening_and_refinement(); + write_vtu(counter); + triangulation.execute_coarsening_and_refinement(); + + if (I_am_host) + { + deallog << counter << " " + << triangulation.n_global_active_cells() << std::endl; + } + n_cell.push_back(triangulation.n_global_active_cells()); + } + + // Output cell center coordinates in a sorted manner to prevent change in + // cell ordering causing failure of this test. + // This is just a local operation, you have to re-implement this part if you + // want to run this test in distributed parallel. + { + deallog << " position of cell centers:" << std::endl; + + for (typename TypeTria::active_cell_iterator cell = triangulation.begin_active(); + cell!= triangulation.end(); + ++cell) + if (cell->is_locally_owned()) + { + Location loc(cell->center()); + position_list.insert(loc); + } + + for (typename std::set >::const_iterator it = position_list.begin(); + it != position_list.end(); + ++it) + { + deallog << *it << std::endl; + } + } + + return; +} + +template +void TriaTest::write_vtu(const unsigned int counter) const +{ +#ifdef __WRITE_VTU__ + // save refine flag + Vector refine_mark; + { + refine_mark.reinit (triangulation.n_active_cells()); + + typename TypeTria::active_cell_iterator + cell = triangulation.begin_active(); + const typename TypeTria::active_cell_iterator + endc = triangulation.end(); + for (; cell != endc; ++cell) + if (cell->is_locally_owned()) + { + refine_mark[cell->active_cell_index()] = + cell->refine_flag_set(); + } + } + + DataOut data_out; + data_out.attach_triangulation (triangulation); + { + const std::string data_name ("refine_flag"); + data_out.add_data_vector (refine_mark, + data_name, + DataOut::type_cell_data); + } + + data_out.build_patches(); + + const std::string output_tag = case_name + + Utilities::int_to_string (counter, 4); + const std::string slot_itag = ".slot-" + Utilities::int_to_string (myid, 4); + + std::ofstream output ((output_tag + slot_itag + ".vtu").c_str()); + data_out.write_vtu (output); + + if (I_am_host) + { + std::vector filenames; + for (unsigned int i=0; + i