From 990855c91fa20e6f22f490b9d7a0a9c251d2387a Mon Sep 17 00:00:00 2001 From: kanschat Date: Tue, 6 Sep 2011 15:43:55 +0000 Subject: [PATCH] compatibility file for last move git-svn-id: https://svn.dealii.org/trunk@24270 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/include/deal.II/meshworker/dof_info.h | 1184 ---------------- .../deal.II/meshworker/integration_info.h | 1211 +++++++++++++++++ .../deal.II/numerics/mesh_worker_info.h | 20 + 3 files changed, 1231 insertions(+), 1184 deletions(-) create mode 100644 deal.II/include/deal.II/meshworker/integration_info.h create mode 100644 deal.II/include/deal.II/numerics/mesh_worker_info.h diff --git a/deal.II/include/deal.II/meshworker/dof_info.h b/deal.II/include/deal.II/meshworker/dof_info.h index 3a8f6d43db..3eb644ec73 100644 --- a/deal.II/include/deal.II/meshworker/dof_info.h +++ b/deal.II/include/deal.II/meshworker/dof_info.h @@ -18,7 +18,6 @@ #include #include #include -//#include #include #include @@ -280,676 +279,6 @@ namespace MeshWorker bool exterior_face_available[GeometryInfo::faces_per_cell]; }; - - - -/** - * Class for objects handed to local integration functions. - * - * Objects of this class contain one or more objects of type FEValues, - * FEFaceValues or FESubfacevalues to be used in local - * integration. They are stored in an array of pointers to the base - * classes FEValuesBase for cells and FEFaceValuesBase for faces and - * subfaces, respectively. The template parameter VECTOR allows the - * use of different data types for the global system. - * - * The @p FEVALUESBASE template parameter should be either - * FEValuesBase or FEFaceValuesBase, depending on whether the object - * is used to integrate over cells or faces. The actual type of @p - * FEVALUES object is fixed in the constructor and only used to - * initialize the pointers in #fevalv. - * - * Additionally, this function containes space to store the values of - * finite element functions stored in #global_data in the - * quadrature points. These vectors are initialized automatically on - * each cell or face. In order to avoid initializing unused vectors, - * you can use initialize_selector() to select the vectors by name - * that you actually want to use. - * - *

Integration models

