Hello all, I am observing some strange behavior with Parallel::distributed::triangulation<dim>::save(const * char) function.
I make a triangulation and then move its vertices and then save it to a file. After saving the triangulation and I load it in a different function. However, this loaded triangulation is not consistent with the one saved. I am attaching a code that demonstrate this strange behavior / bug. Although, I am using parallel::distributed::triangulation, I run it with zero elements and on one processor. The output in the "femesh" file is as follows: original cell number 1 Vertex 0: 0 0 Vertex 1: 1 0 Vertex 2: 0 1 Vertex 3: 1 1 before writing to file cell number 1 Vertex 0: 0 0 Vertex 1: 0.75 0 Vertex 2: 0 1 Vertex 3: 0.75 1 after loading from file cell number 1 Vertex 0: 0 0 Vertex 1: 1 0 Vertex 2: 0 1 Vertex 3: 1 1 Thanks for your help. -- The deal.II project is located at http://www.dealii.org/ For mailing list/forum options, see https://groups.google.com/d/forum/dealii?hl=en --- You received this message because you are subscribed to the Google Groups "deal.II User Group" group. To unsubscribe from this group and stop receiving emails from it, send an email to dealii+unsubscr...@googlegroups.com. For more options, visit https://groups.google.com/d/optout.
#include <deal.II/dofs/dof_handler.h> #include <deal.II/fe/fe_q.h> #include <deal.II/fe/fe_values.h> #include <deal.II/grid/grid_generator.h> #include <deal.II/grid/grid_tools.h> #include <deal.II/grid/tria.h> #include <deal.II/grid/tria_accessor.h> #include <deal.II/grid/tria_iterator.h> #include <iostream> #include <fstream> #include <deal.II/grid/tria_boundary_lib.h> #include <deal.II/dofs/dof_handler.h> #include <deal.II/dofs/dof_accessor.h> #include <deal.II/fe/fe_q.h> #include <deal.II/distributed/tria.h> #include <deal.II/distributed/grid_refinement.h> using namespace dealii; template<int dim> void print_to_file (typename DoFHandler<dim>::active_cell_iterator cell, typename DoFHandler<dim>::active_cell_iterator endc, std::string str) { const std::string filename = "femesh"; std::ofstream out(filename.c_str(), std::ios::app); out << str << std::endl; int id = 0; for(;cell != endc; ++cell) if (cell -> is_locally_owned()) { id++; out << "cell number " << id << std::endl; for (unsigned int v = 0; v < GeometryInfo<dim>::vertices_per_cell; ++v) { out << "Vertex\t" << cell->vertex_index(v) << ":\t" << cell->vertex(v) << std::endl; } } out << std::endl; } template<int dim> void load_and_write_mesh() { MPI_Comm mpi_communicator(MPI_COMM_WORLD); parallel::distributed::Triangulation<dim> tria (mpi_communicator, typename Triangulation<dim>::MeshSmoothing (Triangulation<dim>::smoothing_on_refinement | Triangulation<dim>::smoothing_on_coarsening)); GridGenerator::hyper_cube(tria, 0.0, 1.0); const int degree = 1; FE_Q<dim> fe(degree); QGauss<dim> qf_cell(degree+1); std::string file_name = "mymesh"; tria.load(file_name.c_str()); DoFHandler<dim> dof_handler(tria); dof_handler.distribute_dofs(fe); typename DoFHandler<dim>::active_cell_iterator cell = dof_handler.begin_active(), endc = dof_handler.end(); print_to_file<dim>(cell, endc, std::string("after loading from file")); } int main(int argc, char *argv[]) { dealii::Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); const int dim = 2; MPI_Comm mpi_communicator(MPI_COMM_WORLD); parallel::distributed::Triangulation<dim> tria (mpi_communicator, typename Triangulation<dim>::MeshSmoothing (Triangulation<dim>::smoothing_on_refinement | Triangulation<dim>::smoothing_on_coarsening)); const int degree = 1; FE_Q<dim> fe(degree); QGauss<dim> qf_cell(degree+1); DoFHandler<dim> dof_handler(tria); // Single cell of unit length GridGenerator::hyper_cube(tria, 0.0, 1.0); tria.refine_global(0); dof_handler.distribute_dofs(fe); FEValues<dim> fe_values (fe, qf_cell, update_JxW_values); typename DoFHandler<dim>::active_cell_iterator cell = dof_handler.begin_active(), endc = dof_handler.end(); print_to_file<dim>(cell, endc, std::string("original")); std::string file_name = "mymesh"; // Move one of the grid vertices for (typename Triangulation<dim>::active_cell_iterator cell = tria.begin_active(); cell != tria.end(); ++cell) { // Decrease the x-length of the geometry by 25% // By moving only the right-most vertices for (unsigned int v=0; v<GeometryInfo<dim>::vertices_per_cell; ++v) if (std::abs(cell->vertex(v)[0] - 1.0) < 1e-9) cell->vertex(v)[0] *= 0.75; } print_to_file<dim>(cell, endc, std::string("before writing to file")); tria.save(file_name.c_str()); load_and_write_mesh<dim>(); return 0; }
## # CMake script for the step-1 tutorial program: ## # Set the name of the project and target: SET(TARGET mesh_movement) # Declare all source files the target consists of: SET(TARGET_SRC ${TARGET}.cc # You can specify additional files here! ) # Usually, you will not need to modify anything beyond this point... CMAKE_MINIMUM_REQUIRED(VERSION 2.8.8) FIND_PACKAGE(deal.II 8.3 QUIET HINTS ${deal.II_DIR} ${DEAL_II_DIR} ../ ../../ $ENV{DEAL_II_DIR} ) IF(NOT ${deal.II_FOUND}) MESSAGE(FATAL_ERROR "\n" "*** Could not locate deal.II. ***\n\n" "You may want to either pass a flag -DDEAL_II_DIR=/path/to/deal.II to cmake\n" "or set an environment variable \"DEAL_II_DIR\" that contains this path." ) ENDIF() DEAL_II_INITIALIZE_CACHED_VARIABLES() PROJECT(${TARGET}) DEAL_II_INVOKE_AUTOPILOT()