From f9eeb5d42d623f2cb0bdf459ddc1904d6a40191c Mon Sep 17 00:00:00 2001 From: bangerth Date: Mon, 29 Dec 2008 22:12:33 +0000 Subject: [PATCH] Fix doc markup in a number of places. Move one internal function from the private section of the class to the .cc file. git-svn-id: https://svn.dealii.org/trunk@18044 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/fe/fe_tools.h | 36 +- deal.II/deal.II/source/fe/fe_tools.cc | 468 +++++++++++++------------- 2 files changed, 255 insertions(+), 249 deletions(-) diff --git a/deal.II/deal.II/include/fe/fe_tools.h b/deal.II/deal.II/include/fe/fe_tools.h index 2b96003a85..56f64f6c71 100644 --- a/deal.II/deal.II/include/fe/fe_tools.h +++ b/deal.II/deal.II/include/fe/fe_tools.h @@ -358,27 +358,25 @@ class FETools * FiniteElement::prolongation * matrices. * - * @param fe The finite element + * @param fe The finite element * class for which we compute the * embedding matrices. * - * @param - * matrices A reference to - * RefinementCase::isotropic_refinement + * @param matrices A reference to + * RefinementCase::isotropic_refinement * vectors of FullMatrix * objects. Each vector * corresponds to one * RefinementCase @p * refinement_case and is of the * vector size - * GeometryInfo::n_children(refinement_case). This + * GeometryInfo::n_children(refinement_case). This * is the format used in * FiniteElement, where we want * to use this function mostly. * - * @param - * isotropic_only: set - * this true if you only + * @param isotropic_only Set + * to true if you only * want to compute matrices for * isotropic refinement. */ @@ -848,14 +846,14 @@ class FETools OutVector& u1_interpolated); /** - * Gives $(Id-I_h)z1$ for a given - * @p dof1-function @p z1, where $I_h$ + * Gives $(Id-I_h)z_1$ for a given + * @p dof1-function $z_1$, where $I_h$ * is the interpolation from @p fe1 - * to @p fe2. $(Id-I_h)z1$ is - * denoted by @p z1_difference. + * to @p fe2. The result $(Id-I_h)z_1$ is + * written into @p z1_difference. * * Note, that this function does - * not work on continuous + * not work for continuous * elements at hanging nodes. For * that case use the * @p interpolation_difference @@ -870,11 +868,11 @@ class FETools OutVector &z1_difference); /** - * Gives $(Id-I_h)z1$ for a given - * @p dof1-function @p z1, where $I_h$ + * Gives $(Id-I_h)z_1$ for a given + * @p dof1-function $z_1$, where $I_h$ * is the interpolation from @p fe1 - * to @p fe2. $(Id-I_h)z1$ is - * denoted by @p z1_difference. + * to @p fe2. The result $(Id-I_h)z_1$ is + * written into @p z1_difference. * @p constraints1 and * @p constraints2 are the * hanging node constraints @@ -938,7 +936,7 @@ class FETools * nodes, please use the * @p extrapolate function with * an additional - * @p ConstraintMatrix argument, + * ConstraintMatrix argument, * see below. * * Since this function operates @@ -1081,7 +1079,7 @@ class FETools * * The name must be in the form which * is returned by the - * @p FiniteElement::get_name + * FiniteElement::get_name * function, where a few * modifications are allowed: * diff --git a/deal.II/deal.II/source/fe/fe_tools.cc b/deal.II/deal.II/source/fe/fe_tools.cc index 76cf999ff3..592b4bc017 100644 --- a/deal.II/deal.II/source/fe/fe_tools.cc +++ b/deal.II/deal.II/source/fe/fe_tools.cc @@ -1550,6 +1550,242 @@ FETools::add_fe_name(const std::string& parameter_name, } +namespace internal +{ + namespace + { + template + FiniteElement* + get_fe_from_name (std::string &name) + { + // Extract the name of the + // finite element class, which only + // contains characters, numbers and + // underscores. + unsigned int name_end = + name.find_first_not_of(std::string("abcdefghijklmnopqrstuvwxyzABCDEFGHIJKLMNOPQRSTUVWXYZ0123456789_")); + const std::string name_part(name, 0, name_end); + name.erase(0, name_part.size()); + + // now things get a little more + // complicated: FESystem. it's + // more complicated, since we + // have to figure out what the + // base elements are. this can + // only be done recursively + if (name_part == "FESystem") + { + // next we have to get at the + // base elements. start with + // the first. wrap the whole + // block into try-catch to + // make sure we destroy the + // pointers we got from + // recursive calls if one of + // these calls should throw + // an exception + std::vector*> base_fes; + std::vector base_multiplicities; + try + { + // Now, just the [...] + // part should be left. + if (name.size() == 0 || name[0] != '[') + throw (std::string("Invalid first character in ") + name); + do + { + // Erase the + // leading '[' or '-' + name.erase(0,1); + // Now, the name of the + // first base element is + // first... Let's get it + base_fes.push_back (get_fe_from_name (name)); + // next check whether + // FESystem placed a + // multiplicity after + // the element name + if (name[0] == '^') + { + // yes. Delete the '^' + // and read this + // multiplicity + name.erase(0,1); + + const std::pair tmp + = Utilities::get_integer_at_position (name, 0); + name.erase(0, tmp.second); + // add to length, + // including the '^' + base_multiplicities.push_back (tmp.first); + } + else + // no, so + // multiplicity is + // 1 + base_multiplicities.push_back (1); + + // so that's it for + // this base + // element. base + // elements are + // separated by '-', + // and the list is + // terminated by ']', + // so loop while the + // next character is + // '-' + } + while (name[0] == '-'); + + // so we got to the end + // of the '-' separated + // list. make sure that + // we actually had a ']' + // there + if (name.size() == 0 || name[0] != ']') + throw (std::string("Invalid first character in ") + name); + name.erase(0,1); + // just one more sanity check + Assert ((base_fes.size() == base_multiplicities.size()) + && + (base_fes.size() > 0), + ExcInternalError()); + + // ok, apparently + // everything went ok. so + // generate the composed + // element + FiniteElement *system_element = 0; + switch (base_fes.size()) + { + case 1: + { + system_element = new FESystem(*base_fes[0], + base_multiplicities[0]); + break; + } + + case 2: + { + system_element = new FESystem(*base_fes[0], + base_multiplicities[0], + *base_fes[1], + base_multiplicities[1]); + break; + } + + case 3: + { + system_element = new FESystem(*base_fes[0], + base_multiplicities[0], + *base_fes[1], + base_multiplicities[1], + *base_fes[2], + base_multiplicities[2]); + break; + } + + default: + AssertThrow (false, ExcNotImplemented()); + } + + // now we don't need the + // list of base elements + // any more + for (unsigned int i=0; i(degree+1)) + // part should be left. + if (name.size() == 0 || name[0] != '(') + throw (std::string("Invalid first character in ") + name); + name.erase(0,1); + if (name[0] != 'Q') + { + const std::pair tmp + = Utilities::get_integer_at_position (name, 0); + name.erase(0, tmp.second+1); + return fe_name_map.find(name_part)->second->get(tmp.first); + } + else + { + unsigned int position = name.find('('); + const std::string quadrature_name(name, 0, position-1); + name.erase(0,position); + if (quadrature_name.compare("QGaussLobatto") == 0) + { + const std::pair tmp + = Utilities::get_integer_at_position (name, 0); + name.erase(0, tmp.second+1); +//TODO: Implement a get function taking Quadrature<1> in fe_tools.h. +//return fe_name_map.find(name_part)->second->get(QGaussLobatto<1>(tmp.first)); + AssertThrow (false, ExcNotImplemented()); + } + else + { + AssertThrow (false,ExcNotImplemented()); + } + } + } + + + // hm, if we have come thus far, we + // didn't know what to do with the + // string we got. so do as the docs + // say: raise an exception + AssertThrow (false, FETools::ExcInvalidFEName(name)); + + // make some compilers happy that + // do not realize that we can't get + // here after throwing + return 0; + } + } +} + + + + + template FiniteElement * FETools::get_fe_from_name (const std::string ¶meter_name) @@ -1601,7 +1837,7 @@ FETools::get_fe_from_name (const std::string ¶meter_name) try { - FiniteElement *fe = get_fe_from_name_aux (name); + FiniteElement *fe = internal::get_fe_from_name (name); // Make sure the auxiliary function // ate up all characters of the name. @@ -1625,238 +1861,10 @@ FETools::get_fe_from_name (const std::string ¶meter_name) // FiniteElement * // FETools::get_fe_from_name (const std::string ¶meter_name) // { -// return get_fe_from_name(parameter_name); +// return internal::get_fe_from_name(parameter_name); // } -template -FiniteElement* -FETools::get_fe_from_name_aux (std::string &name) -{ - // Extract the name of the - // finite element class, which only - // contains characters, numbers and - // underscores. - unsigned int name_end = - name.find_first_not_of(std::string("abcdefghijklmnopqrstuvwxyzABCDEFGHIJKLMNOPQRSTUVWXYZ0123456789_")); - const std::string name_part(name, 0, name_end); - name.erase(0, name_part.size()); - - // now things get a little more - // complicated: FESystem. it's - // more complicated, since we - // have to figure out what the - // base elements are. this can - // only be done recursively - if (name_part == "FESystem") - { - // next we have to get at the - // base elements. start with - // the first. wrap the whole - // block into try-catch to - // make sure we destroy the - // pointers we got from - // recursive calls if one of - // these calls should throw - // an exception - std::vector*> base_fes; - std::vector base_multiplicities; - try - { - // Now, just the [...] - // part should be left. - if (name.size() == 0 || name[0] != '[') - throw (std::string("Invalid first character in ") + name); - do - { - // Erase the - // leading '[' or '-' - name.erase(0,1); - // Now, the name of the - // first base element is - // first... Let's get it - base_fes.push_back (get_fe_from_name_aux (name)); - // next check whether - // FESystem placed a - // multiplicity after - // the element name - if (name[0] == '^') - { - // yes. Delete the '^' - // and read this - // multiplicity - name.erase(0,1); - - const std::pair tmp - = Utilities::get_integer_at_position (name, 0); - name.erase(0, tmp.second); - // add to length, - // including the '^' - base_multiplicities.push_back (tmp.first); - } - else - // no, so - // multiplicity is - // 1 - base_multiplicities.push_back (1); - - // so that's it for - // this base - // element. base - // elements are - // separated by '-', - // and the list is - // terminated by ']', - // so loop while the - // next character is - // '-' - } - while (name[0] == '-'); - - // so we got to the end - // of the '-' separated - // list. make sure that - // we actually had a ']' - // there - if (name.size() == 0 || name[0] != ']') - throw (std::string("Invalid first character in ") + name); - name.erase(0,1); - // just one more sanity check - Assert ((base_fes.size() == base_multiplicities.size()) - && - (base_fes.size() > 0), - ExcInternalError()); - - // ok, apparently - // everything went ok. so - // generate the composed - // element - FiniteElement *system_element = 0; - switch (base_fes.size()) - { - case 1: - { - system_element = new FESystem(*base_fes[0], - base_multiplicities[0]); - break; - } - - case 2: - { - system_element = new FESystem(*base_fes[0], - base_multiplicities[0], - *base_fes[1], - base_multiplicities[1]); - break; - } - - case 3: - { - system_element = new FESystem(*base_fes[0], - base_multiplicities[0], - *base_fes[1], - base_multiplicities[1], - *base_fes[2], - base_multiplicities[2]); - break; - } - - default: - AssertThrow (false, ExcNotImplemented()); - } - - // now we don't need the - // list of base elements - // any more - for (unsigned int i=0; i(degree+1)) - // part should be left. - if (name.size() == 0 || name[0] != '(') - throw (std::string("Invalid first character in ") + name); - name.erase(0,1); - if (name[0] != 'Q') - { - const std::pair tmp - = Utilities::get_integer_at_position (name, 0); - name.erase(0, tmp.second+1); - return fe_name_map.find(name_part)->second->get(tmp.first); - } - else - { - unsigned int position = name.find('('); - const std::string quadrature_name(name, 0, position-1); - name.erase(0,position); - if (quadrature_name.compare("QGaussLobatto") == 0) - { - const std::pair tmp - = Utilities::get_integer_at_position (name, 0); - name.erase(0, tmp.second+1); -//TODO: Implement a get function taking Quadrature<1> in fe_tools.h. - //return fe_name_map.find(name_part)->second->get(QGaussLobatto<1>(tmp.first)); - AssertThrow (false, ExcNotImplemented()); - } - else - { - AssertThrow (false,ExcNotImplemented()); - } - } - } - - - // hm, if we have come thus far, we - // didn't know what to do with the - // string we got. so do as the docs - // say: raise an exception - AssertThrow (false, ExcInvalidFEName(name)); - - // make some compilers happy that - // do not realize that we can't get - // here after throwing - return 0; -} - - - template void FETools:: -- 2.39.5