From e6083994f9b578539a62645a22dadc822ffabe6f Mon Sep 17 00:00:00 2001 From: wolf Date: Mon, 24 Sep 2001 13:58:28 +0000 Subject: [PATCH] Add Eulerian mapping. git-svn-id: https://svn.dealii.org/trunk@5038 0785d39b-7218-0410-832d-ea1e28bc413d --- .../deal.II/include/fe/mapping_q1_eulerian.h | 151 ++++++++++++++++++ .../deal.II/source/fe/mapping_q1_eulerian.cc | 109 +++++++++++++ 2 files changed, 260 insertions(+) create mode 100644 deal.II/deal.II/include/fe/mapping_q1_eulerian.h create mode 100644 deal.II/deal.II/source/fe/mapping_q1_eulerian.cc diff --git a/deal.II/deal.II/include/fe/mapping_q1_eulerian.h b/deal.II/deal.II/include/fe/mapping_q1_eulerian.h new file mode 100644 index 0000000000..c9b24d0ae6 --- /dev/null +++ b/deal.II/deal.II/include/fe/mapping_q1_eulerian.h @@ -0,0 +1,151 @@ +//---------------------------- mapping_q1_eulerian.h --------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2001 by the deal.II authors and Michael Stadler +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//---------------------------- mapping_q1_eulerian.h --------------------------- +#ifndef __deal2__mapping_q1_eulerian_h +#define __deal2__mapping_q1_eulerian_h + +#include +#include + + +/** + * Eulerian mapping of general unit cells by d-linear shape + * functions. Each cell is thus shifted in space by values given to + * the mapping through a finite element field. + * + * @sect3{Usage} + * + * The constructor of this class takes two arguments: a reference to + * the vector that defines the mapping from the reference + * configuration to the current configuration and a reference to the + * @ref{DoFHandler}. The vector should then represent a (flattened out + * version of a) vector valued field defined at nodes defined by the + * the @ref{DoFHandler}, where the number of components of the vector + * field equals the number of space dimensions. Thus, the + * @ref{DoFHandler} shall operate on a finite element that has as many + * components as space dimensions. As an additional requirement, we + * impose that it have as many degree of freedom per vertex as there + * are space dimensions; since this object only evaluates the finite + * element field at the vertices, the values + * of all other degrees of freedom (not associated to vertices) are + * ignored. These requirements are met if the finite element which the + * given @ref{DoFHandler} operates on is constructed as a system + * element (@ref{FESystem}) from @p{dim} continuous @ref{FE_Q} + * objects. + * + * In many cases, the shift vector will also be the solution vector of + * the problem under investigation. If this is not the case (i.e. the + * number of components of the solution variable is not equal to the + * space dimension, e.g. for scalar problems in @p{dim>1} where the + * Eulerian coordinates only give a background field) or for coupled + * problems where more variables are computed than just the flow + * field), then a different @ref{DoFHandler} has to be set up on the + * given triangulation, and the shift vector has then to be associated + * to it. + * + * An example is shown below: + * @begin{verbatim} + * FESystem fe(FE_Q(1), dim); + * DoFHandler flowfield_dof_handler(triangulation); + * flowfield_dof_handler.distribute_dofs(fe); + * Vector map_points(flowfield_dof_handler.n_dofs()); + * MappingQ1Eulerian mymapping(map_points, flowfield_dof_handler); + * @end{verbatim} + * + * Note that since the vector of shift values and the dof handler are + * only associated to this object at construction time, you have to + * make sure that whenever you use this object, the given objects + * still represent valid data. + * + * @author Michael Stadler, 2001 + */ +template +class MappingQ1Eulerian : public MappingQ1 +{ + public: + + /** + * Constructor. It takes a + * @p{Vector &} as its + * first argument to specify the + * transformation of the whole + * problem from the reference to + * the current configuration. + * The organization of the + * elements in the @p{Vector} + * must follow the concept how + * deal.II stores solutions that + * are associated to a + * triangulation. This is + * automatically the case if the + * @p{Vector} represents the + * solution of the previous step + * of a nonlinear problem. + * Alternatively, the @p{Vector} + * can be initialized by + * @p{DoFObjectAccessor::set_dof_values()}. + */ + MappingQ1Eulerian ( const Vector &euler_transform_vectors, + const DoFHandler &shiftmap_dof_handler); + + + /** + * Exception + */ + DeclException0 (ExcWrongNoOfComponents); + + /** + * Exception. + */ + DeclException2 (ExcWrongVectorSize, + int, int, + << "Vector has wrong size " << arg1 + << ", expected size " << arg2); + + /** + * Exception. + */ + DeclException0 (ExcInactiveCell); + + + + protected: + + /** + * Reference to the vector of + * shifts. + */ + const Vector &euler_transform_vectors; + + /** + * Pointer to the DoFHandler to + * which the mapping vector is + * associated. + */ + const SmartPointer > shiftmap_dof_handler; + + + private: + /** + * Computes the support points of + * the mapping. For + * @p{MappingQ1Eulerian} these + * are the vertices. + */ + virtual void compute_mapping_support_points( + const typename Triangulation::cell_iterator &cell, + std::vector > &a) const; + +}; + + +#endif diff --git a/deal.II/deal.II/source/fe/mapping_q1_eulerian.cc b/deal.II/deal.II/source/fe/mapping_q1_eulerian.cc new file mode 100644 index 0000000000..2df84da7fe --- /dev/null +++ b/deal.II/deal.II/source/fe/mapping_q1_eulerian.cc @@ -0,0 +1,109 @@ +//---------------------------- mapping_q1_eulerian.cc --------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2001 by the deal.II authors and Michael Stadler +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//---------------------------- mapping_q1_eulerian.cc --------------------------- + + +#include +#include +#include +#include +#include +#include + + + +template +MappingQ1Eulerian::MappingQ1Eulerian ( const Vector &euler_transform_vectors, + const DoFHandler &shiftmap_dof_handler) + : + euler_transform_vectors(euler_transform_vectors), + shiftmap_dof_handler(&shiftmap_dof_handler) +{} + + + +template +void +MappingQ1Eulerian::compute_mapping_support_points( + const typename Triangulation::cell_iterator &cell, + std::vector > &a) const +{ + + // The assertions can not be in the + // constructor, since this would + // require to call + // dof_handler.distribute_dofs(fe) + // *before* the mapping object is + // constructed, which is not + // necessarily what we want. + Assert (dim == shiftmap_dof_handler->get_fe().n_dofs_per_vertex(), + ExcWrongNoOfComponents()); + Assert (shiftmap_dof_handler->get_fe().n_components() == dim, + ExcWrongNoOfComponents()); + + Assert (shiftmap_dof_handler->n_dofs() == euler_transform_vectors.size(), + ExcWrongVectorSize(euler_transform_vectors.size(), + shiftmap_dof_handler->n_dofs())); + + // cast the + // Triangulation::cell_iterator + // into a + // DoFHandler::cell_iterator + // which is necessary for access to + // DoFCellAccessor::get_dof_values() + typename DoFHandler::cell_iterator + dof_cell (const_cast *> (&(cell->get_triangulation())), + cell->level(), + cell->index(), + shiftmap_dof_handler); + + // We require the cell to be + // active. This is determined by + // the user when looping over all + // active cells in the problem + // code. + Assert (dof_cell->active() == true, + ExcInactiveCell()); + + // for Q1 elements, the number of + // support points should equal the + // number of vertices + a.resize(GeometryInfo::vertices_per_cell); + + // now get the values of the shift + // vectors at the vertices + Vector mapping_values (shiftmap_dof_handler->get_fe().dofs_per_cell); + dof_cell->get_dof_values (euler_transform_vectors, mapping_values); + + + for (unsigned int i=0; i::vertices_per_cell; ++i) + { + Point shift_vector; + + // pick out the value of the + // shift vector at the present + // vertex. since vertex dofs + // are always numbered first, + // we can access them easily + for (unsigned int j=0; jvertex(vertex_mapping[i]) + shift_vector; + } +} + + +// explicit instantiation +template class MappingQ1Eulerian; -- 2.39.5