]> https://gitweb.dealii.org/ - dealii.git/commitdiff
add a new class FETools that provides interpolation, back_interpolation, and interpol...
authorRalf Hartmann <Ralf.Hartmann@dlr.de>
Fri, 26 May 2000 16:25:03 +0000 (16:25 +0000)
committerRalf Hartmann <Ralf.Hartmann@dlr.de>
Fri, 26 May 2000 16:25:03 +0000 (16:25 +0000)
git-svn-id: https://svn.dealii.org/trunk@2959 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/fe/fe_tools.h [new file with mode: 0644]
deal.II/deal.II/include/numerics/matrices.h
deal.II/deal.II/source/fe/fe_tools.cc [new file with mode: 0644]
deal.II/deal.II/source/numerics/matrices.cc

diff --git a/deal.II/deal.II/include/fe/fe_tools.h b/deal.II/deal.II/include/fe/fe_tools.h
new file mode 100644 (file)
index 0000000..06ec845
--- /dev/null
@@ -0,0 +1,206 @@
+//----------------------------  fe_tools.h  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 1998, 1999, 2000 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.h  ---------------------------
+#ifndef __deal2__fe_tools_H
+#define __deal2__fe_tools_H
+
+
+
+template <typename number> class FullMatrix;
+template <int dim> class FiniteElement;
+template <int dim> class DoFHandler;
+template <typename number> class Vector;
+
+#include <base/exceptions.h>
+
+/**
+ * This class performs interpolations and extrapolations of discrete
+ * functions of one #FiniteElement# #fe1# to another #FiniteElement#
+ * #fe2#.
+ *
+ * It also provides the local interpolation matrices that interpolate
+ * on each cell. Furthermore it provides the difference matrix
+ * $id-I_h$ that is needed for evaluating $(id-I_h)z$ for e.g. the
+ * dual solution $z$.
+ *
+ * @author Ralf Hartmann, 2000
+ */
+class FETools
+{
+  public:
+                                  /**
+                                   * Gives the interpolation matrix
+                                   * that interpolates a #fe1#-
+                                   * function to a #fe2#-function on
+                                   * each cell. The interpolation_matrix
+                                   * needs to be of size
+                                   * #(fe2.dofs_per_cell, fe1.dofs_per_cell)#.
+                                   *
+                                   * Note, that if the finite element
+                                   * space #fe1# is a subset of
+                                   * the finite element space
+                                   * #fe2# than the #interpolation_matrix#
+                                   * is an embedding matrix.
+                                   */
+  template <int dim, typename number>
+    static void get_interpolation_matrix(const FiniteElement<dim> &fe1,
+                                        const FiniteElement<dim> &fe2,
+                                        FullMatrix<number> &interpolation_matrix);
+    
+                                  /**
+                                   * Gives the interpolation matrix
+                                   * that interpolates a #fe1#-
+                                   * function to a #fe2#-function, and
+                                   * interpolates this to a second
+                                   * #fe1#-function on
+                                   * each cell. The interpolation_matrix
+                                   * needs to be of size
+                                   * #(fe1.dofs_per_cell, fe1.dofs_per_cell)#.
+                                   *
+                                   * Note, that this function only
+                                   * makes sense if the finite element
+                                   * space due to #fe1# is not a subset of
+                                   * the finite element space due to
+                                   * #fe2#, as if it were a subset then
+                                   * the #interpolation_matrix# would be 
+                                   * only the unit matrix.
+                                   */
+  template <int dim, typename number>
+    static void get_back_interpolation_matrix(const FiniteElement<dim> &fe1,
+                                             const FiniteElement<dim> &fe2,
+                                             FullMatrix<number> &interpolation_matrix);
+
+                                  /**
+                                   * Gives the unit matrix minus the
+                                   * back interpolation matrix.
+                                   * The #difference_matrix#
+                                   * needs to be of size
+                                   * #(fe1.dofs_per_cell, fe1.dofs_per_cell)#.
+                                   *
+                                   * This function gives
+                                   * the matrix that transforms a
+                                   * #fe1# function $z$ to $z-I_hz$
+                                   * where $I_h$ denotes the interpolation
+                                   * operator from the #fe1# space to
+                                   * the #fe2# space. This matrix hence
+                                   * is useful to evaluate
+                                   * error-representations where $z$
+                                   * denotes the dual solution.
+                                   */
+  template <int dim, typename number>
+    static void get_interpolation_difference_matrix(
+      const FiniteElement<dim> &fe1,
+      const FiniteElement<dim> &fe2,
+      FullMatrix<number> &difference_matrix);
+
+                                  /**
+                                   * Gives the interpolation of a the
+                                   * #dof1#-function #u1# to a
+                                   * #dof2#-function #u2#. #dof1# and
+                                   * #dof2# need to be #DoFHandler#s
+                                   * based on the same triangulation.
+                                   *
+                                   * If the elements #fe1# and #fe2#
+                                   * are either both continuous or
+                                   * both discontinuous then this
+                                   * interpolation is the usual point
+                                   * interpolation. The same is true
+                                   * if #fe1# is a continuous and
+                                   * #fe2# is a discontinuous finite
+                                   * element. For the case that #fe1#
+                                   * is a discontinuous and #fe2# is
+                                   * a continuous finite element
+                                   * there is no point interpolation
+                                   * defined at the discontinuities.
+                                   * Therefore the meanvalue is taken
+                                   * at the DoF values on the
+                                   * discontinuities.
+                                   */
+  template <int dim, typename number>
+    static void interpolate(const DoFHandler<dim> &dof1,
+                           const Vector<number> &u1,
+                           const DoFHandler<dim> &dof2,
+                           Vector<number> &u2);
+    
+                                  /**
+                                   * Gives the interpolation of the #fe1#-
+                                   * function #u1# to a #fe2#-function, and
+                                   * interpolates this to a second
+                                   * #fe1#-function named #u1_interpolated#.
+                                   *
+                                   * Note, that this function only
+                                   * makes sense if the finite element
+                                   * space due to #fe1# is not a subset of
+                                   * the finite element space due to
+                                   * #fe2#, as if it were a subset then
+                                   * #u1_interpolated# would be equal to #u1#.
+                                   */
+  template <int dim, typename number>
+    static void back_interpolate(const DoFHandler<dim> &dof1,
+                                const Vector<number> &u1,
+                                const FiniteElement<dim> &fe2,
+                                Vector<number> &u1_interpolated);
+
+                                  /**
+                                   * Gives $(Id-I_h)z2$ for a given
+                                   * #fe2#-function #z2#, where $I_h$
+                                   * is the interpolation from #fe2#
+                                   * to #fe1#. $(Id-I_h)z2$ is
+                                   * denoted by #z2_difference#.
+                                   */
+  template <int dim, typename number>
+    static void interpolation_difference(const DoFHandler<dim> &dof1,
+                                        const Vector<number> &z1,
+                                        const FiniteElement<dim> &fe2,
+                                        Vector<number> &z1_difference);    
+    
+                                  /**
+                                   * Gives the patchwise
+                                   * extrapolation of a #dof1#
+                                   * function #z1# to a #dof2#
+                                   * function #z2#.  #dof1# and
+                                   * #dof2# need to be #DoFHandler#
+                                   * based on the same triangulation.
+                                   *
+                                   * This function is interesting for
+                                   * e.g. extrapolating patchwise a
+                                   * piecewise linear dual solution
+                                   * to a piecewise quadratic dual
+                                   * solution.
+                                   */
+  template <int dim, typename number>
+    static void extrapolate(const DoFHandler<dim> &dof1,
+                           const Vector<number> &z1,
+                           const DoFHandler<dim> &dof2,
+                           Vector<number> &z2);    
+  
+  
+                                  /**
+                                   * Exception
+                                   */
+  DeclException0 (ExcTriangulationMismatch);
+  
+                                  /**
+                                   * Exception
+                                   */
+  DeclException4 (ExcMatrixDimensionMismatch,
+                 int, int, int, int,
+                 << "This is a " << arg1 << "x" << arg2 << " matrix, "
+                 << "but should be a " << arg1 << "x" << arg2 << " matrix.");
+};
+
+
+
+/*----------------------------   fe_tools.h     ---------------------------*/
+/* end of #ifndef __deal2__fe_tools_H */
+#endif
+/*----------------------------   fe_tools.h     ---------------------------*/
index c12b39d506146884cd7f6f66fb69bcbadbeee7a7..b8c183e9894b7d6f94a935924894bfd014480163 100644 (file)
@@ -329,27 +329,6 @@ class MatrixCreator
                                       Vector<double>           &rhs_vector,
                                       const Function<dim>      *a = 0);
 
