* @p{FETools::ExcInvalidFEName}
* is thrown.
*
+ * There is one exception,
+ * however, where the names must
+ * not match exactly: while the
+ * finite elements write the
+ * space dimension in the form of
+ * a template argument after the
+ * name of the class (for example
+ * @p{FE_Q<2>}, you can omit the
+ * dimension argument altogether,
+ * or replace it with the string
+ * @p{<dim>}. The reason is that
+ * the dimension argument may be
+ * cumbersome if the name of a
+ * finite element is given in an
+ * input file that may be used to
+ * control operation of the
+ * program in different space
+ * dimensions. Running the
+ * program in another space
+ * dimension would then require
+ * changing the input file as
+ * well. With above exception,
+ * there is a canonical spelling
+ * that doesn't require this.
+ *
* The function returns a pointer
* to a newly create finite
* element. It is in the callers
{
Assert (position < name.size(), ExcInternalError());
- std::string test_string (name.begin()+position,
- name.end());
+ const std::string test_string (name.begin()+position,
+ name.end());
#ifdef HAVE_STD_STRINGSTREAM
std::istringstream str(test_string);
else
return std::make_pair (-1, static_cast<unsigned int>(-1));
}
+
+
+
+ // return how many characters
+ // starting at the given position
+ // of the string match either the
+ // generic string "<dim>" or the
+ // specialized string with "dim"
+ // replaced with the numeric value
+ // of the template argument
+ template <int dim>
+ unsigned int match_dimension (const std::string &name,
+ const unsigned int position)
+ {
+ if (position >= name.size())
+ return 0;
+
+ if ((position+5 < name.size())
+ &&
+ (name[position] == '<')
+ &&
+ (name[position+1] == 'd')
+ &&
+ (name[position+2] == 'i')
+ &&
+ (name[position+3] == 'm')
+ &&
+ (name[position+4] == '>'))
+ return 5;
+
+ Assert (dim<10, ExcNotImplemented());
+ const char dim_char = '0'+dim;
+
+ if ((position+3 < name.size())
+ &&
+ (name[position] == '<')
+ &&
+ (name[position+1] == dim_char)
+ &&
+ (name[position+2] == '>'))
+ return 3;
+
+ // some other string that doesn't
+ // match
+ return 0;
+ }
}
template <int dim>
std::pair<FiniteElement<dim> *, unsigned int>
FETools::get_fe_from_name_aux (const std::string &name)
-{
-#ifdef HAVE_STD_STRINGSTREAM
- std::ostringstream s;
-#else
- std::ostrstream s;
-#endif
-
- s << '<' << dim << '>';
-#ifndef HAVE_STD_STRINGSTREAM
- s << std::ends;
-#endif
-
- const std::string dim_name = s.str();
-
+{
// so, let's see what's at position
// 0 of this string, and create a
// respective finite element
// make sure we don't match FE_Q
// when it's actually a
// FE_Q_Hierarchic
- if (match_at_string_start (name, std::string("FE_Q_Hierarchical")+dim_name))
+ if (match_at_string_start (name, std::string("FE_Q_Hierarchical")))
{
- unsigned int position = (std::string("FE_Q_Hierarchical")+dim_name).size();
+ unsigned int position = std::string("FE_Q_Hierarchical").size();
+ // as described in the
+ // documentation, at this point
+ // we may have either a) no
+ // dimension specification, b)
+ // the correct dimension (like
+ // <2>), or c) a generic
+ // <dim>. check how many of
+ // these characters match, and
+ // advance the cursor by the
+ // respective amount
+ position += match_dimension<dim> (name, position);
+
// make sure the next character
// is an opening parenthesis
AssertThrow (name[position] == '(', ExcInvalidFEName(name));
}
// check other possibilities in
// exactly the same way
- else if (match_at_string_start (name, std::string("FE_RaviartThomas")+dim_name))
+ else if (match_at_string_start (name, std::string("FE_RaviartThomas")))
{
- unsigned int position = (std::string("FE_RaviartThomas")+dim_name).size();
+ unsigned int position = std::string("FE_RaviartThomas").size();
+ position += match_dimension<dim> (name, position);
AssertThrow (name[position] == '(', ExcInvalidFEName(name));
++position;
const std::pair<int,unsigned int> tmp = get_integer (name, position);
return std::make_pair (new FE_RaviartThomas<dim>(tmp.first),
position);
}
- else if (match_at_string_start (name, std::string("FE_Nedelec")+dim_name))
+ else if (match_at_string_start (name, std::string("FE_Nedelec")))
{
- unsigned int position = (std::string("FE_Nedelec")+dim_name).size();
+ unsigned int position = std::string("FE_Nedelec").size();
+ position += match_dimension<dim> (name, position);
AssertThrow (name[position] == '(', ExcInvalidFEName(name));
++position;
const std::pair<int,unsigned int> tmp = get_integer (name, position);
return std::make_pair (new FE_Nedelec<dim>(tmp.first),
position);
}
- else if (match_at_string_start (name, std::string("FE_DGPNonparametric")+dim_name))
+ else if (match_at_string_start (name, std::string("FE_DGPNonparametric")))
{
- unsigned int position = (std::string("FE_DGPNonparametric")+dim_name).size();
+ unsigned int position = std::string("FE_DGPNonparametric").size();
+ position += match_dimension<dim> (name, position);
AssertThrow (name[position] == '(', ExcInvalidFEName(name));
++position;
const std::pair<int,unsigned int> tmp = get_integer (name, position);
return std::make_pair (new FE_DGPNonparametric<dim>(tmp.first),
position);
}
- else if (match_at_string_start (name, std::string("FE_DGP")+dim_name))
+ else if (match_at_string_start (name, std::string("FE_DGP")))
{
- unsigned int position = (std::string("FE_DGP")+dim_name).size();
+ unsigned int position = std::string("FE_DGP").size();
+ position += match_dimension<dim> (name, position);
AssertThrow (name[position] == '(', ExcInvalidFEName(name));
++position;
const std::pair<int,unsigned int> tmp = get_integer (name, position);
return std::make_pair (new FE_DGP<dim>(tmp.first),
position);
}
- else if (match_at_string_start (name, std::string("FE_DGQ")+dim_name))
+ else if (match_at_string_start (name, std::string("FE_DGQ")))
{
- unsigned int position = (std::string("FE_DGQ")+dim_name).size();
+ unsigned int position = std::string("FE_DGQ").size();
+ position += match_dimension<dim> (name, position);
AssertThrow (name[position] == '(', ExcInvalidFEName(name));
++position;
const std::pair<int,unsigned int> tmp = get_integer (name, position);
return std::make_pair (new FE_DGQ<dim>(tmp.first),
position);
}
- else if (match_at_string_start (name, std::string("FE_Q")+dim_name))
+ else if (match_at_string_start (name, std::string("FE_Q")))
{
- unsigned int position = (std::string("FE_Q")+dim_name).size();
+ unsigned int position = std::string("FE_Q").size();
+ position += match_dimension<dim> (name, position);
AssertThrow (name[position] == '(', ExcInvalidFEName(name));
++position;
const std::pair<int,unsigned int> tmp = get_integer (name, position);
// have to figure out what the
// base elements are. this can
// only be done recursively
- if (match_at_string_start (name, std::string("FESystem")+dim_name))
+ if (match_at_string_start (name, std::string("FESystem")))
{
- unsigned int position = (std::string("FESystem")+dim_name).size();
+ unsigned int position = std::string("FESystem").size();
+ position += match_dimension<dim> (name, position);
// FESystem puts the names of
// the basis elements into
--- /dev/null
+//---------------------------- fe_tools_10.cc ---------------------------
+// $Id$
+// Version: $Name$
+//
+// Copyright (C) 2003 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.
+//
+//---------------------------- fe_tools_10.cc ---------------------------
+
+#include "fe_tools_common.cc"
+#include <lac/sparsity_pattern.h>
+
+// check
+// FETools::get_fe_from_name
+// like fe_tools_09, but with the difference that we use the generic
+// "<dim>" marker in the finite element name, instead of the one with
+// the concrete dimension (see the documentation)
+
+
+std::string output_file_name = "fe_tools_10.output";
+
+template <int dim>
+std::string modify_name (const std::string &name)
+{
+ std::string new_name = name;
+ std::string dim_name = std::string("<");
+ dim_name += '0'+dim;
+ dim_name += '>';
+
+ unsigned int pos;
+ while ((pos = new_name.find(dim_name)) != std::string::npos)
+ new_name.replace (pos, 3, "<dim>");
+
+ return new_name;
+}
+
+
+
+template <int dim>
+void
+check_this (const FiniteElement<dim> &fe1,
+ const FiniteElement<dim> &fe2)
+{
+ FiniteElement<dim> *p1, *p2;
+
+ // check that the name of the fe
+ // and the name of the fe that we
+ // re-create from this name are
+ // identitical. this is also a
+ // pretty good indication that the
+ // two FEs are actually the same
+ deallog << modify_name<dim> (fe1.get_name());
+ p1 = FETools::get_fe_from_name<dim> (modify_name<dim> (fe1.get_name()));
+ Assert (fe1.get_name() == p1->get_name(),
+ ExcInternalError());
+ deallog << " ok" << std::endl;
+ delete p1;
+
+ // same for fe2
+ deallog << modify_name<dim> (fe2.get_name());
+ p2 = FETools::get_fe_from_name<dim> (modify_name<dim> (fe2.get_name()));
+ Assert (fe2.get_name() == p2->get_name(),
+ ExcInternalError());
+ deallog << " ok" << std::endl;
+ delete p2;
+}
+
--- /dev/null
+//---------------------------- fe_tools_11.cc ---------------------------
+// $Id$
+// Version: $Name$
+//
+// Copyright (C) 2003 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.
+//
+//---------------------------- fe_tools_11.cc ---------------------------
+
+#include "fe_tools_common.cc"
+#include <lac/sparsity_pattern.h>
+
+// check
+// FETools::get_fe_from_name
+// like fe_tools_09 and fe_tools_10, but this time with no dimension
+// marker at all (see the documentation)
+
+
+std::string output_file_name = "fe_tools_11.output";
+
+template <int dim>
+std::string modify_name (const std::string &name)
+{
+ std::string new_name = name;
+ std::string dim_name = std::string("<");
+ dim_name += '0'+dim;
+ dim_name += '>';
+
+ unsigned int pos;
+ while ((pos = new_name.find(dim_name)) != std::string::npos)
+ new_name.replace (pos, 3, "");
+
+ return new_name;
+}
+
+
+
+template <int dim>
+void
+check_this (const FiniteElement<dim> &fe1,
+ const FiniteElement<dim> &fe2)
+{
+ FiniteElement<dim> *p1, *p2;
+
+ // check that the name of the fe
+ // and the name of the fe that we
+ // re-create from this name are
+ // identitical. this is also a
+ // pretty good indication that the
+ // two FEs are actually the same
+ deallog << modify_name<dim> (fe1.get_name());
+ p1 = FETools::get_fe_from_name<dim> (modify_name<dim> (fe1.get_name()));
+ Assert (fe1.get_name() == p1->get_name(),
+ ExcInternalError());
+ deallog << " ok" << std::endl;
+ delete p1;
+
+ // same for fe2
+ deallog << modify_name<dim> (fe2.get_name());
+ p2 = FETools::get_fe_from_name<dim> (modify_name<dim> (fe2.get_name()));
+ Assert (fe2.get_name() == p2->get_name(),
+ ExcInternalError());
+ deallog << " ok" << std::endl;
+ delete p2;
+}
+