From: Luca Heltai Date: Sat, 16 Aug 2014 16:47:09 +0000 (+0200) Subject: Modified all get_name functions to include the spacedim parameter. X-Git-Tag: v8.2.0-rc1~191^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=0ce2bc450408e46458e1d4cb5a954a1814a54e7c;p=dealii.git Modified all get_name functions to include the spacedim parameter. Added one test for codimension one finite element generation from names. Reintroduce get_fe_from_name() with one template parameter as deprecated New function get_fe_by_name() --- diff --git a/include/deal.II/base/utilities.h b/include/deal.II/base/utilities.h index 5d11071adb..6f926359a8 100644 --- a/include/deal.II/base/utilities.h +++ b/include/deal.II/base/utilities.h @@ -81,6 +81,19 @@ namespace Utilities int string_to_int (const std::string &s); + /** + * Return a string describing the dimensions of the object. Often, + * functions in the deal.II library as well as in user codes need to + * define a string containing the template dimensions of some + * objects defined using two template parameters: dim (the + * topological dimension of the object) and spacedim (the dimension + * of the embedding Euclidean space). Since in all deal.II classes, + * by default spacedim is equal to dimension, the above string is + * usually contracted to , instead of . This + * function returns a string containing "dim" if dim is equal to + * spacedim, otherwhise it returns "dim,spacedim". + */ + std::string dim_string(const int dim, const int spacedim); /** * Given a list of strings, convert it to a diff --git a/include/deal.II/fe/fe.h b/include/deal.II/fe/fe.h index eca8f11029..c47bcb9b90 100644 --- a/include/deal.II/fe/fe.h +++ b/include/deal.II/fe/fe.h @@ -355,10 +355,10 @@ class FiniteElement : public Subscriptor, { public: /** - * The dimension of the image space, corresponding to Triangulation. + * The dimension of the image space, corresponding to Triangulation. */ static const unsigned int space_dimension = spacedim; - + /** * Base class for internal data. Adds data for second derivatives to * Mapping::InternalDataBase() @@ -1007,7 +1007,7 @@ public: const bool face_orientation = true, const bool face_flip = false, const bool face_rotation = false) const; - + /** * For lines with non-standard line_orientation in 3D, the dofs on lines * have to be permuted in order to be combined with the correct shape diff --git a/include/deal.II/fe/fe_tools.h b/include/deal.II/fe/fe_tools.h index 0435498d3b..53724b70de 100644 --- a/include/deal.II/fe/fe_tools.h +++ b/include/deal.II/fe/fe_tools.h @@ -20,6 +20,7 @@ #include +#include #include #include #include @@ -76,10 +77,6 @@ namespace FETools * This class is used in the FETools::get_fe_from_name() and * FETools::add_fe_name() functions. * - * The Subscriptor base class was introduced such that we can have - * common pointers to the base class and then get the actual object - * through a `dynamic_cast`. - * * @author Guido Kanschat, 2006 */ template @@ -814,7 +811,7 @@ namespace FETools template void hierarchic_to_lexicographic_numbering (unsigned int degree, - std::vector &h2l); + std::vector &h2l); template void @@ -879,8 +876,16 @@ namespace FETools * function, use the add_fe_name() function. This function does not work if * one wants to get a codimension 1 finite element. */ - template + template FiniteElement * + get_fe_by_name (const std::string &name); + + + /** + * @deprecated Use get_fe_by_name() with two template parameters instead + */ + template + FiniteElement * get_fe_from_name (const std::string &name); @@ -926,7 +931,7 @@ namespace FETools * then you should call this function for each space dimension for which you * want your finite element added to the map. */ - template + template void add_fe_name (const std::string &name, const FEFactoryBase *factory); diff --git a/source/base/utilities.cc b/source/base/utilities.cc index e7934eda30..6823771b6b 100644 --- a/source/base/utilities.cc +++ b/source/base/utilities.cc @@ -70,7 +70,6 @@ namespace Utilities << "Can't convert the string " << arg1 << " to the desired type"); - std::string int_to_string (const unsigned int i, const unsigned int digits) @@ -124,6 +123,15 @@ namespace Utilities } + std::string + dim_string(const int dim, const int spacedim) + { + if (dim == spacedim) + return int_to_string(dim); + else + return int_to_string(dim)+","+int_to_string(spacedim); + } + unsigned int needed_digits (const unsigned int max_number) diff --git a/source/fe/fe_dg_vector.cc b/source/fe/fe_dg_vector.cc index a22345159f..7fce049e94 100644 --- a/source/fe/fe_dg_vector.cc +++ b/source/fe/fe_dg_vector.cc @@ -41,7 +41,9 @@ FE_DGNedelec::get_name () const // have to be kept in synch std::ostringstream namebuf; - namebuf << "FE_DGNedelec<" << dim << ',' << spacedim << ">(" << this->degree-1 << ")"; + namebuf << "FE_DGNedelec<" + << Utilities::dim_string(dim,spacedim) + << ">(" << this->degree-1 << ")"; return namebuf.str(); } @@ -65,7 +67,9 @@ FE_DGRaviartThomas::get_name () const // have to be kept in synch std::ostringstream namebuf; - namebuf << "FE_DGRaviartThomas<" << dim << ',' << spacedim << ">(" << this->degree-1 << ")"; + namebuf << "FE_DGRaviartThomas<" + << Utilities::dim_string(dim,spacedim) + << ">(" << this->degree-1 << ")"; return namebuf.str(); } @@ -89,7 +93,9 @@ FE_DGBDM::get_name () const // have to be kept in synch std::ostringstream namebuf; - namebuf << "FE_DGBDM<" << dim << ',' << spacedim << ">(" << this->degree-1 << ")"; + namebuf << "FE_DGBDM<" + << Utilities::dim_string(dim,spacedim) + << ">(" << this->degree-1 << ")"; return namebuf.str(); } diff --git a/source/fe/fe_dgp.cc b/source/fe/fe_dgp.cc index daa2b51d3a..61c67289d5 100644 --- a/source/fe/fe_dgp.cc +++ b/source/fe/fe_dgp.cc @@ -54,7 +54,9 @@ FE_DGP::get_name () const // kept in synch std::ostringstream namebuf; - namebuf << "FE_DGP<" << dim << ">(" << this->degree << ")"; + namebuf << "FE_DGP<" + << Utilities::dim_string(dim,spacedim) + << ">(" << this->degree << ")"; return namebuf.str(); } @@ -244,7 +246,7 @@ FE_DGP::get_constant_modes () const Table<2,bool> constant_modes(1, this->dofs_per_cell); constant_modes(0,0) = true; return std::pair, std::vector > - (constant_modes, std::vector(1, 0)); + (constant_modes, std::vector(1, 0)); } diff --git a/source/fe/fe_dgp_nonparametric.cc b/source/fe/fe_dgp_nonparametric.cc index abecac8462..448f5525e3 100644 --- a/source/fe/fe_dgp_nonparametric.cc +++ b/source/fe/fe_dgp_nonparametric.cc @@ -106,7 +106,9 @@ FE_DGPNonparametric::get_name () const // have to be kept in synch std::ostringstream namebuf; - namebuf << "FE_DGPNonparametric<" << dim << ">(" << degree << ")"; + namebuf << "FE_DGPNonparametric<" + << Utilities::dim_string(dim,spacedim) + << ">(" << degree << ")"; return namebuf.str(); } diff --git a/source/fe/fe_dgq.cc b/source/fe/fe_dgq.cc index 0ca4d1c021..b71127a9f6 100644 --- a/source/fe/fe_dgq.cc +++ b/source/fe/fe_dgq.cc @@ -224,7 +224,9 @@ FE_DGQ::get_name () const // have to be kept in synch std::ostringstream namebuf; - namebuf << "FE_DGQ<" << dim << ">(" << this->degree << ")"; + namebuf << "FE_DGQ<" + << Utilities::dim_string(dim,spacedim) + << ">(" << this->degree << ")"; return namebuf.str(); } diff --git a/source/fe/fe_face.cc b/source/fe/fe_face.cc index 3a210bd849..54f2465536 100644 --- a/source/fe/fe_face.cc +++ b/source/fe/fe_face.cc @@ -119,7 +119,9 @@ FE_FaceQ::get_name () const // particular format of the string this function returns, so they have to be // kept in synch std::ostringstream namebuf; - namebuf << "FE_FaceQ<" << dim << ">(" << this->degree << ")"; + namebuf << "FE_FaceQ<" + << Utilities::dim_string(dim,spacedim) + << ">(" << this->degree << ")"; return namebuf.str(); } @@ -296,7 +298,7 @@ FE_FaceQ::get_constant_modes () const for (unsigned int i=0; idofs_per_cell; ++i) constant_modes(0,i) = true; return std::pair, std::vector > - (constant_modes, std::vector(1, 0)); + (constant_modes, std::vector(1, 0)); } @@ -337,7 +339,9 @@ FE_FaceQ<1,spacedim>::get_name () const // particular format of the string this function returns, so they have to be // kept in synch std::ostringstream namebuf; - namebuf << "FE_FaceQ<1>(" << this->degree << ")"; + namebuf << "FE_FaceQ<" + << Utilities::dim_string(1,spacedim) + << ">(" << this->degree << ")"; return namebuf.str(); } @@ -424,7 +428,7 @@ FE_FaceQ<1,spacedim>::get_constant_modes () const for (unsigned int i=0; idofs_per_cell; ++i) constant_modes(0,i) = true; return std::pair, std::vector > - (constant_modes, std::vector(1,0)); + (constant_modes, std::vector(1,0)); } @@ -471,7 +475,7 @@ typename Mapping<1,spacedim>::InternalDataBase * FE_FaceQ<1,spacedim>::get_face_data ( const UpdateFlags update_flags, const Mapping<1,spacedim> &, - const Quadrature<0>& quadrature) const + const Quadrature<0> &quadrature) const { // generate a new data object and initialize some fields typename Mapping<1,spacedim>::InternalDataBase *data = @@ -503,7 +507,7 @@ typename Mapping<1,spacedim>::InternalDataBase * FE_FaceQ<1,spacedim>::get_subface_data ( const UpdateFlags flags, const Mapping<1,spacedim> &mapping, - const Quadrature<0>& quadrature) const + const Quadrature<0> &quadrature) const { return get_face_data (flags, mapping, quadrature); } @@ -532,7 +536,7 @@ FE_FaceQ<1,spacedim>::fill_fe_face_values ( const Mapping<1,spacedim> &, const typename Triangulation<1,spacedim>::cell_iterator &, const unsigned int face, - const Quadrature<0>& quadrature, + const Quadrature<0> &quadrature, typename Mapping<1,spacedim>::InternalDataBase &, typename Mapping<1,spacedim>::InternalDataBase &fedata, FEValuesData<1,spacedim> &data) const @@ -556,7 +560,7 @@ FE_FaceQ<1,spacedim>::fill_fe_subface_values ( const typename Triangulation<1,spacedim>::cell_iterator &, const unsigned int , const unsigned int , - const Quadrature<0>& , + const Quadrature<0> & , typename Mapping<1,spacedim>::InternalDataBase &, typename Mapping<1,spacedim>::InternalDataBase &, FEValuesData<1,spacedim> &) const @@ -596,7 +600,9 @@ FE_FaceP::get_name () const // particular format of the string this function returns, so they have to be // kept in synch std::ostringstream namebuf; - namebuf << "FE_FaceP<" << dim << ">(" << this->degree << ")"; + namebuf << "FE_FaceP<" + << Utilities::dim_string(dim,spacedim) + << ">(" << this->degree << ")"; return namebuf.str(); } @@ -789,7 +795,7 @@ FE_FaceP::get_constant_modes () const for (unsigned int face=0; face::faces_per_cell; ++face) constant_modes(0, face*this->dofs_per_face) = true; return std::pair, std::vector > - (constant_modes, std::vector(1, 0)); + (constant_modes, std::vector(1, 0)); } @@ -810,7 +816,9 @@ FE_FaceP<1,spacedim>::get_name () const // particular format of the string this function returns, so they have to be // kept in synch std::ostringstream namebuf; - namebuf << "FE_FaceP<1>(" << this->degree << ")"; + namebuf << "FE_FaceP<" + << Utilities::dim_string(1,spacedim) + << ">(" << this->degree << ")"; return namebuf.str(); } diff --git a/source/fe/fe_q.cc b/source/fe/fe_q.cc index 717d783d82..4239d40cc9 100644 --- a/source/fe/fe_q.cc +++ b/source/fe/fe_q.cc @@ -93,7 +93,9 @@ FE_Q::get_name () const } if (equidistant == true) - namebuf << "FE_Q<" << dim << ">(" << this->degree << ")"; + namebuf << "FE_Q<" + << Utilities::dim_string(dim,spacedim) + << ">(" << this->degree << ")"; else { // Check whether the support points come from QGaussLobatto. @@ -106,9 +108,13 @@ FE_Q::get_name () const break; } if (gauss_lobatto == true) - namebuf << "FE_Q<" << dim << ">(QGaussLobatto(" << this->degree+1 << "))"; + namebuf << "FE_Q<" + << Utilities::dim_string(dim,spacedim) + << ">(QGaussLobatto(" << this->degree+1 << "))"; else - namebuf << "FE_Q<" << dim << ">(QUnknownNodes(" << this->degree << "))"; + namebuf << "FE_Q<" + << Utilities::dim_string(dim,spacedim) + << ">(QUnknownNodes(" << this->degree << "))"; } return namebuf.str(); } diff --git a/source/fe/fe_q_dg0.cc b/source/fe/fe_q_dg0.cc index 99ba0e946b..fbf453a6f2 100644 --- a/source/fe/fe_q_dg0.cc +++ b/source/fe/fe_q_dg0.cc @@ -133,7 +133,9 @@ FE_Q_DG0::get_name () const } if (type == true) - namebuf << "FE_Q_DG0<" << dim << ">(" << this->degree << ")"; + namebuf << "FE_Q_DG0<" + << Utilities::dim_string(dim,spacedim) + << ">(" << this->degree << ")"; else { @@ -147,9 +149,13 @@ FE_Q_DG0::get_name () const break; } if (type == true) - namebuf << "FE_Q_DG0<" << dim << ">(QGaussLobatto(" << this->degree+1 << "))"; + namebuf << "FE_Q_DG0<" + << Utilities::dim_string(dim,spacedim) + << ">(QGaussLobatto(" << this->degree+1 << "))"; else - namebuf << "FE_Q_DG0<" << dim << ">(QUnknownNodes(" << this->degree << "))"; + namebuf << "FE_Q_DG0<" + << Utilities::dim_string(dim,spacedim) + << ">(QUnknownNodes(" << this->degree << "))"; } return namebuf.str(); } diff --git a/source/fe/fe_q_iso_q1.cc b/source/fe/fe_q_iso_q1.cc index 6aba97cfac..11f1081c7a 100644 --- a/source/fe/fe_q_iso_q1.cc +++ b/source/fe/fe_q_iso_q1.cc @@ -60,7 +60,9 @@ FE_Q_iso_Q1::get_name () const // kept in synch std::ostringstream namebuf; - namebuf << "FE_Q_iso_Q1<" << dim << ">(" << this->degree << ")"; + namebuf << "FE_Q_iso_Q1<" + << Utilities::dim_string(dim,spacedim) + << ">(" << this->degree << ")"; return namebuf.str(); } diff --git a/source/fe/fe_system.cc b/source/fe/fe_system.cc index 78bce6132b..0d1acf413b 100644 --- a/source/fe/fe_system.cc +++ b/source/fe/fe_system.cc @@ -347,7 +347,9 @@ FESystem::get_name () const std::ostringstream namebuf; - namebuf << "FESystem<" << dim << ">["; + namebuf << "FESystem<" + << Utilities::dim_string(dim,spacedim) + << ">["; for (unsigned int i=0; i< this->n_base_elements(); ++i) { namebuf << base_element(i).get_name(); diff --git a/source/fe/fe_tools.cc b/source/fe/fe_tools.cc index 8c3e00c2d5..b9edc22d63 100644 --- a/source/fe/fe_tools.cc +++ b/source/fe/fe_tools.cc @@ -108,17 +108,12 @@ namespace FETools namespace { - // use a shared pointer to factory - // objects, to ensure that they get - // deleted at the end of the - // program run and don't end up as - // apparent memory leaks to - // programs like valgrind. - // a function that returns the - // default set of finite element - // names and factory objects for - // them. used to initialize - // fe_name_map below + // The following three functions serve to fill the maps from element + // names to elements fe_name_map below. The first one exists because + // we have finite elements which are not implemented for nonzero + // codimension. These should be transfered to the second function + // eventually. + template void fill_no_codim_fe_names (std::map >& result) @@ -151,6 +146,9 @@ namespace = FEFactoryPointer(new FETools::FEFactory >); } + // This function fills a map from names to finite elements for any + // dimension and codimension for those elements which support + // nonzero codimension. template void fill_codim_fe_names (std::map >& result) @@ -167,6 +165,10 @@ namespace = FEFactoryPointer(new FETools::FEFactory >); } + // The function filling the vector fe_name_map below. It iterates + // through all legal dimension/spacedimension pairs and fills + // fe_name_map[dimension][spacedimension] with the maps generated + // by the functions above. std::vector > > > @@ -192,44 +194,37 @@ fill_default_map() } - // have a lock that guarantees that - // at most one thread is changing - // and accessing the fe_name_map - // variable. make this lock local - // to this file. + // have a lock that guarantees that at most one thread is changing + // and accessing the fe_name_map variable. make this lock local to + // this file. // - // this and the next variable are - // declared static (even though - // they're in an anonymous - // namespace) in order to make icc - // happy (which otherwise reports a - // multiply defined symbol when - // linking libraries for more than - // one space dimension together + // this and the next variable are declared static (even though + // they're in an anonymous namespace) in order to make icc happy + // (which otherwise reports a multiply defined symbol when linking + // libraries for more than one space dimension together static Threads::Mutex fe_name_map_lock; - // This is the map used by - // FETools::get_fe_from_name and - // FETools::add_fe_name. Since - // FEFactoryBase has a template - // parameter dim, it could not be a - // member variable of FETools. On - // the other hand, it is only - // accessed by functions in this - // file, so it is safe to make it a - // static variable here. It must be - // static so that we can link - // several dimensions together. - // - // it is initialized at program start time - // using the function above. because at - // this time there are no threads running, - // there are no thread-safety issues - // here. since this is compiled for all - // dimensions at once, need to create - // objects for each dimension and then - // separate between them further down + // This is the map used by FETools::get_fe_from_name and + // FETools::add_fe_name. It is only accessed by functions in this + // file, so it is safe to make it a static variable here. It must be + // static so that we can link several dimensions together. + + // The organization of this storage is such that + // fe_name_map[dim][spacedim][name] points to an + // FEFactoryBase with the name given. Since + // all entries of this vector are of different type, we store + // pointers to generic objects and cast them when needed. + + // We use a shared pointer to factory objects, to ensure that they + // get deleted at the end of the program run and don't end up as + // apparent memory leaks to programs like valgrind. + + // This vector is initialized at program start time using the + // function above. because at this time there are no threads + // running, there are no thread-safety issues here. since this is + // compiled for all dimensions at once, need to create objects for + // each dimension and then separate between them further down static std::vector FiniteElement * - get_fe_from_name (const std::string ¶meter_name) + get_fe_by_name (const std::string ¶meter_name) { // Create a version of the name // string where all template @@ -1634,13 +1629,13 @@ namespace FETools } -// template -// FiniteElement * -// get_fe_from_name (const std::string ¶meter_name) -// { -// return internal::get_fe_from_name(parameter_name); -// } - + template + FiniteElement * + get_fe_from_name (const std::string ¶meter_name) + { + return get_fe_by_name (parameter_name); + } + template void diff --git a/source/fe/fe_tools.inst.in b/source/fe/fe_tools.inst.in index fc0ae9bbf2..7954113c65 100644 --- a/source/fe/fe_tools.inst.in +++ b/source/fe/fe_tools.inst.in @@ -21,19 +21,19 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS \{ #if deal_II_dimension <= deal_II_space_dimension template - void compute_block_renumbering ( - const FiniteElement& , - std::vector&, std::vector&, bool); + void compute_block_renumbering ( + const FiniteElement & , + std::vector &, std::vector &, bool); template - void compute_projection_matrices - (const FiniteElement &, - std::vector > >&, bool); + void compute_projection_matrices + (const FiniteElement &, + std::vector > > &, bool); template - void compute_embedding_matrices - (const FiniteElement &, - std::vector > >&,bool); + void compute_embedding_matrices + (const FiniteElement &, + std::vector > > &,bool); #endif \} } @@ -47,16 +47,16 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS { namespace FETools \{ - + #if deal_II_dimension <= deal_II_space_dimension template - void get_interpolation_matrix - (const FiniteElement &, - const FiniteElement &, - FullMatrix &); + void get_interpolation_matrix + (const FiniteElement &, + const FiniteElement &, + FullMatrix &); template FiniteElement * - get_fe_from_name (const std::string &); + get_fe_by_name (const std::string &); #endif @@ -64,120 +64,123 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS template class FEFactoryBase; + template FiniteElement * + get_fe_from_name (const std::string &); + template - void compute_node_matrix( - FullMatrix&, - const FiniteElement&); + void compute_node_matrix( + FullMatrix &, + const FiniteElement &); template - void compute_component_wise( - const FiniteElement& , - std::vector&, std::vector >&); + void compute_component_wise( + const FiniteElement & , + std::vector &, std::vector > &); template - void get_back_interpolation_matrix - (const FiniteElement &, - const FiniteElement &, - FullMatrix &); + void get_back_interpolation_matrix + (const FiniteElement &, + const FiniteElement &, + FullMatrix &); template - void get_interpolation_difference_matrix - (const FiniteElement &, - const FiniteElement &, - FullMatrix &); + void get_interpolation_difference_matrix + (const FiniteElement &, + const FiniteElement &, + FullMatrix &); template - void get_interpolation_matrix - (const FiniteElement &, - const FiniteElement &, - FullMatrix &); + void get_interpolation_matrix + (const FiniteElement &, + const FiniteElement &, + FullMatrix &); template - void get_back_interpolation_matrix - (const FiniteElement &, - const FiniteElement &, - FullMatrix &); + void get_back_interpolation_matrix + (const FiniteElement &, + const FiniteElement &, + FullMatrix &); template - void get_interpolation_difference_matrix - (const FiniteElement &, - const FiniteElement &, - FullMatrix &); + void get_interpolation_difference_matrix + (const FiniteElement &, + const FiniteElement &, + FullMatrix &); template - void get_projection_matrix - (const FiniteElement &, - const FiniteElement &, - FullMatrix &); + void get_projection_matrix + (const FiniteElement &, + const FiniteElement &, + FullMatrix &); template - void compute_face_embedding_matrices - (const FiniteElement &, FullMatrix (&)[GeometryInfo::max_children_per_face], - unsigned int, unsigned int); + void compute_face_embedding_matrices + (const FiniteElement &, FullMatrix ( &)[GeometryInfo::max_children_per_face], + unsigned int, unsigned int); template - void - compute_projection_from_quadrature_points_matrix (const FiniteElement &, - const Quadrature &, - const Quadrature &, - FullMatrix &); + void + compute_projection_from_quadrature_points_matrix (const FiniteElement &, + const Quadrature &, + const Quadrature &, + FullMatrix &); template - void - compute_projection_from_quadrature_points( - const FullMatrix &, - const std::vector< Tensor<1, deal_II_dimension > > &, - std::vector< Tensor<1, deal_II_dimension > > &); + void + compute_projection_from_quadrature_points( + const FullMatrix &, + const std::vector< Tensor<1, deal_II_dimension > > &, + std::vector< Tensor<1, deal_II_dimension > > &); template - void - compute_projection_from_quadrature_points( - const FullMatrix &, - const std::vector > &, - std::vector > &); + void + compute_projection_from_quadrature_points( + const FullMatrix &, + const std::vector > &, + std::vector > &); template - void - compute_interpolation_to_quadrature_points_matrix (const FiniteElement &, - const Quadrature &, - FullMatrix &); + void + compute_interpolation_to_quadrature_points_matrix (const FiniteElement &, + const Quadrature &, + FullMatrix &); #if deal_II_dimension != 1 template - void - compute_projection_from_face_quadrature_points_matrix (const FiniteElement &, - const Quadrature &, - const Quadrature &, - const DoFHandler::active_cell_iterator & , - unsigned int, - FullMatrix &); + void + compute_projection_from_face_quadrature_points_matrix (const FiniteElement &, + const Quadrature &, + const Quadrature &, + const DoFHandler::active_cell_iterator & , + unsigned int, + FullMatrix &); #endif template - void - hierarchic_to_lexicographic_numbering - (unsigned int, - std::vector &); + void + hierarchic_to_lexicographic_numbering + (unsigned int, + std::vector &); template - void - hierarchic_to_lexicographic_numbering - (const FiniteElementData &, - std::vector &); + void + hierarchic_to_lexicographic_numbering + (const FiniteElementData &, + std::vector &); template - void - lexicographic_to_hierarchic_numbering - (const FiniteElementData &, - std::vector &); + void + lexicographic_to_hierarchic_numbering + (const FiniteElementData &, + std::vector &); template - std::vector - hierarchic_to_lexicographic_numbering - (const FiniteElementData &); + std::vector + hierarchic_to_lexicographic_numbering + (const FiniteElementData &); template - std::vector - lexicographic_to_hierarchic_numbering - (const FiniteElementData &); - + std::vector + lexicographic_to_hierarchic_numbering + (const FiniteElementData &); + #endif \} } diff --git a/source/fe/fe_trace.cc b/source/fe/fe_trace.cc index 4906f0d54e..2584e074f5 100644 --- a/source/fe/fe_trace.cc +++ b/source/fe/fe_trace.cc @@ -32,11 +32,11 @@ DEAL_II_NAMESPACE_OPEN template FE_TraceQ::FE_TraceQ (const unsigned int degree) - : - FE_PolyFace, dim, spacedim> ( - TensorProductPolynomials(Polynomials::LagrangeEquidistant::generate_complete_basis(degree)), - FiniteElementData(get_dpo_vector(degree), 1, degree, FiniteElementData::L2), - std::vector(1,true)) + : + FE_PolyFace, dim, spacedim> ( + TensorProductPolynomials(Polynomials::LagrangeEquidistant::generate_complete_basis(degree)), + FiniteElementData(get_dpo_vector(degree), 1, degree, FiniteElementData::L2), + std::vector(1,true)) { Assert (degree > 0, ExcMessage ("FE_Trace can only be used for polynomial degrees " @@ -52,7 +52,7 @@ FE_TraceQ::FE_TraceQ (const unsigned int degree) template -FiniteElement* +FiniteElement * FE_TraceQ::clone() const { return new FE_TraceQ(this->degree); @@ -63,15 +63,17 @@ template std::string FE_TraceQ::get_name () const { - // note that the - // FETools::get_fe_from_name - // function depends on the - // particular format of the string - // this function returns, so they - // have to be kept in synch + // note that the + // FETools::get_fe_from_name + // function depends on the + // particular format of the string + // this function returns, so they + // have to be kept in synch std::ostringstream namebuf; - namebuf << "FE_TraceQ<" << dim << ">(" << this->degree << ")"; + namebuf << "FE_TraceQ<" + << Utilities::dim_string(dim,spacedim) + << ">(" << this->degree << ")"; return namebuf.str(); } @@ -86,17 +88,17 @@ FE_TraceQ::has_support_on_face (const unsigned int shape_index, Assert (face_index < GeometryInfo::faces_per_cell, ExcIndexRange (face_index, 0, GeometryInfo::faces_per_cell)); - // let's see whether this is a - // vertex + // let's see whether this is a + // vertex if (shape_index < this->first_line_index) { - // for Q elements, there is - // one dof per vertex, so - // shape_index==vertex_number. check - // whether this vertex is - // on the given face. thus, - // for each face, give a - // list of vertices + // for Q elements, there is + // one dof per vertex, so + // shape_index==vertex_number. check + // whether this vertex is + // on the given face. thus, + // for each face, give a + // list of vertices const unsigned int vertex_no = shape_index; Assert (vertex_no < GeometryInfo::vertices_per_cell, ExcInternalError()); @@ -108,23 +110,23 @@ FE_TraceQ::has_support_on_face (const unsigned int shape_index, return false; } else if (shape_index < this->first_quad_index) - // ok, dof is on a line + // ok, dof is on a line { const unsigned int line_index = (shape_index - this->first_line_index) / this->dofs_per_line; Assert (line_index < GeometryInfo::lines_per_cell, ExcInternalError()); - // in 2d, the line is the - // face, so get the line - // index + // in 2d, the line is the + // face, so get the line + // index if (dim == 2) return (line_index == face_index); else if (dim == 3) { - // see whether the - // given line is on the - // given face. + // see whether the + // given line is on the + // given face. for (unsigned int l=0; l::lines_per_face; ++l) if (GeometryInfo<3>::face_to_cell_lines(face_index, l) == line_index) return true; @@ -135,7 +137,7 @@ FE_TraceQ::has_support_on_face (const unsigned int shape_index, Assert (false, ExcNotImplemented()); } else if (shape_index < this->first_hex_index) - // dof is on a quad + // dof is on a quad { const unsigned int quad_index = (shape_index - this->first_quad_index) / this->dofs_per_quad; @@ -143,31 +145,31 @@ FE_TraceQ::has_support_on_face (const unsigned int shape_index, static_cast(GeometryInfo::quads_per_cell), ExcInternalError()); - // in 2d, cell bubble are - // zero on all faces. but - // we have treated this - // case above already + // in 2d, cell bubble are + // zero on all faces. but + // we have treated this + // case above already Assert (dim != 2, ExcInternalError()); - // in 3d, - // quad_index=face_index + // in 3d, + // quad_index=face_index if (dim == 3) return (quad_index == face_index); else Assert (false, ExcNotImplemented()); } else - // dof on hex + // dof on hex { - // can only happen in 3d, - // but this case has - // already been covered - // above + // can only happen in 3d, + // but this case has + // already been covered + // above Assert (false, ExcNotImplemented()); return false; } - // we should not have gotten here + // we should not have gotten here Assert (false, ExcInternalError()); return false; } @@ -182,7 +184,7 @@ FE_TraceQ::get_constant_modes () const for (unsigned int i=0; idofs_per_cell; ++i) constant_modes(0,i) = true; return std::pair, std::vector > - (constant_modes, std::vector(1, 0)); + (constant_modes, std::vector(1, 0)); } @@ -195,7 +197,7 @@ FE_TraceQ::get_dpo_vector (const unsigned int deg) std::vector dpo(dim+1, 1U); dpo[dim]=0; dpo[0]=1; - for (unsigned int i=1;i(3) DEAL::FE_FaceP<2>(0) DEAL::FE_FaceP<2>(1) DEAL::FE_FaceP<2>(3) -DEAL::FE_DGRaviartThomas<2,2>(0) -DEAL::FE_DGRaviartThomas<2,2>(1) -DEAL::FE_DGNedelec<2,2>(0) -DEAL::FE_DGNedelec<2,2>(1) +DEAL::FE_DGRaviartThomas<2>(0) +DEAL::FE_DGRaviartThomas<2>(1) +DEAL::FE_DGNedelec<2>(0) +DEAL::FE_DGNedelec<2>(1) DEAL::FE_RaviartThomas<2>(0) DEAL::FE_RaviartThomas<2>(1) DEAL::FE_RaviartThomas<2>(2) DEAL::FESystem<2>[FE_RaviartThomas<2>(1)-FE_DGQ<2>(1)] DEAL::FE_Nedelec<2>(0) DEAL::FE_Nedelec<2>(1) -DEAL::FE_DGBDM<2,2>(1) -DEAL::FE_DGBDM<2,2>(2) +DEAL::FE_DGBDM<2>(1) +DEAL::FE_DGBDM<2>(2) DEAL::FE_BDM<2>(1) DEAL::FE_BDM<2>(2) DEAL::FE_RaviartThomasNodal<2>(0) @@ -1483,7 +1483,7 @@ DEAL::support on face 1: 4 5 6 7 DEAL::support on face 2: 8 9 10 11 DEAL::support on face 3: 12 13 14 15 DEAL:: -DEAL::fe_data[21]:FE_DGRaviartThomas<2,2>(0) +DEAL::fe_data[21]:FE_DGRaviartThomas<2>(0) DEAL::dofs_per_vertex=0 DEAL::dofs_per_line=0 DEAL::dofs_per_quad=4 @@ -1514,7 +1514,7 @@ DEAL::support on face 1: 0 1 2 3 DEAL::support on face 2: 0 1 2 3 DEAL::support on face 3: 0 1 2 3 DEAL:: -DEAL::fe_data[22]:FE_DGRaviartThomas<2,2>(1) +DEAL::fe_data[22]:FE_DGRaviartThomas<2>(1) DEAL::dofs_per_vertex=0 DEAL::dofs_per_line=0 DEAL::dofs_per_quad=12 @@ -1545,7 +1545,7 @@ DEAL::support on face 1: 0 1 2 3 4 5 6 7 8 9 10 11 DEAL::support on face 2: 0 1 2 3 4 5 6 7 8 9 10 11 DEAL::support on face 3: 0 1 2 3 4 5 6 7 8 9 10 11 DEAL:: -DEAL::fe_data[23]:FE_DGNedelec<2,2>(0) +DEAL::fe_data[23]:FE_DGNedelec<2>(0) DEAL::dofs_per_vertex=0 DEAL::dofs_per_line=0 DEAL::dofs_per_quad=4 @@ -1576,7 +1576,7 @@ DEAL::support on face 1: 0 1 2 3 DEAL::support on face 2: 0 1 2 3 DEAL::support on face 3: 0 1 2 3 DEAL:: -DEAL::fe_data[24]:FE_DGNedelec<2,2>(1) +DEAL::fe_data[24]:FE_DGNedelec<2>(1) DEAL::dofs_per_vertex=0 DEAL::dofs_per_line=0 DEAL::dofs_per_quad=12 @@ -1793,7 +1793,7 @@ DEAL::support on face 1: 2 3 4 5 6 7 DEAL::support on face 2: 0 1 2 3 4 5 DEAL::support on face 3: 0 1 2 3 6 7 8 9 10 11 DEAL:: -DEAL::fe_data[31]:FE_DGBDM<2,2>(1) +DEAL::fe_data[31]:FE_DGBDM<2>(1) DEAL::dofs_per_vertex=0 DEAL::dofs_per_line=0 DEAL::dofs_per_quad=8 @@ -1824,7 +1824,7 @@ DEAL::support on face 1: 0 1 2 3 4 5 6 7 DEAL::support on face 2: 0 1 2 3 4 5 6 7 DEAL::support on face 3: 0 1 2 3 4 5 6 7 DEAL:: -DEAL::fe_data[32]:FE_DGBDM<2,2>(2) +DEAL::fe_data[32]:FE_DGBDM<2>(2) DEAL::dofs_per_vertex=0 DEAL::dofs_per_line=0 DEAL::dofs_per_quad=14 @@ -2248,10 +2248,10 @@ DEAL::FE_FaceQ<3>(3) DEAL::FE_FaceP<3>(0) DEAL::FE_FaceP<3>(1) DEAL::FE_FaceP<3>(3) -DEAL::FE_DGRaviartThomas<3,3>(0) -DEAL::FE_DGRaviartThomas<3,3>(1) -DEAL::FE_DGNedelec<3,3>(0) -DEAL::FE_DGNedelec<3,3>(1) +DEAL::FE_DGRaviartThomas<3>(0) +DEAL::FE_DGRaviartThomas<3>(1) +DEAL::FE_DGNedelec<3>(0) +DEAL::FE_DGNedelec<3>(1) DEAL::FE_RaviartThomas<3>(0) DEAL::FE_RaviartThomas<3>(1) DEAL::FE_RaviartThomas<3>(2) @@ -3002,7 +3002,7 @@ DEAL::support on face 3: 30 31 32 33 34 35 36 37 38 39 DEAL::support on face 4: 40 41 42 43 44 45 46 47 48 49 DEAL::support on face 5: 50 51 52 53 54 55 56 57 58 59 DEAL:: -DEAL::fe_data[21]:FE_DGRaviartThomas<3,3>(0) +DEAL::fe_data[21]:FE_DGRaviartThomas<3>(0) DEAL::dofs_per_vertex=0 DEAL::dofs_per_line=0 DEAL::dofs_per_quad=0 @@ -3037,7 +3037,7 @@ DEAL::support on face 3: 0 1 2 3 4 5 DEAL::support on face 4: 0 1 2 3 4 5 DEAL::support on face 5: 0 1 2 3 4 5 DEAL:: -DEAL::fe_data[22]:FE_DGRaviartThomas<3,3>(1) +DEAL::fe_data[22]:FE_DGRaviartThomas<3>(1) DEAL::dofs_per_vertex=0 DEAL::dofs_per_line=0 DEAL::dofs_per_quad=0 @@ -3072,7 +3072,7 @@ DEAL::support on face 3: 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 DEAL::support on face 4: 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 DEAL::support on face 5: 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 DEAL:: -DEAL::fe_data[23]:FE_DGNedelec<3,3>(0) +DEAL::fe_data[23]:FE_DGNedelec<3>(0) DEAL::dofs_per_vertex=0 DEAL::dofs_per_line=0 DEAL::dofs_per_quad=0 @@ -3107,7 +3107,7 @@ DEAL::support on face 3: 0 1 2 3 4 5 6 7 8 9 10 11 DEAL::support on face 4: 0 1 2 3 4 5 6 7 8 9 10 11 DEAL::support on face 5: 0 1 2 3 4 5 6 7 8 9 10 11 DEAL:: -DEAL::fe_data[24]:FE_DGNedelec<3,3>(1) +DEAL::fe_data[24]:FE_DGNedelec<3>(1) DEAL::dofs_per_vertex=0 DEAL::dofs_per_line=0 DEAL::dofs_per_quad=0 diff --git a/tests/fe/get_fe_from_name_01.cc b/tests/fe/get_fe_from_name_01.cc new file mode 100644 index 0000000000..7eb8c1055a --- /dev/null +++ b/tests/fe/get_fe_from_name_01.cc @@ -0,0 +1,112 @@ +// --------------------------------------------------------------------- +// $Id$ +// +// Copyright (C) 2005 - 2013 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + + +#include "../tests.h" +#include +#include +#include +#include + +class Test +{ +public: + + void generate(const char *myname) + { + std::string name = myname; + + fe1 = FETools::get_fe_by_name<1,1>(name); + fe2 = FETools::get_fe_by_name<2,2>(name); + fe3 = FETools::get_fe_by_name<3,3>(name); + fe4 = FETools::get_fe_by_name<1,2>(name); + fe5 = FETools::get_fe_by_name<1,3>(name); + fe6 = FETools::get_fe_by_name<2,3>(name); + + deallog << "Read " << name << std::endl; + deallog << "Generated :" << std::endl; + deallog << fe1->get_name() << std::endl; + deallog << fe2->get_name() << std::endl; + deallog << fe3->get_name() << std::endl; + deallog << fe4->get_name() << std::endl; + deallog << fe5->get_name() << std::endl; + deallog << fe6->get_name() << std::endl; + + delete fe1; + delete fe2; + delete fe3; + delete fe4; + delete fe5; + delete fe6; + } + + FiniteElement<1,1> *fe1; + FiniteElement<2,2> *fe2; + FiniteElement<3,3> *fe3; + FiniteElement<1,2> *fe4; + FiniteElement<1,3> *fe5; + FiniteElement<2,3> *fe6; + +}; + +int main () +{ + std::ofstream logfile("output"); + deallog.attach(logfile); + + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + Test gen; + // Some of the finite element types are commented out, since their + // implementation or instantiation for codimension != 0 is still + // missing. This was opened as issue #92. + + // gen.generate("FE_Q_Hierarchical(1)"); + // gen.generate("FE_DGPNonparametric(1)"); + gen.generate("FE_DGQ(1)"); + gen.generate("FE_Q(1)"); + + // gen.generate("FE_Q_Hierarchical(2)"); + // gen.generate("FE_DGPNonparametric(2)"); + gen.generate("FE_DGQ(2)"); + gen.generate("FE_Q(2)"); + + // gen.generate("FESystem[FE_Q_Hierarchical(1)^2-FE_Q_Hierarchical(1)]"); + // gen.generate("FESystem[FE_DGPNonparametric(1)^2-FE_Q_Hierarchical(1)]"); + // gen.generate("FESystem[FE_DGQ(1)^2-FE_Q_Hierarchical(1)]"); + // gen.generate("FESystem[FE_Q(1)^2-FE_Q_Hierarchical(1)]"); + + // gen.generate("FESystem[FE_Q_Hierarchical(1)^2-FE_DGPNonparametric(1)]"); + // gen.generate("FESystem[FE_DGPNonparametric(1)^2-FE_DGPNonparametric(1)]"); + // gen.generate("FESystem[FE_DGQ(1)^2-FE_DGPNonparametric(1)]"); + // gen.generate("FESystem[FE_Q(1)^2-FE_DGPNonparametric(1)]"); + + // gen.generate("FESystem[FE_Q_Hierarchical(1)^2-FE_DGQ(1)]"); + // gen.generate("FESystem[FE_DGPNonparametric(1)^2-FE_DGQ(1)]"); + gen.generate("FESystem[FE_DGQ(1)^2-FE_DGQ(1)]"); + gen.generate("FESystem[FE_Q(1)^2-FE_DGQ(1)]"); + + + // gen.generate("FESystem[FE_Q_Hierarchical(1)^2-FE_Q(1)]"); + // gen.generate("FESystem[FE_DGPNonparametric(1)^2-FE_Q(1)]"); + gen.generate("FESystem[FE_DGQ(1)^2-FE_Q(1)]"); + gen.generate("FESystem[FE_Q(1)^2-FE_Q(1)]"); + + + return 0; +} + diff --git a/tests/fe/get_fe_from_name_01.output b/tests/fe/get_fe_from_name_01.output new file mode 100644 index 0000000000..f2427fe774 --- /dev/null +++ b/tests/fe/get_fe_from_name_01.output @@ -0,0 +1,65 @@ + +DEAL::Read FE_DGQ(1) +DEAL::Generated : +DEAL::FE_DGQ<1>(1) +DEAL::FE_DGQ<2>(1) +DEAL::FE_DGQ<3>(1) +DEAL::FE_DGQ<1,2>(1) +DEAL::FE_DGQ<1,3>(1) +DEAL::FE_DGQ<2,3>(1) +DEAL::Read FE_Q(1) +DEAL::Generated : +DEAL::FE_Q<1>(1) +DEAL::FE_Q<2>(1) +DEAL::FE_Q<3>(1) +DEAL::FE_Q<1,2>(1) +DEAL::FE_Q<1,3>(1) +DEAL::FE_Q<2,3>(1) +DEAL::Read FE_DGQ(2) +DEAL::Generated : +DEAL::FE_DGQ<1>(2) +DEAL::FE_DGQ<2>(2) +DEAL::FE_DGQ<3>(2) +DEAL::FE_DGQ<1,2>(2) +DEAL::FE_DGQ<1,3>(2) +DEAL::FE_DGQ<2,3>(2) +DEAL::Read FE_Q(2) +DEAL::Generated : +DEAL::FE_Q<1>(2) +DEAL::FE_Q<2>(2) +DEAL::FE_Q<3>(2) +DEAL::FE_Q<1,2>(2) +DEAL::FE_Q<1,3>(2) +DEAL::FE_Q<2,3>(2) +DEAL::Read FESystem[FE_DGQ(1)^2-FE_DGQ(1)] +DEAL::Generated : +DEAL::FESystem<1>[FE_DGQ<1>(1)^2-FE_DGQ<1>(1)] +DEAL::FESystem<2>[FE_DGQ<2>(1)^2-FE_DGQ<2>(1)] +DEAL::FESystem<3>[FE_DGQ<3>(1)^2-FE_DGQ<3>(1)] +DEAL::FESystem<1,2>[FE_DGQ<1,2>(1)^2-FE_DGQ<1,2>(1)] +DEAL::FESystem<1,3>[FE_DGQ<1,3>(1)^2-FE_DGQ<1,3>(1)] +DEAL::FESystem<2,3>[FE_DGQ<2,3>(1)^2-FE_DGQ<2,3>(1)] +DEAL::Read FESystem[FE_Q(1)^2-FE_DGQ(1)] +DEAL::Generated : +DEAL::FESystem<1>[FE_Q<1>(1)^2-FE_DGQ<1>(1)] +DEAL::FESystem<2>[FE_Q<2>(1)^2-FE_DGQ<2>(1)] +DEAL::FESystem<3>[FE_Q<3>(1)^2-FE_DGQ<3>(1)] +DEAL::FESystem<1,2>[FE_Q<1,2>(1)^2-FE_DGQ<1,2>(1)] +DEAL::FESystem<1,3>[FE_Q<1,3>(1)^2-FE_DGQ<1,3>(1)] +DEAL::FESystem<2,3>[FE_Q<2,3>(1)^2-FE_DGQ<2,3>(1)] +DEAL::Read FESystem[FE_DGQ(1)^2-FE_Q(1)] +DEAL::Generated : +DEAL::FESystem<1>[FE_DGQ<1>(1)^2-FE_Q<1>(1)] +DEAL::FESystem<2>[FE_DGQ<2>(1)^2-FE_Q<2>(1)] +DEAL::FESystem<3>[FE_DGQ<3>(1)^2-FE_Q<3>(1)] +DEAL::FESystem<1,2>[FE_DGQ<1,2>(1)^2-FE_Q<1,2>(1)] +DEAL::FESystem<1,3>[FE_DGQ<1,3>(1)^2-FE_Q<1,3>(1)] +DEAL::FESystem<2,3>[FE_DGQ<2,3>(1)^2-FE_Q<2,3>(1)] +DEAL::Read FESystem[FE_Q(1)^2-FE_Q(1)] +DEAL::Generated : +DEAL::FESystem<1>[FE_Q<1>(1)^2-FE_Q<1>(1)] +DEAL::FESystem<2>[FE_Q<2>(1)^2-FE_Q<2>(1)] +DEAL::FESystem<3>[FE_Q<3>(1)^2-FE_Q<3>(1)] +DEAL::FESystem<1,2>[FE_Q<1,2>(1)^2-FE_Q<1,2>(1)] +DEAL::FESystem<1,3>[FE_Q<1,3>(1)^2-FE_Q<1,3>(1)] +DEAL::FESystem<2,3>[FE_Q<2,3>(1)^2-FE_Q<2,3>(1)] diff --git a/tests/fe/get_name_02.cc b/tests/fe/get_name_02.cc new file mode 100644 index 0000000000..69774a0d3b --- /dev/null +++ b/tests/fe/get_name_02.cc @@ -0,0 +1,76 @@ +// --------------------------------------------------------------------- +// $Id$ +// +// Copyright (C) 2003 - 2013 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + + + +// test get_name() + +#include "../tests.h" +#include +#include +#include +#include +#include +#include +#include + +#include +#include + + +template +void test(const FiniteElement &fe) +{ + deallog << fe.get_name() << std::endl; +} + + +int +main() +{ + initlog(); + + { + FE_Q<2,3> fe(1); + test(fe); + } + { + FE_Q<2,3> fe(QGaussLobatto<1>(4)); + test(fe); + } + { + QGauss<1> quadrature_g(5); + FE_DGQArbitraryNodes<2,3> fe(quadrature_g); + test(fe); + } + { + QGaussLobatto<1> quadrature_gl(5); + FE_DGQArbitraryNodes<2,3> fe(quadrature_gl); + test(fe); + } + { + QGaussLog<1> quadrature(3); + FE_DGQArbitraryNodes<2,3> fe(quadrature); + test(fe); + } + + + + return 0; +} + + + diff --git a/tests/fe/get_name_02.output b/tests/fe/get_name_02.output new file mode 100644 index 0000000000..8cec2e7b53 --- /dev/null +++ b/tests/fe/get_name_02.output @@ -0,0 +1,6 @@ + +DEAL::FE_Q<2,3>(1) +DEAL::FE_Q<2,3>(QGaussLobatto(4)) +DEAL::FE_DGQArbitraryNodes<2>(QUnknownNodes(4)) +DEAL::FE_DGQArbitraryNodes<2>(QGaussLobatto(5)) +DEAL::FE_DGQArbitraryNodes<2>(QUnknownNodes(2)) diff --git a/tests/serialization/dof_handler_01.with_64bit_indices=off.output b/tests/serialization/dof_handler_01.with_64bit_indices=off.output index e9b71d8425..8d46a72f32 100644 --- a/tests/serialization/dof_handler_01.with_64bit_indices=off.output +++ b/tests/serialization/dof_handler_01.with_64bit_indices=off.output @@ -7,7 +7,7 @@ DEAL::0 0 0 0 0 0 0 0 0 0 0 0 0 0 10 0 0 1 11 12 6 7 3 4 9 10 0 0 14 14 0 0 0 0 DEAL::0 0 0 0 0 0 0 0 0 0 0 0 0 0 10 0 0 1 11 12 6 7 3 4 9 10 0 0 14 14 0 0 0 0 1 0 0 0 0 14 0 1 14 0 1 0 14 0 0 1 0 1 0 0 14 0 1 14 0 0 0 3 0 11 1 0 0 5 0 4294967295 4294967295 4294967295 4294967295 4294967295 0 0 1 0 4294967295 11 1 10 0 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 2 0 4294967295 4294967295 11 -2 20 0 0 1 3 4 2 3 4 6 7 5 6 7 9 10 8 9 10 11 12 13 4 0 2 5 8 13 -1 7 34 FESystem<1>[FE_Q<1>(2)-FE_Q<1>(1)] 23 Policy::Sequential<1,2> +2 20 0 0 1 3 4 2 3 4 6 7 5 6 7 9 10 8 9 10 11 12 13 4 0 2 5 8 13 -1 7 40 FESystem<1,2>[FE_Q<1,2>(2)-FE_Q<1,2>(1)] 23 Policy::Sequential<1,2> DEAL::0 0 0 0 0 0 0 0 0 0 0 0 0 0 75 0 0 1 2 114 115 116 127 128 129 148 149 150 32 33 34 39 40 41 162 163 164 172 173 174 56 57 58 9 10 11 73 74 75 12 13 14 84 85 86 157 158 159 177 178 179 167 168 169 182 183 184 50 51 52 108 109 110 53 54 55 102 103 104 19 20 21 99 100 101 105 106 107 111 112 113 0 0 187 187 0 0 0 0 1 0 0 0 0 187 0 1 187 0 1 0 187 0 0 1 0 1 0 0 187 0 1 187 0 0 0 3 0 11 1 0 0 22 0 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 0 0 2 0 4294967295 4294967295 11 @@ -19,7 +19,7 @@ DEAL::0 0 0 0 0 0 0 0 0 0 0 0 0 0 75 0 0 1 2 114 115 116 127 128 129 148 149 150 0 22 0 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 0 0 2 0 4294967295 4294967295 11 1 88 0 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 8 0 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 11 2 352 0 0 1 2 9 10 11 12 13 14 19 20 21 3 4 15 16 5 6 17 18 7 8 9 10 11 32 33 34 19 20 21 50 51 52 15 16 35 36 22 23 37 38 24 25 12 13 14 19 20 21 39 40 41 53 54 55 26 27 42 43 17 18 44 45 28 29 19 20 21 50 51 52 53 54 55 56 57 58 42 43 46 47 37 38 48 49 30 31 32 33 34 73 74 75 50 51 52 99 100 101 35 36 76 77 59 60 78 79 61 62 73 74 75 114 115 116 99 100 101 157 158 159 76 77 117 118 119 120 160 161 121 122 50 51 52 99 100 101 56 57 58 102 103 104 46 47 80 81 78 79 82 83 63 64 99 100 101 157 158 159 102 103 104 162 163 164 80 81 123 124 160 161 165 166 125 126 39 40 41 53 54 55 84 85 86 105 106 107 65 66 87 88 44 45 89 90 67 68 53 54 55 56 57 58 105 106 107 108 109 110 87 88 91 92 48 49 93 94 69 70 84 85 86 105 106 107 127 128 129 167 168 169 130 131 170 171 89 90 132 133 134 135 105 106 107 108 109 110 167 168 169 172 173 174 170 171 175 176 93 94 136 137 138 139 56 57 58 102 103 104 108 109 110 111 112 113 91 92 95 96 82 83 97 98 71 72 102 103 104 162 163 164 111 112 113 177 178 179 95 96 140 141 165 166 180 181 142 143 108 109 110 111 112 113 172 173 174 182 183 184 175 176 185 186 97 98 144 145 146 147 111 112 113 177 178 179 182 183 184 148 149 150 185 186 151 152 180 181 153 154 155 156 32 0 7 8 24 25 28 29 30 31 61 62 121 122 63 64 125 126 67 68 69 70 134 135 138 139 71 72 142 143 146 147 155 156 13 1 0 -3 0 0 112 0 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 5 6 22 23 59 60 119 120 3 4 26 27 65 66 130 131 117 118 123 124 140 141 151 152 132 133 136 137 144 145 153 154 35 36 46 47 91 92 175 176 44 45 48 49 82 83 165 166 15 16 42 43 17 18 37 38 76 77 80 81 78 79 160 161 87 88 170 171 89 90 93 94 95 96 185 186 97 98 180 181 21 36 FESystem<2>[FE_Q<2>(2)^2-FE_Q<2>(1)] 23 Policy::Sequential<2,3> +3 0 0 112 0 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 5 6 22 23 59 60 119 120 3 4 26 27 65 66 130 131 117 118 123 124 140 141 151 152 132 133 136 137 144 145 153 154 35 36 46 47 91 92 175 176 44 45 48 49 82 83 165 166 15 16 42 43 17 18 37 38 76 77 80 81 78 79 160 161 87 88 170 171 89 90 93 94 95 96 185 186 97 98 180 181 21 42 FESystem<2,3>[FE_Q<2,3>(2)^2-FE_Q<2,3>(1)] 23 Policy::Sequential<2,3> DEAL::0 0 0 0 0 0 0 0 0 0 0 0 0 0 500 0 0 1 2 3 1093 1094 1095 1096 1148 1149 1150 1151 1227 1228 1229 1230 1276 1277 1278 1279 1361 1362 1363 1364 1440 1441 1442 1443 1531 1532 1533 1534 146 147 148 149 171 172 173 174 211 212 213 214 1582 1583 1584 1585 1601 1602 1603 1604 1658 1659 1660 1661 1677 1678 1679 1680 1783 1784 1785 1786 1822 1823 1824 1825 1841 1842 1843 1844 1947 1948 1949 1950 2035 2036 2037 2038 339 340 341 342 309 310 311 312 352 353 354 355 2144 2145 2146 2147 2172 2173 2174 2175 2228 2229 2230 2231 398 399 400 401 25 26 27 28 513 514 515 516 29 30 31 32 577 578 579 580 33 34 35 36 665 666 667 668 1556 1557 1558 1559 1708 1709 1710 1711 1560 1561 1562 1563 1872 1873 1874 1875 1632 1633 1634 1635 1727 1728 1729 1730 1636 1637 1638 1639 1960 1961 1962 1963 1746 1747 1748 1749 2072 2073 2074 2075 1796 1797 1798 1799 1903 1904 1905 1906 1800 1801 1802 1803 1991 1992 1993 1994 1922 1923 1924 1925 2097 2098 2099 2100 2022 2023 2024 2025 2110 2111 2112 2113 322 323 324 325 816 817 818 819 279 280 281 282 929 930 931 932 275 276 277 278 855 856 857 858 292 293 294 295 803 804 805 806 296 297 298 299 942 943 944 945 326 327 328 329 868 869 870 871 2130 2131 2132 2133 2242 2243 2244 2245 2137 2138 2139 2140 2193 2194 2195 2196 2165 2166 2167 2168 2200 2201 2202 2203 2158 2159 2160 2161 2270 2271 2272 2273 2214 2215 2216 2217 2284 2285 2286 2287 2221 2222 2223 2224 2256 2257 2258 2259 1069 1070 1071 1072 382 383 384 385 1029 1030 1031 1032 394 395 396 397 1045 1046 1047 1048 390 391 392 393 68 69 70 71 912 913 914 915 790 791 792 793 964 965 966 967 64 65 66 67 786 787 788 789 838 839 840 841 890 891 892 893 72 73 74 75 842 843 844 845 916 917 918 919 986 987 988 989 2123 2124 2125 2126 2179 2180 2181 2182 2235 2236 2237 2238 2291 2292 2293 2294 2151 2152 2153 2154 2263 2264 2265 2266 2186 2187 2188 2189 2298 2299 2300 2301 2207 2208 2209 2210 2249 2250 2251 2252 2277 2278 2279 2280 2305 2306 2307 2308 1053 1054 1055 1056 1041 1042 1043 1044 1025 1026 1027 1028 386 387 388 389 1077 1078 1079 1080 1021 1022 1023 1024 1065 1066 1067 1068 378 379 380 381 1085 1086 1087 1088 1061 1062 1063 1064 1037 1038 1039 1040 374 375 376 377 85 86 87 88 1017 1018 1019 1020 1033 1034 1035 1036 1049 1050 1051 1052 1057 1058 1059 1060 1073 1074 1075 1076 1081 1082 1083 1084 1089 1090 1091 1092 0 0 2312 2312 0 0 0 0 1 0 0 0 0 2312 0 1 2312 0 1 0 2312 0 0 1 0 1 0 0 2312 0 1 2312 0 0 0 3 0 11 1 0 0 89 0 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 4294967295 0 0 3 0 4294967295 4294967295 4294967295 11