-                                    /**
-                                     * Lagrange interpolation
-                                     * matrix for different
-                                     * elements.
-                                     *
-                                     * This function builds a matrix
-                                     * $A$ such that a function
-                                     * $u_{high}$ is interpolated to
-                                     * a function of lower order
-                                     * $u_{low}$ by cell-wise
-                                     * multiplication
-                                     * $u_{low} = A u_{high}$.
-                                     */
-    static void create_interpolation_matrix(const FiniteElement<dim> &high,
-                                           const FiniteElement<dim> &low,
-                                           FullMatrix<double>& result);
-
-                                    /**
-                                     * Exception
-                                     */
-    DeclException0 (ExcInvalidFE);
                                     /**
                                      * Exception
                                      */
diff --git a/deal.II/deal.II/source/fe/fe_tools.cc b/deal.II/deal.II/source/fe/fe_tools.cc
new file mode 100644 (file)
index 0000000..f698b4f
--- /dev/null
@@ -0,0 +1,384 @@
+//----------------------------  fe_tools.cc  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 1998, 1999, 2000 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.cc  ---------------------------
+
+
+#include <base/quadrature.h>
+#include <lac/full_matrix.h>
+#include <lac/vector.h>
+#include <grid/persistent_tria.h>
+#include <grid/tria_iterator.h>
+#include <fe/fe_tools.h>
+#include <fe/fe.h>
+#include <fe/fe_values.h>
+#include <dofs/dof_handler.h>
+#include <dofs/dof_accessor.h>
+
+
+
+template <int dim, typename number>
+void FETools::get_interpolation_matrix(const FiniteElement<dim> &fe1,
+                                      const FiniteElement<dim> &fe2,
+                                      FullMatrix<number> &interpolation_matrix)
+{
+  Assert (fe1.n_components() == fe2.n_components(),
+         ExcDimensionMismatch(fe1.n_components(), fe2.n_components()));
+  Assert(interpolation_matrix.m()==fe2.dofs_per_cell &&
+        interpolation_matrix.n()==fe1.dofs_per_cell,
+        ExcMatrixDimensionMismatch(interpolation_matrix.m(),
+                                   interpolation_matrix.n(),
+                                   fe2.dofs_per_cell,
+                                   fe1.dofs_per_cell));
+
+                                  // Initialize FEValues for fe1 at
+                                  // the unit support points of the
+                                  // fe2 element.
+  vector<double> phantom_weights(fe2.dofs_per_cell,1.);
+  vector<Point<dim> > fe2_support_points (fe2.dofs_per_cell);
+  fe2.get_unit_support_points (fe2_support_points);
+  Quadrature<dim> fe2_support_points_quadrature(fe2_support_points,
+                                               phantom_weights);
+  
+  FEValues<dim> fe_values(
+    fe1, fe2_support_points_quadrature, update_values);
+
+  
+  for (unsigned int i=0; i<fe2.dofs_per_cell; ++i)     
+    {
+      const unsigned int i1
+       = fe2.system_to_component_index(i).first;
+      for (unsigned int j=0; j<fe1.dofs_per_cell; ++j)
+       {
+         const unsigned int j1
+           = fe1.system_to_component_index(j).first;
+         if (i1==j1)
+           interpolation_matrix(i,j) =
+             fe_values.shape_value (j,i);
+         else
+           interpolation_matrix(i,j)=0.;
+       }  
+    }
+  
+/*
+  vector<Point<dim> > unit_support_points (fe2.dofs_per_cell);
+  fe2.get_unit_support_points (unit_support_points);
+  
+  for (unsigned int i=0; i<fe2.dofs_per_cell; ++i)     
+    {
+      const unsigned int i1
+       = fe2.system_to_component_index(i).first;
+      for (unsigned int j=0; j<fe1.dofs_per_cell; ++j)
+       {
+         const unsigned int j1
+           = fe1.system_to_component_index(j).first;
+         if (i1==j1)
+           interpolation_matrix(i,j) =
+             fe1.shape_value (j, unit_support_points[i]);
+         else
+           interpolation_matrix(i,j)=0.;
+       }
+    }
+*/
+}
+
+
+template <int dim, typename number>
+void FETools::get_back_interpolation_matrix(const FiniteElement<dim> &fe1,
+                                           const FiniteElement<dim> &fe2,
+                                           FullMatrix<number> &interpolation_matrix)
+{
+  Assert (fe1.n_components() == fe2.n_components(),
+         ExcDimensionMismatch(fe1.n_components(), fe2.n_components()));
+  Assert(interpolation_matrix.m()==fe1.dofs_per_cell &&
+        interpolation_matrix.n()==fe1.dofs_per_cell, 
+        ExcMatrixDimensionMismatch(interpolation_matrix.m(),
+                                   interpolation_matrix.n(),
+                                   fe1.dofs_per_cell,
+                                   fe1.dofs_per_cell));
+    
+  FullMatrix<number> first_matrix (fe2.dofs_per_cell, fe1.dofs_per_cell);
+  FullMatrix<number> second_matrix(fe1.dofs_per_cell, fe2.dofs_per_cell);
+  
+  get_interpolation_matrix(fe1, fe2, first_matrix);
+  get_interpolation_matrix(fe2, fe1, second_matrix);
+
+                                  // int_matrix=second_matrix*first_matrix
+  second_matrix.mmult(interpolation_matrix, first_matrix);
+}
+
+
+
+template <int dim, typename number>
+void FETools::get_interpolation_difference_matrix(const FiniteElement<dim> &fe1,
+                                                 const FiniteElement<dim> &fe2,
+                                                 FullMatrix<number> &difference_matrix)
+{
+  Assert (fe1.n_components() == fe2.n_components(),
+         ExcDimensionMismatch(fe1.n_components(), fe2.n_components()));
+  Assert(difference_matrix.m()==fe1.dofs_per_cell &&
+        difference_matrix.n()==fe1.dofs_per_cell, 
+        ExcMatrixDimensionMismatch(difference_matrix.m(),
+                                   difference_matrix.n(),
+                                   fe1.dofs_per_cell,
+                                   fe1.dofs_per_cell));
+   
+  FullMatrix<number> interpolation_matrix(fe1.dofs_per_cell);
+  get_back_interpolation_matrix(fe1, fe2, interpolation_matrix);
+  
+  for (unsigned int i=0; i<fe1.dofs_per_cell; ++i)
+    difference_matrix(i,i) = 1.;
+  
+                                    // compute difference
+  difference_matrix.add (-1, interpolation_matrix);
+}
+
+
+
+template <int dim, typename number>
+void FETools::interpolate(const DoFHandler<dim> &dof1,
+                         const Vector<number> &u1,
+                         const DoFHandler<dim> &dof2,
+                         Vector<number> &u2)
+{
+  Assert(&dof1.get_tria()==&dof2.get_tria(), ExcTriangulationMismatch());
+  Assert(dof1.get_fe().n_components() == dof2.get_fe().n_components(),
+        ExcDimensionMismatch(dof1.get_fe().n_components(), dof2.get_fe().n_components()));
+  Assert(u1.size()==dof1.n_dofs(), ExcDimensionMismatch(u1.size(), dof1.n_dofs()));  
+  Assert(u2.size()==dof2.n_dofs(), ExcDimensionMismatch(u2.size(), dof2.n_dofs()));
+
+  const unsigned int dofs_per_cell1=dof1.get_fe().dofs_per_cell;
+  const unsigned int dofs_per_cell2=dof2.get_fe().dofs_per_cell;
+  
+  Vector<number> u1_local(dofs_per_cell1);
+  Vector<number> u2_local(dofs_per_cell2);
+
+  FullMatrix<number> interpolation_matrix(dofs_per_cell2,
+                                         dofs_per_cell1);
+  FETools::get_interpolation_matrix(dof1.get_fe(), dof2.get_fe(),
+                                   interpolation_matrix);
+  
+  DoFHandler<dim>::active_cell_iterator cell1 = dof1.begin_active(),
+                                       endc1 = dof1.end(),
+                                       cell2 = dof2.begin_active(),
+                                       endc2 = dof2.end();
+
+  vector<unsigned int> index_multiplicity(dof2.n_dofs(),0);
+  vector<unsigned int> dofs (dofs_per_cell2);
+  u2=0;
+  
+  for (; cell1!=endc1, cell2!=endc2; ++cell1, ++cell2) 
+    {
+      cell1->get_dof_values(u1, u1_local);
+      interpolation_matrix.vmult(u2_local, u1_local);
+      cell2->get_dof_indices(dofs);
+      for (unsigned int i=0; i<dofs_per_cell2; ++i)
+       {
+         u2(dofs[i])+=u2_local(i);
+         ++index_multiplicity[dofs[i]];
+       }
+    }
+  for (unsigned int i=0; i<dof2.n_dofs(); ++i)
+    {
+      Assert(index_multiplicity[i]!=0, ExcInternalError());
+      u2(i)/=static_cast<double>(index_multiplicity[i]);
+    }
+}
+
+
+
+template <int dim, typename number>
+void FETools::back_interpolate(const DoFHandler<dim> &dof1,
+                              const Vector<number> &u1,
+                              const FiniteElement<dim> &fe2,
+                              Vector<number> &u1_interpolated)
+{
+  Assert(dof1.get_fe().n_components() == fe2.n_components(),
+        ExcDimensionMismatch(dof1.get_fe().n_components(), fe2.n_components()));
+  Assert(u1.size()==dof1.n_dofs(), ExcDimensionMismatch(u1.size(), dof1.n_dofs()));  
+  Assert(u1_interpolated.size()==dof1.n_dofs(),
+        ExcDimensionMismatch(u1_interpolated.size(), dof1.n_dofs()));  
+
+  const unsigned int dofs_per_cell1=dof1.get_fe().dofs_per_cell;
+
+  Vector<number> u1_local(dofs_per_cell1);
+  Vector<number> u1_int_local(dofs_per_cell1);
+  
+  FullMatrix<number> interpolation_matrix(dofs_per_cell1, dofs_per_cell1);
+  FETools::get_back_interpolation_matrix(dof1.get_fe(), fe2,
+                                        interpolation_matrix);
+
+  DoFHandler<dim>::active_cell_iterator cell1 = dof1.begin_active(),
+                                       endc1 = dof1.end();
+  
+  for (; cell1!=endc1; ++cell1) 
+    {
+      cell1->get_dof_values(u1, u1_local);
+      interpolation_matrix.vmult(u1_int_local, u1_local);
+      cell1->set_dof_values(u1_int_local, u1_interpolated);
+    }
+}
+
+
+
+template <int dim, typename number>
+void FETools::interpolation_difference(const DoFHandler<dim> &dof1,
+                                      const Vector<number> &u1,
+                                      const FiniteElement<dim> &fe2,
+                                      Vector<number> &u1_difference)
+{
+  Assert(dof1.get_fe().n_components() == fe2.n_components(),
+        ExcDimensionMismatch(dof1.get_fe().n_components(), fe2.n_components()));
+  Assert(u1.size()==dof1.n_dofs(), ExcDimensionMismatch(u1.size(), dof1.n_dofs()));  
+  Assert(u1_difference.size()==dof1.n_dofs(),
+        ExcDimensionMismatch(u1_difference.size(), dof1.n_dofs()));  
+
+  const unsigned int dofs_per_cell1=dof1.get_fe().dofs_per_cell;
+
+  Vector<number> u1_local(dofs_per_cell1);
+  Vector<number> u1_diff_local(dofs_per_cell1);
+  
+  FullMatrix<number> difference_matrix(dofs_per_cell1, dofs_per_cell1);
+  FETools::get_interpolation_difference_matrix(dof1.get_fe(), fe2,
+                                              difference_matrix);
+  
+  DoFHandler<dim>::active_cell_iterator cell1 = dof1.begin_active(),
+                                       endc1 = dof1.end();
+  
+  for (; cell1!=endc1; ++cell1) 
+    {
+      cell1->get_dof_values(u1, u1_local);
+      difference_matrix.vmult(u1_diff_local, u1_local);
+      cell1->set_dof_values(u1_diff_local, u1_difference);
+    }
+}
+
+
+template <int dim, typename number>
+void FETools::extrapolate(const DoFHandler<dim> &dof1,
+                         const Vector<number> &u1,
+                         const DoFHandler<dim> &dof2,
+                         Vector<number> &u2)
+{
+  Assert(dof1.get_fe().n_components() == dof2.get_fe().n_components(),
+        ExcDimensionMismatch(dof1.get_fe().n_components(), dof2.get_fe().n_components()));
+  Assert(&dof1.get_tria()==&dof2.get_tria(), ExcTriangulationMismatch());
+  Assert(u1.size()==dof1.n_dofs(), ExcDimensionMismatch(u1.size(), dof1.n_dofs()));
+  Assert(u2.size()==dof2.n_dofs(), ExcDimensionMismatch(u2.size(), dof2.n_dofs()));
+
+  Vector<number> u3(dof2.n_dofs());
+  interpolate(dof1, u1, dof2, u3);
+
+  const unsigned int dofs_per_cell  = dof2.get_fe().dofs_per_cell;
+  const Triangulation<dim> &tria=dof1.get_tria();
+  Vector<number> dof_values(dofs_per_cell);
+  
+  for (unsigned int level=0; level<tria.n_levels()-1; ++level)
+    {
+      DoFHandler<dim>::cell_iterator cell=dof2.begin(level),
+                                    endc=dof2.end(level);
+
+      for (; cell!=endc; ++cell)
+       {
+         if (!cell->active())
+           {
+             bool active_children=false;
+             for (unsigned int child_n=0;
+                  child_n<GeometryInfo<dim>::children_per_cell; ++child_n)
+               if (cell->child(child_n)->active())
+                 {
+                   active_children=true;
+                   break;
+                 }
+             
+             if (active_children)
+               {
+                 cell->get_interpolated_dof_values(u3, dof_values);
+                 cell->set_dof_values_by_interpolation(dof_values, u2);
+               }
+           }
+       }
+    }  
+}
+
+
+
+
+/*-------------- Explicit Instantiations -------------------------------*/
+
+template
+void FETools::get_interpolation_matrix(const FiniteElement<deal_II_dimension> &,
+                                      const FiniteElement<deal_II_dimension> &,
+                                      FullMatrix<double> &);
+template
+void FETools::get_back_interpolation_matrix(const FiniteElement<deal_II_dimension> &,
+                                           const FiniteElement<deal_II_dimension> &,
+                                           FullMatrix<double> &);
+template
+void FETools::get_interpolation_difference_matrix(const FiniteElement<deal_II_dimension> &,
+                                                 const FiniteElement<deal_II_dimension> &,
+                                                 FullMatrix<double> &);
+template
+void FETools::interpolate(const DoFHandler<deal_II_dimension> &,
+                         const Vector<double> &,
+                         const DoFHandler<deal_II_dimension> &,
+                         Vector<double> &);
+template
+void FETools::back_interpolate(const DoFHandler<deal_II_dimension> &,
+                              const Vector<double> &,
+                              const FiniteElement<deal_II_dimension> &,
+                              Vector<double> &);
+template
+void FETools::interpolation_difference(const DoFHandler<deal_II_dimension> &,
+                                      const Vector<double> &,
+                                      const FiniteElement<deal_II_dimension> &,
+                                      Vector<double> &);
+template
+void FETools::extrapolate(const DoFHandler<deal_II_dimension> &,
+                         const Vector<double> &,
+                         const DoFHandler<deal_II_dimension> &,
+                         Vector<double> &);
+
+
+template
+void FETools::get_interpolation_matrix(const FiniteElement<deal_II_dimension> &,
+                                      const FiniteElement<deal_II_dimension> &,
+                                      FullMatrix<float> &);
+template
+void FETools::get_back_interpolation_matrix(const FiniteElement<deal_II_dimension> &,
+                                           const FiniteElement<deal_II_dimension> &,
+                                           FullMatrix<float> &);
+template
+void FETools::get_interpolation_difference_matrix(const FiniteElement<deal_II_dimension> &,
+                                                 const FiniteElement<deal_II_dimension> &,
+                                                 FullMatrix<float> &);
+template
+void FETools::interpolate(const DoFHandler<deal_II_dimension> &,
+                         const Vector<float> &,
+                         const DoFHandler<deal_II_dimension> &,
+                         Vector<float> &);
+template
+void FETools::back_interpolate(const DoFHandler<deal_II_dimension> &,
+                              const Vector<float> &,
+                              const FiniteElement<deal_II_dimension> &,
+                              Vector<float> &);
+template
+void FETools::interpolation_difference(const DoFHandler<deal_II_dimension> &,
+                                      const Vector<float> &,
+                                      const FiniteElement<deal_II_dimension> &,
+                                      Vector<float> &);
+template
+void FETools::extrapolate(const DoFHandler<deal_II_dimension> &,
+                         const Vector<float> &,
+                         const DoFHandler<deal_II_dimension> &,
+                         Vector<float> &);
+
+/*----------------------------   fe_tools.cc     ---------------------------*/
index aeebcbefac1d2311f3eb63ee93bea94ff4d60bb4..6a9648598b3598d7e90ff5a04fae3f6635f03839 100644 (file)
@@ -1368,42 +1368,6 @@ void LaplaceMatrix<dim>::assemble (Vector<double>      &rhs,
 
 
 
-template<int dim>
-void
-MatrixCreator<dim>::create_interpolation_matrix(const FiniteElement<dim> &high,
-                                               const FiniteElement<dim> &low,
-                                               FullMatrix<double>& result)
-{
-  Assert (high.n_components() == low.n_components(),
-         ExcInvalidFE());
-  
-  Assert (result.m() == low.dofs_per_cell,
-         ExcDimensionMismatch(result.m(), low.dofs_per_cell));
-  Assert (result.n() == high.dofs_per_cell,
-         ExcDimensionMismatch(result.n(), high.dofs_per_cell));
-
-
-                                  // Initialize FEValues at the support points
-                                  // of the low element.
-  vector<double> phantom_weights(low.dofs_per_cell,1.);
-  vector<Point<dim> > support_points(low.dofs_per_cell);
-  low.get_unit_support_points(support_points);
-  Quadrature<dim> low_points(support_points,
-                            phantom_weights);
-
-  FEValues<dim> fe(high, low_points, update_values);
-  
-  for (unsigned int i=0; i<low.dofs_per_cell; ++i)
-    for (unsigned int j=0; j<high.dofs_per_cell; ++j)
-                                      // shape functions need to belong
-                                      // to the same component
-      if (low.system_to_component_index(i).first ==
-         high.system_to_component_index(j).first)
-       result(i,j) = fe.shape_value(j,i);
-}
-
-
-
 // explicit instantiations
 
 template class MatrixCreator<deal_II_dimension>;

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.