From 646610a6a34101fd972a598fd3de56cf743a3b42 Mon Sep 17 00:00:00 2001 From: bangerth Date: Wed, 14 Sep 2011 23:56:36 +0000 Subject: [PATCH] Make GridTools, VectorTools, MGTools namespaces, rather than a class with all static member functions. Now only MatrixTools remains. git-svn-id: https://svn.dealii.org/trunk@24324 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/doc/news/changes.h | 14 +- deal.II/include/deal.II/grid/grid_tools.h | 2061 +++-- deal.II/include/deal.II/multigrid/mg_tools.h | 488 +- deal.II/include/deal.II/numerics/vectors.h | 3319 ++++--- .../deal.II/numerics/vectors.templates.h | 8080 ++++++++--------- deal.II/source/grid/grid_tools.cc | 2370 ++--- deal.II/source/grid/grid_tools.inst.in | 87 +- deal.II/source/multigrid/mg_tools.cc | 3475 ++++--- deal.II/source/multigrid/mg_tools.inst.in | 161 +- .../source/multigrid/mg_transfer_component.cc | 2 +- deal.II/source/numerics/vectors.cc | 53 +- deal.II/source/numerics/vectors.inst.in | 204 +- 12 files changed, 10117 insertions(+), 10197 deletions(-) diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h index a5248326e6..19a097f0f9 100644 --- a/deal.II/doc/news/changes.h +++ b/deal.II/doc/news/changes.h @@ -21,6 +21,13 @@ inconvenience this causes.

    +
  1. Changed: GridTools, DoFTools, MGTools and VectorTools are now namespaces. They have long +been classes that had only public, static member functions, making +the end result semantically exactly equivalent to a namespace, which is +also how it was used. This is now also reflected in the actual code. +
    +(Wolfgang Bangerth, 2011/04/27, 2011/09/14) +
  2. Changed: The PETScWrapper::VectorBase and PETScWrapper::MatrixBase classes tried to keep track of whether the last operation done on a vector was to add to an element or to @@ -508,13 +515,6 @@ TrilinosWrappers::SparsityPattern::exists(). This is now fixed.
    (Habib Talavatifard, Wolfgang Bangerth 2011/05/09, 2011/05/27) -
  3. Changed: DoFTools is now a namespace. It has long been a class that -had only public, static member functions, making the end result semantically -exactly equivalent to a namespace, which is also how it was used. This is -now also reflected in the actual code. -
    -(Wolfgang Bangerth, 2011/04/27) -
  4. New: The version of DoFTools::make_flux_sparsity_pattern that takes the coupling masks is now also available for hp::DoFHandler objects.
    diff --git a/deal.II/include/deal.II/grid/grid_tools.h b/deal.II/include/deal.II/grid/grid_tools.h index 917736ad39..f17b2c0f3f 100644 --- a/deal.II/include/deal.II/grid/grid_tools.h +++ b/deal.II/include/deal.II/grid/grid_tools.h @@ -37,1082 +37,1051 @@ class SparsityPattern; /** - * This class is a collection of algorithms working on triangulations, + * This namespace is a collection of algorithms working on triangulations, * such as shifting or rotating triangulations, but also finding a * cell that contains a given point. See the descriptions of the * individual functions for more information. * * @ingroup grid */ -class GridTools +namespace GridTools { - public: - /** - * Return the diameter of a - * triangulation. The diameter is - * computed using only the - * vertices, i.e. if the diameter - * should be larger than the - * maximal distance between - * boundary vertices due to a - * higher order mapping, then - * this function will not catch - * this. - */ - template - static - double diameter (const Triangulation &tria); - - /** - * Compute the volume (i.e. the dim-dimensional measure) of the - * triangulation. We compute the measure using the integral - * $\int 1 \; dx$. The integral approximated is approximated - * via quadrature for which we need the mapping argument. - * - * This function also works for objects of type - * parallel::distributed::Triangulation, in which case the - * function is a collective operation. - */ - template - static - double volume (const Triangulation &tria, - const Mapping &mapping = (StaticMappingQ1::mapping)); - - /** - * Given a list of vertices (typically - * obtained using - * Triangulation::get_vertices) as the - * first, and a list of vertex indices - * that characterize a single cell as the - * second argument, return the measure - * (area, volume) of this cell. If this - * is a real cell, then you can get the - * same result using - * cell-@>measure(), but - * this function also works for cells - * that do not exist except that you make - * it up by naming its vertices from the - * list. - */ - template - static - double cell_measure (const std::vector > &all_vertices, - const unsigned int (&vertex_indices)[GeometryInfo::vertices_per_cell]); - - /** - * Remove vertices that are not - * referenced by any of the - * cells. This function is called - * by all GridIn::read_* - * functions to eliminate - * vertices that are listed in - * the input files but are not - * used by the cells in the input - * file. While these vertices - * should not be in the input - * from the beginning, they - * sometimes are, most often when - * some cells have been removed - * by hand without wanting to - * update the vertex lists, as - * they might be lengthy. - * - * This function is called by all - * GridIn::read_* - * functions as the triangulation - * class requires them to be - * called with used vertices - * only. This is so, since the - * vertices are copied verbatim - * by that class, so we have to - * eliminate unused vertices - * beforehand. - * - * Not implemented for the - * codimension one case. - */ - template - static - void delete_unused_vertices (std::vector > &vertices, - std::vector > &cells, - SubCellData &subcelldata); - - /** - * Remove vertices that are duplicated, - * due to the input of a structured grid, - * for example. If these vertices are not - * removed, the faces bounded by these - * vertices become part of the boundary, - * even if they are in the interior of - * the mesh. - * - * This function is called by some - * GridIn::read_* functions. Only - * the vertices with indices in @p - * considered_vertices are tested for - * equality. This speeds up the algorithm, - * which is quadratic and thus quite slow - * to begin with. However, if you wish to - * consider all vertices, simply pass an - * empty vector. - * - * Two vertices are considered equal if - * their difference in each coordinate - * direction is less than @p tol. - */ - template - static - void delete_duplicated_vertices (std::vector > &all_vertices, - std::vector > &cells, - SubCellData &subcelldata, - std::vector &considered_vertices, - const double tol=1e-12); - - /** - * Transform the vertices of the given - * triangulation by applying the - * predicate to all its vertices. Since - * the internal consistency of a - * triangulation can only be guaranteed - * if the transformation is applied to - * the vertices of only one level of a - * hierarchically refined cells, this - * function may only be used if all cells - * of the triangulation are on the same - * refinement level. - * - * The predicate given as - * argument is used to transform - * each vertex. Its respective - * type has to offer a - * function-like syntax, i.e. the - * predicate is either an object - * of a type that has an - * operator(), or it is a - * pointer to the function. In - * either case, argument and - * return value have to be of - * type Point. - * - * This function is used in the - * "Possibilities for extensions" section - * of step-38. - */ - template - static - void transform (const Predicate &predicate, - Triangulation &triangulation); - - /** - * Shift each vertex of the - * triangulation by the given - * shift vector. This function - * uses the transform() - * function above, so the - * requirements on the - * triangulation stated there - * hold for this function as - * well. - */ - template - static - void shift (const Point &shift_vector, - Triangulation &triangulation); - - - /** - * Rotate all vertices of the - * given two-dimensional - * triangulation in - * counter-clockwise sense around - * the origin of the coordinate - * system by the given angle - * (given in radians, rather than - * degrees). This function uses - * the transform() function - * above, so the requirements on - * the triangulation stated there - * hold for this function as - * well. - */ - static - void rotate (const double angle, - Triangulation<2> &triangulation); - - /** - * Scale the entire triangulation - * by the given factor. To - * preserve the orientation of - * the triangulation, the factor - * must be positive. - * - * This function uses the - * transform() function - * above, so the requirements on - * the triangulation stated there - * hold for this function as - * well. - */ - template - static - void scale (const double scaling_factor, - Triangulation &triangulation); - - /** - * Find and return the number of - * the used vertex in a given - * Container that is located closest - * to a given point @p p. The - * type of the first parameter - * may be either Triangulation, - * DoFHandler, hp::DoFHandler, or - * MGDoFHandler. - * - * @author Ralf B. Schulz, 2006 - */ - template class Container, int spacedim> - static - unsigned int - find_closest_vertex (const Container &container, - const Point &p); - - /** - * Find and return a vector of - * iterators to active cells that - * surround a given vertex @p vertex. - * The type of the first parameter - * may be either Triangulation, - * DoFHandler, hp::DoFHandler, or - * MGDoFHandler. - * - * For locally refined grids, the - * vertex itself might not be a vertex - * of all adjacent cells, but will - * always be located on a face or an - * edge of the adjacent cells returned. - * - * @author Ralf B. Schulz, - * Wolfgang Bangerth, 2006 - */ - template class Container, int spacedim> - static - std::vector::active_cell_iterator> - find_cells_adjacent_to_vertex (const Container &container, - const unsigned int vertex); - - - /** - * Find and return an iterator to - * the active cell that surrounds - * a given point @p ref. The - * type of the first parameter - * may be either - * Triangulation, - * DoFHandler, or - * MGDoFHandler, i.e. we - * can find the cell around a - * point for iterators into each - * of these classes. - * - * This is solely a wrapper function - * for the @p interpolate function - * given below, - * providing backward compatibility. - * A Q1 mapping is used for the - * boundary, and the iterator to - * the cell in which the point - * resides is returned. - * - * It is recommended to use the - * other version of this function, - * as it simultaneously delivers the - * local coordinate of the given point - * without additional computational cost. - */ - template class Container, int spacedim> - static - typename Container::active_cell_iterator - find_active_cell_around_point (const Container &container, - const Point &p); - - /** - * Find and return an iterator to - * the active cell that surrounds - * a given point @p p. The - * type of the first parameter - * may be either - * Triangulation, - * DoFHandler, hp::DoFHandler, or - * MGDoFHandler, i.e., we - * can find the cell around a - * point for iterators into each - * of these classes. - * - * The algorithm used in this - * function proceeds by first - * looking for vertex located - * closest to the given point, see - * find_closest_vertex(). Secondly, - * all adjacent cells to this point - * are found in the mesh, see - * find_cells_adjacent_to_vertex(). - * Lastly, for each of these cells, - * it is tested whether the point is - * inside. This check is performed - * using arbitrary boundary mappings. - * Still, it is possible that due - * to roundoff errors, the point - * cannot be located exactly inside - * the unit cell. In this case, - * even points at a very small - * distance outside the unit cell - * are allowed. - * - * If a point lies on the - * boundary of two or more cells, - * then the algorithm tries to identify - * the cell that is of highest - * refinement level. - * - * The function returns an - * iterator to the cell, as well - * as the local position of the - * point inside the unit - * cell. This local position - * might be located slightly - * outside an actual unit cell, - * due to numerical roundoff. - * Therefore, the point returned - * by this function should - * be projected onto the unit cell, - * using GeometryInfo::project_to_unit_cell. - * This is not automatically performed - * by the algorithm. - */ - template class Container, int spacedim> - static - std::pair::active_cell_iterator, Point > - find_active_cell_around_point (const Mapping &mapping, - const Container &container, - const Point &p); - - /** - * A version of the previous function - * where we use that mapping on a given - * cell that corresponds to the active - * finite element index of that - * cell. This is obviously only useful - * for hp problems, since the active - * finite element index for all other DoF - * handlers is always zero. - */ - template - static - std::pair::active_cell_iterator, Point > - find_active_cell_around_point (const hp::MappingCollection &mapping, - const hp::DoFHandler &container, - const Point &p); - - /** - * Return a list of all descendents of - * the given cell that are active. For - * example, if the current cell is once - * refined but none of its children are - * any further refined, then the returned - * list will contain all its children. - * - * If the current cell is already active, - * then the returned list is empty - * (because the cell has no children that - * may be active). - * - * Since in C++ the type of the Container - * template argument (which can be - * Triangulation, DoFHandler, - * MGDoFHandler, or hp::DoFHandler) can - * not be deduced from a function call, - * you will have to specify it after the - * function name, as for example in - * GridTools::get_active_child_cells@ - * @> (cell). - */ - template - static - std::vector - get_active_child_cells (const typename Container::cell_iterator &cell); - - /** - * Extract the active cells around a given - * cell @p cell and return them in the - * vector @p active_neighbors. - */ - template - static void - get_active_neighbors (const typename Container::active_cell_iterator &cell, - std::vector &active_neighbors); - - /** - * Produce a sparsity pattern in which - * nonzero entries indicate that two - * cells are connected via a common - * face. The diagonal entries of the - * sparsity pattern are also set. - * - * The rows and columns refer to the - * cells as they are traversed in their - * natural order using cell iterators. - */ - template - static void - get_face_connectivity_of_cells (const Triangulation &triangulation, - SparsityPattern &connectivity); - - /** - * Use the METIS partitioner to generate - * a partitioning of the active cells - * making up the entire domain. After - * calling this function, the subdomain - * ids of all active cells will have - * values between zero and - * @p n_partitions-1. You can access the - * subdomain id of a cell by using - * cell-@>subdomain_id(). - * - * This function will generate an error - * if METIS is not installed unless - * @p n_partitions is one. I.e., you can - * write a program so that it runs in the - * single-processor single-partition case - * without METIS installed, and only - * requires METIS when multiple - * partitions are required. - */ - template - static - void - partition_triangulation (const unsigned int n_partitions, - Triangulation &triangulation); - - /** - * This function does the same as the - * previous one, i.e. it partitions a - * triangulation using METIS into a - * number of subdomains identified by the - * cell-@>subdomain_id() - * flag. - * - * The difference to the previous - * function is the second argument, a - * sparsity pattern that represents the - * connectivity pattern between cells. - * - * While the function above builds it - * directly from the triangulation by - * considering which cells neighbor each - * other, this function can take a more - * refined connectivity graph. The - * sparsity pattern needs to be of size - * $N\times N$, where $N$ is the number - * of active cells in the - * triangulation. If the sparsity pattern - * contains an entry at position $(i,j)$, - * then this means that cells $i$ and $j$ - * (in the order in which they are - * traversed by active cell iterators) - * are to be considered connected; METIS - * will then try to partition the domain - * in such a way that (i) the subdomains - * are of roughly equal size, and (ii) a - * minimal number of connections are - * broken. - * - * This function is mainly useful in - * cases where connections between cells - * exist that are not present in the - * triangulation alone (otherwise the - * previous function would be the simpler - * one to use). Such connections may - * include that certain parts of the - * boundary of a domain are coupled - * through symmetric boundary conditions - * or integrals (e.g. friction contact - * between the two sides of a crack in - * the domain), or if a numerical scheme - * is used that not only connects - * immediate neighbors but a larger - * neighborhood of cells (e.g. when - * solving integral equations). - * - * In addition, this function may be - * useful in cases where the default - * sparsity pattern is not entirely - * sufficient. This can happen because - * the default is to just consider face - * neighbors, not neighboring cells that - * are connected by edges or - * vertices. While the latter couple when - * using continuous finite elements, they - * are typically still closely connected - * in the neighborship graph, and METIS - * will not usually cut important - * connections in this case. However, if - * there are vertices in the mesh where - * many cells (many more than the common - * 4 or 6 in 2d and 3d, respectively) - * come together, then there will be a - * significant number of cells that are - * connected across a vertex, but several - * degrees removed in the connectivity - * graph built only using face - * neighbors. In a case like this, METIS - * may sometimes make bad decisions and - * you may want to build your own - * connectivity graph. - */ - template - static - void - partition_triangulation (const unsigned int n_partitions, - const SparsityPattern &cell_connection_graph, - Triangulation &triangulation); - - /** - * For each active cell, return in the - * output array to which subdomain (as - * given by the cell->subdomain_id() - * function) it belongs. The output array - * is supposed to have the right size - * already when calling this function. - * - * This function returns the association - * of each cell with one subdomain. If - * you are looking for the association of - * each @em DoF with a subdomain, use the - * DoFTools::get_subdomain_association - * function. - */ - template - static void - get_subdomain_association (const Triangulation &triangulation, - std::vector &subdomain); - - /** - * Count how many cells are uniquely - * associated with the given @p subdomain - * index. - * - * This function may return zero - * if there are no cells with the - * given @p subdomain index. This - * can happen, for example, if - * you try to partition a coarse - * mesh into more partitions (one - * for each processor) than there - * are cells in the mesh. - * - * This function returns the number of - * cells associated with one - * subdomain. If you are looking for the - * association of @em DoFs with this - * subdomain, use the - * DoFTools::count_dofs_with_subdomain_association - * function. - */ - template - static unsigned int - count_cells_with_subdomain_association (const Triangulation &triangulation, - const types::subdomain_id_t subdomain); - - /** - * Given two mesh containers - * (i.e. objects of type - * Triangulation, DoFHandler, - * hp::DoFHandler, or - * MGDoFHandler) that are based - * on the same coarse mesh, this - * function figures out a set of - * cells that are matched between - * the two meshes and where at - * most one of the meshes is more - * refined on this cell. In other - * words, it finds the smallest - * cells that are common to both - * meshes, and that together - * completely cover the domain. - * - * This function is useful, for - * example, in time-dependent or - * nonlinear application, where - * one has to integrate a - * solution defined on one mesh - * (e.g., the one from the - * previous time step or - * nonlinear iteration) against - * the shape functions of another - * mesh (the next time step, the - * next nonlinear iteration). If, - * for example, the new mesh is - * finer, then one has to obtain - * the solution on the coarse - * mesh (mesh_1) and interpolate - * it to the children of the - * corresponding cell of - * mesh_2. Conversely, if the new - * mesh is coarser, one has to - * express the coarse cell shape - * function by a linear - * combination of fine cell shape - * functions. In either case, one - * needs to loop over the finest - * cells that are common to both - * triangulations. This function - * returns a list of pairs of - * matching iterators to cells in - * the two meshes that can be - * used to this end. - * - * Note that the list of these - * iterators is not necessarily - * order, and does also not - * necessarily coincide with the - * order in which cells are - * traversed in one, or both, of - * the meshes given as arguments. - */ - template - static - std::list > - get_finest_common_cells (const Container &mesh_1, - const Container &mesh_2); - - /** - * Return true if the two - * triangulations are based on - * the same coarse mesh. This is - * determined by checking whether - * they have the same number of - * cells on the coarsest level, - * and then checking that they - * have the same vertices. - * - * The two meshes may have - * different refinement histories - * beyond the coarse mesh. - */ - template - static - bool - have_same_coarse_mesh (const Triangulation &mesh_1, - const Triangulation &mesh_2); - - /** - * The same function as above, - * but working on arguments of - * type DoFHandler, - * hp::DoFHandler, or - * MGDoFHandler. This function is - * provided to allow calling - * have_same_coarse_mesh for all - * types of containers - * representing triangulations or - * the classes built on - * triangulations. - */ - template - static - bool - have_same_coarse_mesh (const Container &mesh_1, - const Container &mesh_2); - - /** - * Return the diamater of the smallest - * active cell of a triangulation. See - * step-24 for an example - * of use of this function. - */ - template - static - double - minimal_cell_diameter (const Triangulation &triangulation); - - /** - * Return the diamater of the largest - * active cell of a triangulation. - */ - template - static - double - maximal_cell_diameter (const Triangulation &triangulation); - - /** - * Given the two triangulations - * specified as the first two - * arguments, create the - * triangulation that contains - * the finest cells of both - * triangulation and store it in - * the third parameter. Previous - * content of @p result will be - * deleted. - * - * @note This function is intended - * to create an adaptively refined - * triangulation that contains the - * most refined cells from - * two input triangulations that - * were derived from the same - * coarse grid by adaptive refinement. - * This is an operation sometimes - * needed when one solves for two - * variables of a coupled problem - * on separately refined meshes on - * the same domain (for example - * because these variables have - * boundary layers in different places) - * but then needs to compute something - * that involves both variables or - * wants to output the result into a - * single file. In both cases, in - * order not to lose information, - * the two solutions can not be - * interpolated onto the respectively - * other mesh because that may be - * coarser than the ones on which - * the variable was computed. Rather, - * one needs to have a mesh for the - * domain that is at least as fine - * as each of the two initial meshes. - * This function computes such a mesh. - * - * @note If you want to create - * a mesh that is the merger of - * two other coarse meshes, for - * example in order to compose a mesh - * for a complicated geometry from - * meshes for simpler geometries, - * take a look at - * GridGenerator::merge_triangulations . - */ - template - static - void - create_union_triangulation (const Triangulation &triangulation_1, - const Triangulation &triangulation_2, - Triangulation &result); - - /** - * Given a triangulation and a - * list of cells whose children - * have become distorted as a - * result of mesh refinement, try - * to fix these cells up by - * moving the center node around. - * - * The function returns a list of - * cells with distorted children - * that couldn't be fixed up for - * whatever reason. The returned - * list is therefore a subset of - * the input argument. - * - * For a definition of the - * concept of distorted cells, - * see the - * @ref GlossDistorted "glossary entry". - * The first argument passed to the - * current function is typically - * the exception thrown by the - * Triangulation::execute_coarsening_and_refinement - * function. - */ - template - static - typename Triangulation::DistortedCellList - fix_up_distorted_child_cells (const typename Triangulation::DistortedCellList &distorted_cells, - Triangulation &triangulation); - - /** - * This function implements a boundary - * subgrid extraction. Given a - * -Triangulation (the - * "volume mesh") the function extracts a - * subset of its boundary (the "surface - * mesh"). The boundary to be extracted - * is specified by a list of - * boundary_ids. If none is specified - * the whole boundary will be - * extracted. The function is used in - * step-38. - * - * It also builds a mapping linking the - * cells on the surface mesh to the - * corresponding faces on the volume - * one. This mapping is the return value - * of the function. - * - * @note The function builds the surface - * mesh by creating a coarse mesh from - * the selected faces of the coarse cells - * of the volume mesh. It copies the - * boundary indicators of these faces to - * the cells of the coarse surface - * mesh. The surface mesh is then refined - * in the same way as the faces of the - * volume mesh are. In order to ensure - * that the surface mesh has the same - * vertices as the volume mesh, it is - * therefore important that you assign - * appropriate boundary objects through - * Triangulation::set_boundary to the - * surface mesh object before calling - * this function. If you don't, the - * refinement will happen under the - * assumption that all faces are straight - * (i.e using the StraightBoundary class) - * rather than any curved boundary object - * you may want to use to determine the - * location of new vertices. - * - * @note Oftentimes, the - * Container - * template type will be of kind - * Triangulation; in that case, - * the map that is returned will - * be between Triangulation cell - * iterators of the surface mesh - * and Triangulation face - * iterators of the volume - * mesh. However, one often needs - * to have this mapping between - * DoFHandler (or hp::DoFHandler) - * iterators. In that case, you - * can pass DoFHandler arguments - * as first and second parameter; - * the function will in that case - * re-build the triangulation - * underlying the second argument - * and return a map between - * DoFHandler iterators. However, - * the function will not actually - * distribute degrees of freedom - * on this newly created surface - * mesh. - * - * @note The algorithm outlined - * above assumes that all faces - * on higher refinement levels - * always have exactly the same - * boundary indicator as their - * parent face. Consequently, we - * can start with coarse level - * faces and build the surface - * mesh based on that. It would - * not be very difficult to - * extend the function to also - * copy boundary indicators from - * finer level faces to their - * corresponding surface mesh - * cells, for example to - * accomodate different geometry - * descriptions in the case of - * curved boundaries. - */ - template