]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Be more lenient about template dimension parameters when parsing finite element names.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Wed, 9 Jul 2003 16:20:44 +0000 (16:20 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Wed, 9 Jul 2003 16:20:44 +0000 (16:20 +0000)
git-svn-id: https://svn.dealii.org/trunk@7867 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/fe/fe_tools.h
deal.II/deal.II/source/fe/fe_tools.cc
tests/bits/Makefile
tests/bits/fe_tools_10.cc [new file with mode: 0644]
tests/bits/fe_tools_11.cc [new file with mode: 0644]

index 7c4ff25437fcc8e107a611c287d32d8fffb3d429..e204ea7a148fb6e7e14f210e6865ce7579b493f3 100644 (file)
@@ -484,6 +484,31 @@ class FETools
                                      * @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
index 81ceaebc0bb083fb36366b1c8b39c7c473ce8e4e..461ffb1fde698d0db9917f3ea07de19a038f47f6 100644 (file)
@@ -113,8 +113,8 @@ namespace
   {
     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);
@@ -148,6 +148,52 @@ namespace
     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;
+  }
 }
 
 
@@ -1040,20 +1086,7 @@ FETools::get_fe_from_name (const std::string &name)
 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
@@ -1062,9 +1095,21 @@ FETools::get_fe_from_name_aux (const std::string &name)
                                   // 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));
@@ -1089,9 +1134,10 @@ FETools::get_fe_from_name_aux (const std::string &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);
@@ -1102,9 +1148,10 @@ FETools::get_fe_from_name_aux (const std::string &name)
       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);
@@ -1115,9 +1162,10 @@ FETools::get_fe_from_name_aux (const std::string &name)
       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);
@@ -1128,9 +1176,10 @@ FETools::get_fe_from_name_aux (const std::string &name)
       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);
@@ -1141,9 +1190,10 @@ FETools::get_fe_from_name_aux (const std::string &name)
       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);
@@ -1154,9 +1204,10 @@ FETools::get_fe_from_name_aux (const std::string &name)
       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);
@@ -1175,9 +1226,10 @@ FETools::get_fe_from_name_aux (const std::string &name)
                                     // 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
index 5b40a59f90d88b92cc805359bb842cc93c0697df..5135c7425079299251405f8f61e0b3a2e9c2a1a1 100644 (file)
@@ -84,9 +84,12 @@ fe_tools_06.exe         : fe_tools_06.g.$(OBJEXT)          $(libraries)
 fe_tools_07.exe         : fe_tools_07.g.$(OBJEXT)          $(libraries)
 fe_tools_08.exe         : fe_tools_08.g.$(OBJEXT)          $(libraries)
 fe_tools_09.exe         : fe_tools_09.g.$(OBJEXT)          $(libraries)
+fe_tools_10.exe         : fe_tools_10.g.$(OBJEXT)          $(libraries)
+fe_tools_11.exe         : fe_tools_11.g.$(OBJEXT)          $(libraries)
 gerold_1.exe            : gerold_1.g.$(OBJEXT)             $(libraries)
 roy_1.exe               : roy_1.g.$(OBJEXT)                $(libraries)
 denis_1.exe             : denis_1.g.$(OBJEXT)              $(libraries)
+point_inside_1.exe      : point_inside_1.g.$(OBJEXT)      $(libraries)
 unit_support_points.exe : unit_support_points.g.$(OBJEXT)  $(libraries)
 parameter_handler_1.exe : parameter_handler_1.g.$(OBJEXT)  $(libraries)
 parameter_handler_2.exe : parameter_handler_2.g.$(OBJEXT)  $(libraries)
@@ -95,7 +98,7 @@ block_sparse_matrix_1.exe:block_sparse_matrix_1.g.$(OBJEXT) $(libraries)
 
 
 tests = anna_1 anna_2 anna_3 anna_4 anna_5 anna_6 \
-       geometry_info_1 \
+       geometry_info_1 point_inside_1 \
        data_out_01 data_out_02 data_out_faces_01 data_out_faces_02 \
        data_out_rotation_01 data_out_rotation_02 data_out_stack_01 \
        dof_tools_00a \
@@ -111,6 +114,7 @@ tests = anna_1 anna_2 anna_3 anna_4 anna_5 anna_6 \
        fe_tools_01a fe_tools_01b fe_tools_01c \
        fe_tools_02 fe_tools_03 fe_tools_04 \
        fe_tools_05 fe_tools_06 fe_tools_07 fe_tools_08 fe_tools_09 \
+       fe_tools_10 fe_tools_11\
        roy_1 \
         denis_1 \
        unit_support_points parameter_handler_1 parameter_handler_2 \
diff --git a/tests/bits/fe_tools_10.cc b/tests/bits/fe_tools_10.cc
new file mode 100644 (file)
index 0000000..b2c6c1a
--- /dev/null
@@ -0,0 +1,71 @@
+//----------------------------  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;
+}
+
diff --git a/tests/bits/fe_tools_11.cc b/tests/bits/fe_tools_11.cc
new file mode 100644 (file)
index 0000000..24f59b2
--- /dev/null
@@ -0,0 +1,70 @@
+//----------------------------  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;
+}
+

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.