- * - * This class supports two local integration models, corresponding to - * the data models in the documentation of the Assembler namespace. - * One is the - * standard model suggested by the use of FESystem. Namely, there is - * one FEValuseBase object in this class, containing all shape - * functions of the whole system, and having as many components as the - * system. Using this model involves loops over all system shape - * functions. It requires to identify the system components - * for each shape function and to select the correct bilinear form, - * usually in an @p if or @p switch statement. - * - * The second integration model builds one FEValuesBase object per - * base element of the system. The degrees of freedom on each cell are - * renumbered by block, such that they represent the same block - * structure as the global system. Objects performing the integration - * can then process each block separately, which improves reusability - * of code considerably. - * - * @note As described in DoFInfo, the use of the local block model is - * triggered by calling BlockInfo::initialize_local() before - * using initialize() in this class. - * - * @ingroup MeshWorker - * @author Guido Kanschat, 2009 - */ - template - class IntegrationInfo - { - private: - /// vector of FEValues objects - std::vector > > fevalv; - public: - static const unsigned int dimension = dim; - static const unsigned int space_dimension = spacedim; - - /** - * Constructor. - */ - IntegrationInfo(); - - /** - * Copy constructor, creating a - * clone to be used by - * WorksTream::run(). - */ - IntegrationInfo(const IntegrationInfo& other); - - /** - * Build all internal - * structures, in particular - * the FEValuesBase objects - * and allocate space for - * data vectors. - * - * @param el is the finite - * element of the DoFHandler. - * - * @param mapping is the Mapping - * object used to map the - * mesh cells. - * - * @param quadrature is a - * Quadrature formula used in - * the constructor of the - * FEVALUES objects. - * - * @param flags are the - * UpdateFlags used in - * the constructor of the - * FEVALUES objects. - * - * @param local_block_info is - * an optional parameter for - * systems of PDE. If it is - * provided with reasonable - * data, then the degrees of - * freedom on the cells will be - * re-ordered to reflect the - * block structure of the system. - */ - template - void initialize(const FiniteElement& el, - const Mapping& mapping, - const Quadrature& quadrature, - const UpdateFlags flags, - const BlockInfo* local_block_info = 0); - - /** - * Initialize the data - * vector and cache the - * selector. - */ - void initialize_data(const std_cxx1x::shared_ptr > &data); - - /** - * Delete the data created by initialize(). - */ - void clear(); - - /// This is true if we are assembling for multigrid - bool multigrid; - /// Access to finite element - /** - * This is the access - * function being used, if - * the constructor for a - * single element was - * used. It throws an - * exception, if applied to a - * vector of elements. - */ - const FEValuesBase& fe_values () const; - - /// Access to finite elements - /** - * This access function must - * be used if the constructor - * for a group of elements - * was used. - * - * @see DGBlockSplitApplication - */ - const FEValuesBase& fe_values (const unsigned int i) const; - - /** - * The vector containing the - * values of finite element - * functions in the quadrature - * points. - * - * There is one vector per - * selected finite element - * function, containing one - * vector for each component, - * containing vectors with - * values for each quadrature - * point. - */ - std::vector > > values; - - /** - * The vector containing the - * derivatives of finite - * element functions in the - * quadrature points. - * - * There is one vector per - * selected finite element - * function, containing one - * vector for each component, - * containing vectors with - * values for each quadrature - * point. - */ - std::vector > > > gradients; - - /** - * The vector containing the - * second derivatives of finite - * element functions in the - * quadrature points. - * - * There is one vector per - * selected finite element - * function, containing one - * vector for each component, - * containing vectors with - * values for each quadrature - * point. - */ - std::vector > > > hessians; - - /** - * Reinitialize internal data - * structures for use on a cell. - */ - void reinit(const DoFInfo& i); - - /** - * Use the finite element - * functions in #global_data - * and fill the vectors - * #values, #gradients and - * #hessians. - */ - template - void fill_local_data(const DoFInfo& info, bool split_fevalues); - - /** - * The global data vector - * used to compute function - * values in quadrature - * points. - */ - std_cxx1x::shared_ptr > global_data; - - /** - * The memory used by this object. - */ - std::size_t memory_consumption () const; - - private: - /** - * Use the finite element - * functions in #global_data - * and fill the vectors - * #values, #gradients and - * #hessians with values - * according to the - * selector. - */ - template - void fill_local_data( - std::vector > >& data, - VectorSelector& selector, - bool split_fevalues) const; - /** - * Cache the number of - * components of the system element. - */ - unsigned int n_components; - }; - -/** - * The object holding the scratch data for integrating over cells and - * faces. IntegrationInfoBox serves three main purposes: - * - *
    - *
  1. It provides the interface needed by MeshWorker::loop(), namely - * the two functions post_cell() and post_faces(), as well as - * the data members #cell, #boundary, #face, - * #subface, and #neighbor. - * - *
  2. It contains all information needed to initialize the FEValues - * and FEFaceValues objects in the IntegrationInfo data members. - * - *
  3. It stores information on finite element vectors and whether - * their data should be used to compute values or derivatives of - * functions at quadrature points. - * - *
  4. It makes educated guesses on quadrature rules and update - * flags, so that minimal code has to be written when default - * parameters are sufficient. - *
- * - * In order to allow for sufficient generality, a few steps have to be - * undertaken to use this class. - * - * First, you should consider if you need values from any vectors in a - * NamedData object. If so, fill the VectorSelector objects - * #cell_selector, #boundary_selector and #face_selector with their names - * and the data type (value, gradient, Hessian) to be extracted. - * - * Afterwards, you will need to consider UpdateFlags for FEValues - * objects. A good start is initialize_update_flags(), which looks at - * the selectors filled before and adds all the flags needed to get - * the selection. Additional flags can be set with add_update_flags(). - * - * Finally, we need to choose quadrature formulas. In the simplest - * case, you might be happy with the default settings, which are - * n-point Gauss formulas. If only derivatives of the shape - * functions are used (#update_values is not set) n equals the - * highest polynomial degree in the FiniteElement, if #update_values - * is set, n is one higher than this degree. If you choose to - * use Gauss formulas of other size, use initialize_gauss_quadrature() - * with appropriate values. Otherwise, you can fill the variables - * #cell_quadrature, #boundary_quadrature and #face_quadrature - * directly. - * - * In order to save time, you can set the variables boundary_fluxes - * and interior_fluxes of the base class to false, thus telling the - * Meshworker::loop() not to loop over those faces. - * - * All the information in here is used to set up IntegrationInfo - * objects correctly, typically in an IntegrationInfoBox. - * - * @ingroup MeshWorker - * @author Guido Kanschat, 2009 - */ - template - class IntegrationInfoBox - { - public: - - /** - * The type of the info object - * for cells. - */ - typedef IntegrationInfo CellInfo; - - /** - * Default constructor. - */ - IntegrationInfoBox (); - - /** - * Initialize the - * IntegrationInfo objects - * contained. - * - * Before doing so, add update - * flags necessary to produce - * the data needed and also - * set uninitialized quadrature - * rules to Gauss formulas, - * which integrate polynomial - * bilinear forms exactly. - */ - void initialize(const FiniteElement& el, - const Mapping& mapping, - const BlockInfo* block_info = 0); - - /** - * Initialize the - * IntegrationInfo objects - * contained. - * - * Before doing so, add update - * flags necessary to produce - * the data needed and also - * set uninitialized quadrature - * rules to Gauss formulas, - * which integrate polynomial - * bilinear forms exactly. - */ - template - void initialize(const FiniteElement& el, - const Mapping& mapping, - const NamedData& data, - const BlockInfo* block_info = 0); - - /** - * Initialize the - * IntegrationInfo objects - * contained. - * - * Before doing so, add update - * flags necessary to produce - * the data needed and also - * set uninitialized quadrature - * rules to Gauss formulas, - * which integrate polynomial - * bilinear forms exactly. - */ - template - void initialize(const FiniteElement& el, - const Mapping& mapping, - const NamedData*>& data, - const BlockInfo* block_info = 0); - /** - * @name FEValues setup - */ - /* @{ */ - void initialize_update_flags(); - - /** - * Add FEValues UpdateFlags for - * integration on all objects - * (cells, boundary faces and - * all interior faces). - */ - void add_update_flags_all (const UpdateFlags flags); - - /** - * Add FEValues UpdateFlags for - * integration on cells. - */ - void add_update_flags_cell(const UpdateFlags flags); - - /** - * Add FEValues UpdateFlags for - * integration on boundary faces. - */ - void add_update_flags_boundary(const UpdateFlags flags); - - /** - * Add FEValues UpdateFlags for - * integration on interior faces. - */ - void add_update_flags_face(const UpdateFlags flags); - - /** - * Add additional update flags - * to the ones already set in - * this program. The four - * boolean flags indicate - * wether the additional flags - * should be set for cell, - * boundary, interelement face - * for the cell itself - * or neighbor cell, or - * any combination thereof. - */ - void add_update_flags(const UpdateFlags flags, - const bool cell = true, - const bool boundary = true, - const bool face = true, - const bool neighbor = true); - - /** - * Assign n-point Gauss - * quadratures to each of the - * quadrature rules. Here, a - * size of zero points means - * that no loop over these grid - * entities should be - * performed. - * - * If the parameter - * force is true, then - * all quadrature sets are - * filled with new quadrature - * ruels. If it is false, then - * only empty rules are changed. - */ - void initialize_gauss_quadrature(unsigned int n_cell_points, - unsigned int n_boundary_points, - unsigned int n_face_points, - const bool force = true); - - /** - * The memory used by this object. - */ - std::size_t memory_consumption () const; - - /** - * The set of update flags - * for boundary cell integration. - * - * Defaults to - * #update_JxW_values. - */ - UpdateFlags cell_flags; - /** - * The set of update flags - * for boundary face integration. - * - * Defaults to - * #update_JxW_values and - * #update_normal_vectors. - */ - UpdateFlags boundary_flags; - - /** - * The set of update flags - * for interior face integration. - * - * Defaults to - * #update_JxW_values and - * #update_normal_vectors. - */ - UpdateFlags face_flags; - - /** - * The set of update flags - * for interior face integration. - * - * Defaults to - * #update_default, since - * quadrature weights are - * taken from the other cell. - */ - UpdateFlags neighbor_flags; - - /** - * The quadrature rule used - * on cells. - */ - Quadrature cell_quadrature; - - /** - * The quadrature rule used - * on boundary faces. - */ - Quadrature boundary_quadrature; - - /** - * The quadrature rule used - * on interior faces. - */ - Quadrature face_quadrature; - /* @} */ - - /** - * @name Data vectors - */ - /* @{ */ - - /** - * Initialize the - * VectorSelector objects - * #cell_selector, - * #boundary_selector and - * #face_selector in order to - * save computational - * eeffort. If no selectors - * are used, then values for - * all named vectors in - * DoFInfo::global_data will be - * computed in all quadrature - * points. - * - * This function will also - * add UpdateFlags to the - * flags stored in this class. - */ - /** - * Select the vectors from - * DoFInfo::global_data - * that should be computed in - * the quadrature points on cells. - */ - MeshWorker::VectorSelector cell_selector; - - /** - * Select the vectors from - * DoFInfo::global_data - * that should be computed in - * the quadrature points on - * boundary faces. - */ - MeshWorker::VectorSelector boundary_selector; - - /** - * Select the vectors from - * DoFInfo::global_data - * that should be computed in - * the quadrature points on - * interior faces. - */ - MeshWorker::VectorSelector face_selector; - - std_cxx1x::shared_ptr > cell_data; - std_cxx1x::shared_ptr > boundary_data; - std_cxx1x::shared_ptr > face_data; - /* @} */ - - /** - * @name Interface for MeshWorker::loop() - */ - /* @{ */ - /** - * A callback function which is - * called in the loop over all - * cells, after the action on a - * cell has been performed and - * before the faces are dealt - * with. - * - * In order for this function - * to have this effect, - * at least either of the - * arguments - * boundary_worker or - * face_worker - * arguments of loop() should - * be nonzero. Additionally, - * cells_first should - * be true. If - * cells_first is - * false, this function is - * called before any action on - * a cell is taken. - * - * And empty function in this - * class, but can be replaced - * in other classes given to - * loop() instead. - * - * See loop() and cell_action() - * for more details of how this - * function can be used. - */ - template - void post_cell(const DoFInfoBox&); - - /** - * A callback function which is - * called in the loop over all - * cells, after the action on - * the faces of a cell has been - * performed and before the - * cell itself is dealt with - * (assumes - * cells_first is false). - * - * In order for this function - * to have a reasonable effect, - * at least either of the - * arguments - * boundary_worker or - * face_worker - * arguments of loop() should - * be nonzero. Additionally, - * cells_first should - * be false. - * - * And empty function in this - * class, but can be replaced - * in other classes given to - * loop() instead. - * - * See loop() and cell_action() - * for more details of how this - * function can be used. - */ - template - void post_faces(const DoFInfoBox&); - - /** - * The info object for a cell. - */ - CellInfo cell; - /** - * The info object for a - * boundary face. - */ - CellInfo boundary; - /** - * The info object for a - * regular interior face, seen - * from the first cell. - */ - CellInfo face; - /** - * The info object for the - * refined side of an interior - * face seen from the first - * cell. - */ - CellInfo subface; - /** - * The info object for an - * interior face, seen from the - * other cell. - */ - CellInfo neighbor; - - /* @} */ - }; - - //----------------------------------------------------------------------// template @@ -1124,519 +453,6 @@ namespace MeshWorker } } } - - -//----------------------------------------------------------------------// - - template - inline - IntegrationInfo::IntegrationInfo() - : - fevalv(0), - multigrid(false), - global_data(std_cxx1x::shared_ptr >(new VectorDataBase)) - {} - - - template - inline - IntegrationInfo::IntegrationInfo(const IntegrationInfo& other) - : - multigrid(other.multigrid), - values(other.values), - gradients(other.gradients), - hessians(other.hessians), - global_data(other.global_data), - n_components(other.n_components) - { - fevalv.resize(other.fevalv.size()); - for (unsigned int i=0;i& p = *other.fevalv[i]; - const FEValues* pc = dynamic_cast*>(&p); - const FEFaceValues* pf = dynamic_cast*>(&p); - const FESubfaceValues* ps = dynamic_cast*>(&p); - - if (pc != 0) - fevalv[i] = std_cxx1x::shared_ptr > ( - reinterpret_cast*>( - new FEValues (pc->get_mapping(), pc->get_fe(), - pc->get_quadrature(), pc->get_update_flags()))); - else if (pf != 0) - fevalv[i] = std_cxx1x::shared_ptr > ( - new FEFaceValues (pf->get_mapping(), pf->get_fe(), pf->get_quadrature(), pf->get_update_flags())); - else if (ps != 0) - fevalv[i] = std_cxx1x::shared_ptr > ( - new FESubfaceValues (ps->get_mapping(), ps->get_fe(), ps->get_quadrature(), ps->get_update_flags())); - else - Assert(false, ExcInternalError()); - } - } - - - template <> - inline - IntegrationInfo<1,1>::IntegrationInfo(const IntegrationInfo<1,1>& other) - : - multigrid(other.multigrid), - values(other.values), - gradients(other.gradients), - hessians(other.hessians), - global_data(other.global_data), - n_components(other.n_components) - { - const int dim = 1; - const int sdim = 1; - - fevalv.resize(other.fevalv.size()); - for (unsigned int i=0;i& p = *other.fevalv[i]; - const FEValues* pc = dynamic_cast*>(&p); - const FEFaceValues* pf = dynamic_cast*>(&p); - const FESubfaceValues* ps = dynamic_cast*>(&p); - - if (pc != 0) - fevalv[i] = std_cxx1x::shared_ptr > ( - reinterpret_cast*>( - new FEValues (pc->get_mapping(), pc->get_fe(), - pc->get_quadrature(), pc->get_update_flags()))); - else if (pf != 0) - { - Assert (false, ExcImpossibleInDim(1)); - fevalv[i].reset (); - } - else if (ps != 0) - { - Assert (false, ExcImpossibleInDim(1)); - fevalv[i].reset(); - } - else - Assert(false, ExcInternalError()); - } - } - - - template <> - inline - IntegrationInfo<1,2>::IntegrationInfo(const IntegrationInfo<1,2>& other) - : - multigrid(other.multigrid), - values(other.values), - gradients(other.gradients), - hessians(other.hessians), - global_data(other.global_data), - n_components(other.n_components) - { - const int dim = 1; - const int sdim = 2; - - fevalv.resize(other.fevalv.size()); - for (unsigned int i=0;i& p = *other.fevalv[i]; - const FEValues* pc = dynamic_cast*>(&p); - const FEFaceValues* pf = dynamic_cast*>(&p); - const FESubfaceValues* ps = dynamic_cast*>(&p); - - if (pc != 0) - fevalv[i] = std_cxx1x::shared_ptr > ( - reinterpret_cast*>( - new FEValues (pc->get_mapping(), pc->get_fe(), - pc->get_quadrature(), pc->get_update_flags()))); - else if (pf != 0) - { - Assert (false, ExcImpossibleInDim(1)); - fevalv[i].reset(); - } - else if (ps != 0) - { - Assert (false, ExcImpossibleInDim(1)); - fevalv[i].reset(); - } - else - Assert(false, ExcInternalError()); - } - } - - - - template - template - inline void - IntegrationInfo::initialize( - const FiniteElement& el, - const Mapping& mapping, - const Quadrature& quadrature, - const UpdateFlags flags, - const BlockInfo* block_info) - { - if (block_info == 0 || block_info->local().size() == 0) - { - fevalv.resize(1); - fevalv[0] = std_cxx1x::shared_ptr > ( - new FEVALUES (mapping, el, quadrature, flags)); - } - else - { - fevalv.resize(el.n_base_elements()); - for (unsigned int i=0;i > ( - new FEVALUES (mapping, el.base_element(i), quadrature, flags)); - } - } - n_components = el.n_components(); - } - - - template - inline const FEValuesBase& - IntegrationInfo::fe_values() const - { - AssertDimension(fevalv.size(), 1); - return *fevalv[0]; - } - - - template - inline const FEValuesBase& - IntegrationInfo::fe_values(unsigned int i) const - { - Assert (i - inline void - IntegrationInfo::reinit(const DoFInfo& info) - { - for (unsigned int i=0;i& febase = *fevalv[i]; - if (info.sub_number != deal_II_numbers::invalid_unsigned_int) - { - // This is a subface - FESubfaceValues& fe = dynamic_cast&> (febase); - fe.reinit(info.cell, info.face_number, info.sub_number); - } - else if (info.face_number != deal_II_numbers::invalid_unsigned_int) - { - // This is a face - FEFaceValues& fe = dynamic_cast&> (febase); - fe.reinit(info.cell, info.face_number); - } - else - { - // This is a cell - FEValues& fe = dynamic_cast&> (febase); - fe.reinit(info.cell); - } - } - - const bool split_fevalues = info.block_info != 0; - if (!global_data->empty()) - fill_local_data(info, split_fevalues); - } - - - template <> - inline void - IntegrationInfo<1,1>::reinit(const DoFInfo<1,1>& info) - { - const int dim = 1; - const int spacedim = 1; - - for (unsigned int i=0;i& febase = *fevalv[i]; - if (info.sub_number != deal_II_numbers::invalid_unsigned_int) - { - // This is a subface - Assert (false, ExcImpossibleInDim(1)); - } - else if (info.face_number != deal_II_numbers::invalid_unsigned_int) - { - // This is a face - Assert (false, ExcImpossibleInDim(1)); - } - else - { - // This is a cell - FEValues& fe = dynamic_cast&> (febase); - fe.reinit(info.cell); - } - } - - const bool split_fevalues = info.block_info != 0; - if (!global_data->empty()) - fill_local_data(info, split_fevalues); - } - - - template <> - inline void - IntegrationInfo<1,2>::reinit(const DoFInfo<1,2>& info) - { - const int dim = 1; - const int spacedim = 2; - - for (unsigned int i=0;i& febase = *fevalv[i]; - if (info.sub_number != deal_II_numbers::invalid_unsigned_int) - { - // This is a subface - Assert (false, ExcImpossibleInDim(1)); - } - else if (info.face_number != deal_II_numbers::invalid_unsigned_int) - { - // This is a face - Assert (false, ExcImpossibleInDim(1)); - } - else - { - // This is a cell - FEValues& fe = dynamic_cast&> (febase); - fe.reinit(info.cell); - } - } - - const bool split_fevalues = info.block_info != 0; - if (!global_data->empty()) - fill_local_data(info, split_fevalues); - } - - -//----------------------------------------------------------------------// - - template <> - inline - void - IntegrationInfoBox<1,1>::initialize_gauss_quadrature( - const unsigned int cp, - const unsigned int, - const unsigned int, - bool force) - { - if (force || cell_quadrature.size() == 0) - cell_quadrature = QGauss<1>(cp); - } - - - template <> - inline - void - IntegrationInfoBox<1,2>::initialize_gauss_quadrature( - const unsigned int cp, - const unsigned int, - const unsigned int, - bool force) - { - if (force || cell_quadrature.size() == 0) - cell_quadrature = QGauss<1>(cp); - } - - - template - inline - void - IntegrationInfoBox::initialize_gauss_quadrature( - unsigned int cp, - unsigned int bp, - unsigned int fp, - bool force) - { - if (force || cell_quadrature.size() == 0) - cell_quadrature = QGauss(cp); - if (force || boundary_quadrature.size() == 0) - boundary_quadrature = QGauss(bp); - if (force || face_quadrature.size() == 0) - face_quadrature = QGauss(fp); - } - - - template - inline - void - IntegrationInfoBox::add_update_flags_all (const UpdateFlags flags) - { - add_update_flags(flags, true, true, true, true); - } - - - template - inline - void - IntegrationInfoBox::add_update_flags_cell (const UpdateFlags flags) - { - add_update_flags(flags, true, false, false, false); - } - - - template - inline - void - IntegrationInfoBox::add_update_flags_boundary (const UpdateFlags flags) - { - add_update_flags(flags, false, true, false, false); - } - - - template - inline - void - IntegrationInfoBox::add_update_flags_face (const UpdateFlags flags) - { - add_update_flags(flags, false, false, true, true); - } - - - template <> - inline - void - IntegrationInfoBox<1,1>::initialize( - const FiniteElement<1,1>& el, - const Mapping<1,1>& mapping, - const BlockInfo* block_info) - { - initialize_update_flags(); - initialize_gauss_quadrature( - (cell_flags & update_values) ? (el.tensor_degree()+1) : el.tensor_degree(), 1, 1, false); - - const int dim = 1; - const int sdim = 1; - - cell.initialize >(el, mapping, cell_quadrature, - cell_flags, block_info); - } - - - - template <> - inline - void - IntegrationInfoBox<1,2>::initialize( - const FiniteElement<1,2>& el, - const Mapping<1,2>& mapping, - const BlockInfo* block_info) - { - initialize_update_flags(); - initialize_gauss_quadrature( - (cell_flags & update_values) ? (el.tensor_degree()+1) : el.tensor_degree(), 1, 1, false); - - const int dim = 1; - const int sdim = 2; - - cell.initialize >(el, mapping, cell_quadrature, - cell_flags, block_info); - } - - - template - inline - void - IntegrationInfoBox::initialize( - const FiniteElement& el, - const Mapping& mapping, - const BlockInfo* block_info) - { - initialize_update_flags(); - initialize_gauss_quadrature( - (cell_flags & update_values) ? (el.tensor_degree()+1) : el.tensor_degree(), - (boundary_flags & update_values) ? (el.tensor_degree()+1) : el.tensor_degree(), - (face_flags & update_values) ? (el.tensor_degree()+1) : el.tensor_degree(), false); - - cell.template initialize >(el, mapping, cell_quadrature, - cell_flags, block_info); - boundary.template initialize >(el, mapping, boundary_quadrature, - boundary_flags, block_info); - face.template initialize >(el, mapping, face_quadrature, - face_flags, block_info); - subface.template initialize >(el, mapping, face_quadrature, - face_flags, block_info); - neighbor.template initialize >(el, mapping, face_quadrature, - neighbor_flags, block_info); - } - - - template - template - void - IntegrationInfoBox::initialize( - const FiniteElement& el, - const Mapping& mapping, - const NamedData& data, - const BlockInfo* block_info) - { - initialize(el, mapping, block_info); - std_cxx1x::shared_ptr > p; - - p = std_cxx1x::shared_ptr >(new VectorData (cell_selector)); - p->initialize(data); - cell_data = p; - cell.initialize_data(p); - - p = std_cxx1x::shared_ptr >(new VectorData (boundary_selector)); - p->initialize(data); - boundary_data = p; - boundary.initialize_data(p); - - p = std_cxx1x::shared_ptr >(new VectorData (face_selector)); - p->initialize(data); - face_data = p; - face.initialize_data(p); - subface.initialize_data(p); - neighbor.initialize_data(p); - } - - - template - template - void - IntegrationInfoBox::initialize( - const FiniteElement& el, - const Mapping& mapping, - const NamedData*>& data, - const BlockInfo* block_info) - { - initialize(el, mapping, block_info); - std_cxx1x::shared_ptr > p; - - p = std_cxx1x::shared_ptr >(new MGVectorData (cell_selector)); - p->initialize(data); - cell_data = p; - cell.initialize_data(p); - - p = std_cxx1x::shared_ptr >(new MGVectorData (boundary_selector)); - p->initialize(data); - boundary_data = p; - boundary.initialize_data(p); - - p = std_cxx1x::shared_ptr >(new MGVectorData (face_selector)); - p->initialize(data); - face_data = p; - face.initialize_data(p); - subface.initialize_data(p); - neighbor.initialize_data(p); - } - - - template - template - void - IntegrationInfoBox::post_cell(const DoFInfoBox&) - {} - - - template - template - void - IntegrationInfoBox::post_faces(const DoFInfoBox&) - {} - - } DEAL_II_NAMESPACE_CLOSE diff --git a/deal.II/include/deal.II/meshworker/integration_info.h b/deal.II/include/deal.II/meshworker/integration_info.h new file mode 100644 index 0000000000..3cadf27adf --- /dev/null +++ b/deal.II/include/deal.II/meshworker/integration_info.h @@ -0,0 +1,1211 @@ +//--------------------------------------------------------------------------- +// $Id$ +// +// Copyright (C) 2006, 2007, 2008, 2009, 2010, 2011 by the deal.II authors +// +// 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. +// +//--------------------------------------------------------------------------- + +#ifndef __deal2__mesh_worker_info_h +#define __deal2__mesh_worker_info_h + +#include +#include +#include +#include +#include +#include +#include +#include + +DEAL_II_NAMESPACE_OPEN + +namespace MeshWorker +{ +/** + * Class for objects handed to local integration functions. + * + * Objects of this class contain one or more objects of type FEValues, + * FEFaceValues or FESubfacevalues to be used in local + * integration. They are stored in an array of pointers to the base + * classes FEValuesBase for cells and FEFaceValuesBase for faces and + * subfaces, respectively. The template parameter VECTOR allows the + * use of different data types for the global system. + * + * The @p FEVALUESBASE template parameter should be either + * FEValuesBase or FEFaceValuesBase, depending on whether the object + * is used to integrate over cells or faces. The actual type of @p + * FEVALUES object is fixed in the constructor and only used to + * initialize the pointers in #fevalv. + * + * Additionally, this function containes space to store the values of + * finite element functions stored in #global_data in the + * quadrature points. These vectors are initialized automatically on + * each cell or face. In order to avoid initializing unused vectors, + * you can use initialize_selector() to select the vectors by name + * that you actually want to use. + * + *

Integration models

+ * + * This class supports two local integration models, corresponding to + * the data models in the documentation of the Assembler namespace. + * One is the + * standard model suggested by the use of FESystem. Namely, there is + * one FEValuseBase object in this class, containing all shape + * functions of the whole system, and having as many components as the + * system. Using this model involves loops over all system shape + * functions. It requires to identify the system components + * for each shape function and to select the correct bilinear form, + * usually in an @p if or @p switch statement. + * + * The second integration model builds one FEValuesBase object per + * base element of the system. The degrees of freedom on each cell are + * renumbered by block, such that they represent the same block + * structure as the global system. Objects performing the integration + * can then process each block separately, which improves reusability + * of code considerably. + * + * @note As described in DoFInfo, the use of the local block model is + * triggered by calling BlockInfo::initialize_local() before + * using initialize() in this class. + * + * @ingroup MeshWorker + * @author Guido Kanschat, 2009 + */ + template + class IntegrationInfo + { + private: + /// vector of FEValues objects + std::vector > > fevalv; + public: + static const unsigned int dimension = dim; + static const unsigned int space_dimension = spacedim; + + /** + * Constructor. + */ + IntegrationInfo(); + + /** + * Copy constructor, creating a + * clone to be used by + * WorksTream::run(). + */ + IntegrationInfo(const IntegrationInfo& other); + + /** + * Build all internal + * structures, in particular + * the FEValuesBase objects + * and allocate space for + * data vectors. + * + * @param el is the finite + * element of the DoFHandler. + * + * @param mapping is the Mapping + * object used to map the + * mesh cells. + * + * @param quadrature is a + * Quadrature formula used in + * the constructor of the + * FEVALUES objects. + * + * @param flags are the + * UpdateFlags used in + * the constructor of the + * FEVALUES objects. + * + * @param local_block_info is + * an optional parameter for + * systems of PDE. If it is + * provided with reasonable + * data, then the degrees of + * freedom on the cells will be + * re-ordered to reflect the + * block structure of the system. + */ + template + void initialize(const FiniteElement& el, + const Mapping& mapping, + const Quadrature& quadrature, + const UpdateFlags flags, + const BlockInfo* local_block_info = 0); + + /** + * Initialize the data + * vector and cache the + * selector. + */ + void initialize_data(const std_cxx1x::shared_ptr > &data); + + /** + * Delete the data created by initialize(). + */ + void clear(); + + /// This is true if we are assembling for multigrid + bool multigrid; + /// Access to finite element + /** + * This is the access + * function being used, if + * the constructor for a + * single element was + * used. It throws an + * exception, if applied to a + * vector of elements. + */ + const FEValuesBase& fe_values () const; + + /// Access to finite elements + /** + * This access function must + * be used if the constructor + * for a group of elements + * was used. + * + * @see DGBlockSplitApplication + */ + const FEValuesBase& fe_values (const unsigned int i) const; + + /** + * The vector containing the + * values of finite element + * functions in the quadrature + * points. + * + * There is one vector per + * selected finite element + * function, containing one + * vector for each component, + * containing vectors with + * values for each quadrature + * point. + */ + std::vector > > values; + + /** + * The vector containing the + * derivatives of finite + * element functions in the + * quadrature points. + * + * There is one vector per + * selected finite element + * function, containing one + * vector for each component, + * containing vectors with + * values for each quadrature + * point. + */ + std::vector > > > gradients; + + /** + * The vector containing the + * second derivatives of finite + * element functions in the + * quadrature points. + * + * There is one vector per + * selected finite element + * function, containing one + * vector for each component, + * containing vectors with + * values for each quadrature + * point. + */ + std::vector > > > hessians; + + /** + * Reinitialize internal data + * structures for use on a cell. + */ + void reinit(const DoFInfo& i); + + /** + * Use the finite element + * functions in #global_data + * and fill the vectors + * #values, #gradients and + * #hessians. + */ + template + void fill_local_data(const DoFInfo& info, bool split_fevalues); + + /** + * The global data vector + * used to compute function + * values in quadrature + * points. + */ + std_cxx1x::shared_ptr > global_data; + + /** + * The memory used by this object. + */ + std::size_t memory_consumption () const; + + private: + /** + * Use the finite element + * functions in #global_data + * and fill the vectors + * #values, #gradients and + * #hessians with values + * according to the + * selector. + */ + template + void fill_local_data( + std::vector > >& data, + VectorSelector& selector, + bool split_fevalues) const; + /** + * Cache the number of + * components of the system element. + */ + unsigned int n_components; + }; + +/** + * The object holding the scratch data for integrating over cells and + * faces. IntegrationInfoBox serves three main purposes: + * + *
    + *
  1. It provides the interface needed by MeshWorker::loop(), namely + * the two functions post_cell() and post_faces(), as well as + * the data members #cell, #boundary, #face, + * #subface, and #neighbor. + * + *
  2. It contains all information needed to initialize the FEValues + * and FEFaceValues objects in the IntegrationInfo data members. + * + *
  3. It stores information on finite element vectors and whether + * their data should be used to compute values or derivatives of + * functions at quadrature points. + * + *
  4. It makes educated guesses on quadrature rules and update + * flags, so that minimal code has to be written when default + * parameters are sufficient. + *
+ * + * In order to allow for sufficient generality, a few steps have to be + * undertaken to use this class. + * + * First, you should consider if you need values from any vectors in a + * NamedData object. If so, fill the VectorSelector objects + * #cell_selector, #boundary_selector and #face_selector with their names + * and the data type (value, gradient, Hessian) to be extracted. + * + * Afterwards, you will need to consider UpdateFlags for FEValues + * objects. A good start is initialize_update_flags(), which looks at + * the selectors filled before and adds all the flags needed to get + * the selection. Additional flags can be set with add_update_flags(). + * + * Finally, we need to choose quadrature formulas. In the simplest + * case, you might be happy with the default settings, which are + * n-point Gauss formulas. If only derivatives of the shape + * functions are used (#update_values is not set) n equals the + * highest polynomial degree in the FiniteElement, if #update_values + * is set, n is one higher than this degree. If you choose to + * use Gauss formulas of other size, use initialize_gauss_quadrature() + * with appropriate values. Otherwise, you can fill the variables + * #cell_quadrature, #boundary_quadrature and #face_quadrature + * directly. + * + * In order to save time, you can set the variables boundary_fluxes + * and interior_fluxes of the base class to false, thus telling the + * Meshworker::loop() not to loop over those faces. + * + * All the information in here is used to set up IntegrationInfo + * objects correctly, typically in an IntegrationInfoBox. + * + * @ingroup MeshWorker + * @author Guido Kanschat, 2009 + */ + template + class IntegrationInfoBox + { + public: + + /** + * The type of the info object + * for cells. + */ + typedef IntegrationInfo CellInfo; + + /** + * Default constructor. + */ + IntegrationInfoBox (); + + /** + * Initialize the + * IntegrationInfo objects + * contained. + * + * Before doing so, add update + * flags necessary to produce + * the data needed and also + * set uninitialized quadrature + * rules to Gauss formulas, + * which integrate polynomial + * bilinear forms exactly. + */ + void initialize(const FiniteElement& el, + const Mapping& mapping, + const BlockInfo* block_info = 0); + + /** + * Initialize the + * IntegrationInfo objects + * contained. + * + * Before doing so, add update + * flags necessary to produce + * the data needed and also + * set uninitialized quadrature + * rules to Gauss formulas, + * which integrate polynomial + * bilinear forms exactly. + */ + template + void initialize(const FiniteElement& el, + const Mapping& mapping, + const NamedData& data, + const BlockInfo* block_info = 0); + + /** + * Initialize the + * IntegrationInfo objects + * contained. + * + * Before doing so, add update + * flags necessary to produce + * the data needed and also + * set uninitialized quadrature + * rules to Gauss formulas, + * which integrate polynomial + * bilinear forms exactly. + */ + template + void initialize(const FiniteElement& el, + const Mapping& mapping, + const NamedData*>& data, + const BlockInfo* block_info = 0); + /** + * @name FEValues setup + */ + /* @{ */ + void initialize_update_flags(); + + /** + * Add FEValues UpdateFlags for + * integration on all objects + * (cells, boundary faces and + * all interior faces). + */ + void add_update_flags_all (const UpdateFlags flags); + + /** + * Add FEValues UpdateFlags for + * integration on cells. + */ + void add_update_flags_cell(const UpdateFlags flags); + + /** + * Add FEValues UpdateFlags for + * integration on boundary faces. + */ + void add_update_flags_boundary(const UpdateFlags flags); + + /** + * Add FEValues UpdateFlags for + * integration on interior faces. + */ + void add_update_flags_face(const UpdateFlags flags); + + /** + * Add additional update flags + * to the ones already set in + * this program. The four + * boolean flags indicate + * wether the additional flags + * should be set for cell, + * boundary, interelement face + * for the cell itself + * or neighbor cell, or + * any combination thereof. + */ + void add_update_flags(const UpdateFlags flags, + const bool cell = true, + const bool boundary = true, + const bool face = true, + const bool neighbor = true); + + /** + * Assign n-point Gauss + * quadratures to each of the + * quadrature rules. Here, a + * size of zero points means + * that no loop over these grid + * entities should be + * performed. + * + * If the parameter + * force is true, then + * all quadrature sets are + * filled with new quadrature + * ruels. If it is false, then + * only empty rules are changed. + */ + void initialize_gauss_quadrature(unsigned int n_cell_points, + unsigned int n_boundary_points, + unsigned int n_face_points, + const bool force = true); + + /** + * The memory used by this object. + */ + std::size_t memory_consumption () const; + + /** + * The set of update flags + * for boundary cell integration. + * + * Defaults to + * #update_JxW_values. + */ + UpdateFlags cell_flags; + /** + * The set of update flags + * for boundary face integration. + * + * Defaults to + * #update_JxW_values and + * #update_normal_vectors. + */ + UpdateFlags boundary_flags; + + /** + * The set of update flags + * for interior face integration. + * + * Defaults to + * #update_JxW_values and + * #update_normal_vectors. + */ + UpdateFlags face_flags; + + /** + * The set of update flags + * for interior face integration. + * + * Defaults to + * #update_default, since + * quadrature weights are + * taken from the other cell. + */ + UpdateFlags neighbor_flags; + + /** + * The quadrature rule used + * on cells. + */ + Quadrature cell_quadrature; + + /** + * The quadrature rule used + * on boundary faces. + */ + Quadrature boundary_quadrature; + + /** + * The quadrature rule used + * on interior faces. + */ + Quadrature face_quadrature; + /* @} */ + + /** + * @name Data vectors + */ + /* @{ */ + + /** + * Initialize the + * VectorSelector objects + * #cell_selector, + * #boundary_selector and + * #face_selector in order to + * save computational + * eeffort. If no selectors + * are used, then values for + * all named vectors in + * DoFInfo::global_data will be + * computed in all quadrature + * points. + * + * This function will also + * add UpdateFlags to the + * flags stored in this class. + */ + /** + * Select the vectors from + * DoFInfo::global_data + * that should be computed in + * the quadrature points on cells. + */ + MeshWorker::VectorSelector cell_selector; + + /** + * Select the vectors from + * DoFInfo::global_data + * that should be computed in + * the quadrature points on + * boundary faces. + */ + MeshWorker::VectorSelector boundary_selector; + + /** + * Select the vectors from + * DoFInfo::global_data + * that should be computed in + * the quadrature points on + * interior faces. + */ + MeshWorker::VectorSelector face_selector; + + std_cxx1x::shared_ptr > cell_data; + std_cxx1x::shared_ptr > boundary_data; + std_cxx1x::shared_ptr > face_data; + /* @} */ + + /** + * @name Interface for MeshWorker::loop() + */ + /* @{ */ + /** + * A callback function which is + * called in the loop over all + * cells, after the action on a + * cell has been performed and + * before the faces are dealt + * with. + * + * In order for this function + * to have this effect, + * at least either of the + * arguments + * boundary_worker or + * face_worker + * arguments of loop() should + * be nonzero. Additionally, + * cells_first should + * be true. If + * cells_first is + * false, this function is + * called before any action on + * a cell is taken. + * + * And empty function in this + * class, but can be replaced + * in other classes given to + * loop() instead. + * + * See loop() and cell_action() + * for more details of how this + * function can be used. + */ + template + void post_cell(const DoFInfoBox&); + + /** + * A callback function which is + * called in the loop over all + * cells, after the action on + * the faces of a cell has been + * performed and before the + * cell itself is dealt with + * (assumes + * cells_first is false). + * + * In order for this function + * to have a reasonable effect, + * at least either of the + * arguments + * boundary_worker or + * face_worker + * arguments of loop() should + * be nonzero. Additionally, + * cells_first should + * be false. + * + * And empty function in this + * class, but can be replaced + * in other classes given to + * loop() instead. + * + * See loop() and cell_action() + * for more details of how this + * function can be used. + */ + template + void post_faces(const DoFInfoBox&); + + /** + * The info object for a cell. + */ + CellInfo cell; + /** + * The info object for a + * boundary face. + */ + CellInfo boundary; + /** + * The info object for a + * regular interior face, seen + * from the first cell. + */ + CellInfo face; + /** + * The info object for the + * refined side of an interior + * face seen from the first + * cell. + */ + CellInfo subface; + /** + * The info object for an + * interior face, seen from the + * other cell. + */ + CellInfo neighbor; + + /* @} */ + }; + + +//----------------------------------------------------------------------// + + template + inline + IntegrationInfo::IntegrationInfo() + : + fevalv(0), + multigrid(false), + global_data(std_cxx1x::shared_ptr >(new VectorDataBase)) + {} + + + template + inline + IntegrationInfo::IntegrationInfo(const IntegrationInfo& other) + : + multigrid(other.multigrid), + values(other.values), + gradients(other.gradients), + hessians(other.hessians), + global_data(other.global_data), + n_components(other.n_components) + { + fevalv.resize(other.fevalv.size()); + for (unsigned int i=0;i& p = *other.fevalv[i]; + const FEValues* pc = dynamic_cast*>(&p); + const FEFaceValues* pf = dynamic_cast*>(&p); + const FESubfaceValues* ps = dynamic_cast*>(&p); + + if (pc != 0) + fevalv[i] = std_cxx1x::shared_ptr > ( + reinterpret_cast*>( + new FEValues (pc->get_mapping(), pc->get_fe(), + pc->get_quadrature(), pc->get_update_flags()))); + else if (pf != 0) + fevalv[i] = std_cxx1x::shared_ptr > ( + new FEFaceValues (pf->get_mapping(), pf->get_fe(), pf->get_quadrature(), pf->get_update_flags())); + else if (ps != 0) + fevalv[i] = std_cxx1x::shared_ptr > ( + new FESubfaceValues (ps->get_mapping(), ps->get_fe(), ps->get_quadrature(), ps->get_update_flags())); + else + Assert(false, ExcInternalError()); + } + } + + + template <> + inline + IntegrationInfo<1,1>::IntegrationInfo(const IntegrationInfo<1,1>& other) + : + multigrid(other.multigrid), + values(other.values), + gradients(other.gradients), + hessians(other.hessians), + global_data(other.global_data), + n_components(other.n_components) + { + const int dim = 1; + const int sdim = 1; + + fevalv.resize(other.fevalv.size()); + for (unsigned int i=0;i& p = *other.fevalv[i]; + const FEValues* pc = dynamic_cast*>(&p); + const FEFaceValues* pf = dynamic_cast*>(&p); + const FESubfaceValues* ps = dynamic_cast*>(&p); + + if (pc != 0) + fevalv[i] = std_cxx1x::shared_ptr > ( + reinterpret_cast*>( + new FEValues (pc->get_mapping(), pc->get_fe(), + pc->get_quadrature(), pc->get_update_flags()))); + else if (pf != 0) + { + Assert (false, ExcImpossibleInDim(1)); + fevalv[i].reset (); + } + else if (ps != 0) + { + Assert (false, ExcImpossibleInDim(1)); + fevalv[i].reset(); + } + else + Assert(false, ExcInternalError()); + } + } + + + template <> + inline + IntegrationInfo<1,2>::IntegrationInfo(const IntegrationInfo<1,2>& other) + : + multigrid(other.multigrid), + values(other.values), + gradients(other.gradients), + hessians(other.hessians), + global_data(other.global_data), + n_components(other.n_components) + { + const int dim = 1; + const int sdim = 2; + + fevalv.resize(other.fevalv.size()); + for (unsigned int i=0;i& p = *other.fevalv[i]; + const FEValues* pc = dynamic_cast*>(&p); + const FEFaceValues* pf = dynamic_cast*>(&p); + const FESubfaceValues* ps = dynamic_cast*>(&p); + + if (pc != 0) + fevalv[i] = std_cxx1x::shared_ptr > ( + reinterpret_cast*>( + new FEValues (pc->get_mapping(), pc->get_fe(), + pc->get_quadrature(), pc->get_update_flags()))); + else if (pf != 0) + { + Assert (false, ExcImpossibleInDim(1)); + fevalv[i].reset(); + } + else if (ps != 0) + { + Assert (false, ExcImpossibleInDim(1)); + fevalv[i].reset(); + } + else + Assert(false, ExcInternalError()); + } + } + + + + template + template + inline void + IntegrationInfo::initialize( + const FiniteElement& el, + const Mapping& mapping, + const Quadrature& quadrature, + const UpdateFlags flags, + const BlockInfo* block_info) + { + if (block_info == 0 || block_info->local().size() == 0) + { + fevalv.resize(1); + fevalv[0] = std_cxx1x::shared_ptr > ( + new FEVALUES (mapping, el, quadrature, flags)); + } + else + { + fevalv.resize(el.n_base_elements()); + for (unsigned int i=0;i > ( + new FEVALUES (mapping, el.base_element(i), quadrature, flags)); + } + } + n_components = el.n_components(); + } + + + template + inline const FEValuesBase& + IntegrationInfo::fe_values() const + { + AssertDimension(fevalv.size(), 1); + return *fevalv[0]; + } + + + template + inline const FEValuesBase& + IntegrationInfo::fe_values(unsigned int i) const + { + Assert (i + inline void + IntegrationInfo::reinit(const DoFInfo& info) + { + for (unsigned int i=0;i& febase = *fevalv[i]; + if (info.sub_number != deal_II_numbers::invalid_unsigned_int) + { + // This is a subface + FESubfaceValues& fe = dynamic_cast&> (febase); + fe.reinit(info.cell, info.face_number, info.sub_number); + } + else if (info.face_number != deal_II_numbers::invalid_unsigned_int) + { + // This is a face + FEFaceValues& fe = dynamic_cast&> (febase); + fe.reinit(info.cell, info.face_number); + } + else + { + // This is a cell + FEValues& fe = dynamic_cast&> (febase); + fe.reinit(info.cell); + } + } + + const bool split_fevalues = info.block_info != 0; + if (!global_data->empty()) + fill_local_data(info, split_fevalues); + } + + + template <> + inline void + IntegrationInfo<1,1>::reinit(const DoFInfo<1,1>& info) + { + const int dim = 1; + const int spacedim = 1; + + for (unsigned int i=0;i& febase = *fevalv[i]; + if (info.sub_number != deal_II_numbers::invalid_unsigned_int) + { + // This is a subface + Assert (false, ExcImpossibleInDim(1)); + } + else if (info.face_number != deal_II_numbers::invalid_unsigned_int) + { + // This is a face + Assert (false, ExcImpossibleInDim(1)); + } + else + { + // This is a cell + FEValues& fe = dynamic_cast&> (febase); + fe.reinit(info.cell); + } + } + + const bool split_fevalues = info.block_info != 0; + if (!global_data->empty()) + fill_local_data(info, split_fevalues); + } + + + template <> + inline void + IntegrationInfo<1,2>::reinit(const DoFInfo<1,2>& info) + { + const int dim = 1; + const int spacedim = 2; + + for (unsigned int i=0;i& febase = *fevalv[i]; + if (info.sub_number != deal_II_numbers::invalid_unsigned_int) + { + // This is a subface + Assert (false, ExcImpossibleInDim(1)); + } + else if (info.face_number != deal_II_numbers::invalid_unsigned_int) + { + // This is a face + Assert (false, ExcImpossibleInDim(1)); + } + else + { + // This is a cell + FEValues& fe = dynamic_cast&> (febase); + fe.reinit(info.cell); + } + } + + const bool split_fevalues = info.block_info != 0; + if (!global_data->empty()) + fill_local_data(info, split_fevalues); + } + + +//----------------------------------------------------------------------// + + template <> + inline + void + IntegrationInfoBox<1,1>::initialize_gauss_quadrature( + const unsigned int cp, + const unsigned int, + const unsigned int, + bool force) + { + if (force || cell_quadrature.size() == 0) + cell_quadrature = QGauss<1>(cp); + } + + + template <> + inline + void + IntegrationInfoBox<1,2>::initialize_gauss_quadrature( + const unsigned int cp, + const unsigned int, + const unsigned int, + bool force) + { + if (force || cell_quadrature.size() == 0) + cell_quadrature = QGauss<1>(cp); + } + + + template + inline + void + IntegrationInfoBox::initialize_gauss_quadrature( + unsigned int cp, + unsigned int bp, + unsigned int fp, + bool force) + { + if (force || cell_quadrature.size() == 0) + cell_quadrature = QGauss(cp); + if (force || boundary_quadrature.size() == 0) + boundary_quadrature = QGauss(bp); + if (force || face_quadrature.size() == 0) + face_quadrature = QGauss(fp); + } + + + template + inline + void + IntegrationInfoBox::add_update_flags_all (const UpdateFlags flags) + { + add_update_flags(flags, true, true, true, true); + } + + + template + inline + void + IntegrationInfoBox::add_update_flags_cell (const UpdateFlags flags) + { + add_update_flags(flags, true, false, false, false); + } + + + template + inline + void + IntegrationInfoBox::add_update_flags_boundary (const UpdateFlags flags) + { + add_update_flags(flags, false, true, false, false); + } + + + template + inline + void + IntegrationInfoBox::add_update_flags_face (const UpdateFlags flags) + { + add_update_flags(flags, false, false, true, true); + } + + + template <> + inline + void + IntegrationInfoBox<1,1>::initialize( + const FiniteElement<1,1>& el, + const Mapping<1,1>& mapping, + const BlockInfo* block_info) + { + initialize_update_flags(); + initialize_gauss_quadrature( + (cell_flags & update_values) ? (el.tensor_degree()+1) : el.tensor_degree(), 1, 1, false); + + const int dim = 1; + const int sdim = 1; + + cell.initialize >(el, mapping, cell_quadrature, + cell_flags, block_info); + } + + + + template <> + inline + void + IntegrationInfoBox<1,2>::initialize( + const FiniteElement<1,2>& el, + const Mapping<1,2>& mapping, + const BlockInfo* block_info) + { + initialize_update_flags(); + initialize_gauss_quadrature( + (cell_flags & update_values) ? (el.tensor_degree()+1) : el.tensor_degree(), 1, 1, false); + + const int dim = 1; + const int sdim = 2; + + cell.initialize >(el, mapping, cell_quadrature, + cell_flags, block_info); + } + + + template + inline + void + IntegrationInfoBox::initialize( + const FiniteElement& el, + const Mapping& mapping, + const BlockInfo* block_info) + { + initialize_update_flags(); + initialize_gauss_quadrature( + (cell_flags & update_values) ? (el.tensor_degree()+1) : el.tensor_degree(), + (boundary_flags & update_values) ? (el.tensor_degree()+1) : el.tensor_degree(), + (face_flags & update_values) ? (el.tensor_degree()+1) : el.tensor_degree(), false); + + cell.template initialize >(el, mapping, cell_quadrature, + cell_flags, block_info); + boundary.template initialize >(el, mapping, boundary_quadrature, + boundary_flags, block_info); + face.template initialize >(el, mapping, face_quadrature, + face_flags, block_info); + subface.template initialize >(el, mapping, face_quadrature, + face_flags, block_info); + neighbor.template initialize >(el, mapping, face_quadrature, + neighbor_flags, block_info); + } + + + template + template + void + IntegrationInfoBox::initialize( + const FiniteElement& el, + const Mapping& mapping, + const NamedData& data, + const BlockInfo* block_info) + { + initialize(el, mapping, block_info); + std_cxx1x::shared_ptr > p; + + p = std_cxx1x::shared_ptr >(new VectorData (cell_selector)); + p->initialize(data); + cell_data = p; + cell.initialize_data(p); + + p = std_cxx1x::shared_ptr >(new VectorData (boundary_selector)); + p->initialize(data); + boundary_data = p; + boundary.initialize_data(p); + + p = std_cxx1x::shared_ptr >(new VectorData (face_selector)); + p->initialize(data); + face_data = p; + face.initialize_data(p); + subface.initialize_data(p); + neighbor.initialize_data(p); + } + + + template + template + void + IntegrationInfoBox::initialize( + const FiniteElement& el, + const Mapping& mapping, + const NamedData*>& data, + const BlockInfo* block_info) + { + initialize(el, mapping, block_info); + std_cxx1x::shared_ptr > p; + + p = std_cxx1x::shared_ptr >(new MGVectorData (cell_selector)); + p->initialize(data); + cell_data = p; + cell.initialize_data(p); + + p = std_cxx1x::shared_ptr >(new MGVectorData (boundary_selector)); + p->initialize(data); + boundary_data = p; + boundary.initialize_data(p); + + p = std_cxx1x::shared_ptr >(new MGVectorData (face_selector)); + p->initialize(data); + face_data = p; + face.initialize_data(p); + subface.initialize_data(p); + neighbor.initialize_data(p); + } + + + template + template + void + IntegrationInfoBox::post_cell(const DoFInfoBox&) + {} + + + template + template + void + IntegrationInfoBox::post_faces(const DoFInfoBox&) + {} + + +} + +DEAL_II_NAMESPACE_CLOSE + +#endif diff --git a/deal.II/include/deal.II/numerics/mesh_worker_info.h b/deal.II/include/deal.II/numerics/mesh_worker_info.h new file mode 100644 index 0000000000..082afe94b3 --- /dev/null +++ b/deal.II/include/deal.II/numerics/mesh_worker_info.h @@ -0,0 +1,20 @@ +//--------------------------------------------------------------------------- +// $Id$ +// +// Copyright (C) 2010, 2011 by the deal.II authors +// +// 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. +// +//--------------------------------------------------------------------------- + +/** + * @file + * @deprecated The contents of this file has been split into + * meshworker/dof_info.h and meshworker/integration_info.h + */ + +#include +#include -- 2.39.5