nedelec.exe : nedelec.g.$(OBJEXT) $(libraries)
nedelec_2.exe : nedelec_2.g.$(OBJEXT) $(libraries)
nedelec_3.exe : nedelec_3.g.$(OBJEXT) $(libraries)
+non_primitive_1.exe : non_primitive_1.g.$(OBJEXT) $(libraries)
+non_primitive_2.exe : non_primitive_2.g.$(OBJEXT) $(libraries)
+numbering.exe : numbering.g.$(OBJEXT) $(libraries)
+q_1.exe : q_1.g.$(OBJEXT) $(libraries)
+q_2.exe : q_2.g.$(OBJEXT) $(libraries)
+q_3.exe : q_3.g.$(OBJEXT) $(libraries)
rt_1.exe : rt_1.g.$(OBJEXT) $(libraries)
rt_2.exe : rt_2.g.$(OBJEXT) $(libraries)
rt_3.exe : rt_3.g.$(OBJEXT) $(libraries)
rt_4.exe : rt_4.g.$(OBJEXT) $(libraries)
-non_primitive_1.exe : non_primitive_1.g.$(OBJEXT) $(libraries)
-non_primitive_2.exe : non_primitive_2.g.$(OBJEXT) $(libraries)
-numbering.exe : numbering.g.$(OBJEXT) $(libraries)
shapes.exe : shapes.g.$(OBJEXT) $(libraries)
transfer.exe : transfer.g.$(OBJEXT) $(libraries)
up_and_down.exe : up_and_down.g.$(OBJEXT) $(libraries)
tests = fe_data_test fe_traits fe_tools_test mapping \
mapping_c1 shapes derivatives numbering mapping_q1_eulerian \
transfer internals \
+ q_1 q_2 q_3 \
non_primitive_1 non_primitive_2 \
nedelec nedelec_2 nedelec_3 up_and_down
--- /dev/null
+//---------------------------- q_1.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
+// fuqher information on this license.
+//
+//---------------------------- q_1.cc ---------------------------
+
+
+// Just output the constraint matrices of the FE_Q element. Test
+// introduced when we started to compute them on the fly, rather than
+// precomputing them for a number of elements and storing them in a
+// table
+
+#include <base/logstream.h>
+#include <fe/fe_q.h>
+
+#include <fstream>
+#include <string>
+
+#define PRECISION 5
+
+
+
+template<int dim>
+void
+test(const unsigned int degree)
+{
+ deallog << "FE_Q<" << dim << "> (" << degree << ")"
+ << std::endl;
+
+ FE_Q<dim> fe_q(degree);
+ const FullMatrix<double> & constraints = fe_q.constraints();
+
+ for (unsigned int i=0; i<constraints.m(); ++i)
+ {
+ for (unsigned int j=0; j<constraints.n(); ++j)
+ deallog << constraints(i,j) << ' ';
+ deallog << std::endl;
+ }
+
+ deallog << std::endl;
+}
+
+
+int
+main()
+{
+ std::ofstream logfile ("q_1.output");
+ logfile.precision (PRECISION);
+ logfile.setf(std::ios::fixed);
+ deallog.attach(logfile);
+ deallog.depth_console(0);
+
+ // no constraints in 1d, but we had
+ // the matrices precomputed up to
+ // Q4 for 2d and Q2 for 3d
+ for (unsigned int degree=1; degree<=4; ++degree)
+ test<2>(degree);
+
+ for (unsigned int degree=1; degree<=2; ++degree)
+ test<3>(degree);
+
+ return 0;
+}
+
+
+
--- /dev/null
+//---------------------------- q_2.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
+// fuqher information on this license.
+//
+//---------------------------- q_2.cc ---------------------------
+
+
+// Just output the embedding matrices of the FE_Q element. Test
+// introduced when we started to compute them on the fly, rather than
+// precomputing them for a number of elements and storing them in a
+// table
+
+#include <base/logstream.h>
+#include <fe/fe_q.h>
+
+#include <fstream>
+#include <string>
+
+#define PRECISION 5
+
+
+
+template<int dim>
+void
+test(const unsigned int degree)
+{
+ deallog << "FE_Q<" << dim << "> (" << degree << ")"
+ << std::endl;
+
+ FE_Q<dim> fe_q(degree);
+
+ for (unsigned int c=0; c<GeometryInfo<dim>::children_per_cell; ++c)
+ {
+ const FullMatrix<double> & m = fe_q.prolongate(c);
+
+ for (unsigned int i=0; i<m.m(); ++i)
+ {
+ for (unsigned int j=0; j<m.n(); ++j)
+ deallog << m(i,j) << ' ';
+ deallog << std::endl;
+ }
+
+ deallog << std::endl;
+ }
+
+ deallog << std::endl;
+}
+
+
+int
+main()
+{
+ std::ofstream logfile ("q_2.output");
+ logfile.precision (PRECISION);
+ logfile.setf(std::ios::fixed);
+ deallog.attach(logfile);
+ deallog.depth_console(0);
+
+ // we had the matrices precomputed
+ // up to Q4 for 1d, Q3 for 2d and
+ // Q2 for 3d
+ for (unsigned int degree=1; degree<=4; ++degree)
+ test<1>(degree);
+
+ for (unsigned int degree=1; degree<=3; ++degree)
+ test<2>(degree);
+
+ for (unsigned int degree=1; degree<=2; ++degree)
+ test<3>(degree);
+
+ return 0;
+}
+
+
+
--- /dev/null
+//---------------------------- q_3.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
+// fuqher information on this license.
+//
+//---------------------------- q_3.cc ---------------------------
+
+
+// Just output the restriction matrices of the FE_Q element. Test
+// introduced when we started to compute them on the fly, rather than
+// precomputing them for a number of elements and storing them in a
+// table
+
+#include <base/logstream.h>
+#include <fe/fe_q.h>
+
+#include <fstream>
+#include <string>
+
+#define PRECISION 5
+
+
+
+template<int dim>
+void
+test(const unsigned int degree)
+{
+ deallog << "FE_Q<" << dim << "> (" << degree << ")"
+ << std::endl;
+
+ FE_Q<dim> fe_q(degree);
+
+ for (unsigned int c=0; c<GeometryInfo<dim>::children_per_cell; ++c)
+ {
+ const FullMatrix<double> & m = fe_q.restrict(c);
+
+ for (unsigned int i=0; i<m.m(); ++i)
+ for (unsigned int j=0; j<m.n(); ++j)
+ if (m(i,j)!=0)
+ deallog << '[' << i << ',' << j << ',' << m(i,j) << ']';
+
+ deallog << std::endl;
+ }
+
+ deallog << std::endl;
+}
+
+
+int
+main()
+{
+ std::ofstream logfile ("q_3.output");
+ logfile.precision (PRECISION);
+ logfile.setf(std::ios::fixed);
+ deallog.attach(logfile);
+ deallog.depth_console(0);
+
+ // we had the matrices precomputed
+ // up to Q4 for 1d, 2d and 3d
+ for (unsigned int degree=1; degree<=4; ++degree)
+ test<1>(degree);
+
+ for (unsigned int degree=1; degree<=4; ++degree)
+ test<2>(degree);
+
+ for (unsigned int degree=1; degree<=4; ++degree)
+ test<3>(degree);
+
+ return 0;
+}
+
+
+