]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Templatize on vector types in many functions, to support block vectors and Vector...
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 16 May 2000 09:01:16 +0000 (09:01 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 16 May 2000 09:01:16 +0000 (09:01 +0000)
git-svn-id: https://svn.dealii.org/trunk@2861 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/dofs/dof_accessor.h
deal.II/deal.II/include/fe/fe_values.h
deal.II/deal.II/source/dofs/dof_accessor.cc
deal.II/deal.II/source/fe/fe_values.cc
deal.II/deal.II/source/fe/q1_mapping.jacobians.cc

index dc68f03a7958359feb8b33c370d396d65dccdd94..574117bab09fff6490a58bc9512fadb197838657 100644 (file)
@@ -300,12 +300,23 @@ class DoFObjectAccessor : public DoFAccessor<dim>,
                                      * The vector has to have the
                                      * right size before being passed
                                      * to this function. This
-                                     * function is only callable for active
-                                     * cells.
-                                     */
-    template <typename number>
-    void get_dof_values (const Vector<number> &values,
-                        Vector<number>       &local_values) const;
+                                     * function is only callable for
+                                     * active cells.
+                                     *
+                                     * The input vector may be either
+                                     * a #Vector<float>#,
+                                     * #Vector<double>#, or a
+                                     * #BlockVector<...,double>#. It
+                                     * is in the responsibility of
+                                     * the caller to assure that the
+                                     * types of the numbers stored in
+                                     * input and output vectors are
+                                     * compatible and with similar
+                                     * accuracy.
+                                     */
+    template <class InputVector, typename number>
+    void get_dof_values (const InputVector &values,
+                        Vector<number>    &local_values) const;
 
                                     /**
                                      * This function is the counterpart to
@@ -329,10 +340,21 @@ class DoFObjectAccessor : public DoFAccessor<dim>,
                                      * The vector has to have the
                                      * right size before being passed
                                      * to this function.
-                                     */
-    template <typename number>
+                                     *
+                                     * The output vector may be either
+                                     * a #Vector<float>#,
+                                     * #Vector<double>#, or a
+                                     * #BlockVector<...,double>#. It
+                                     * is in the responsibility of
+                                     * the caller to assure that the
+                                     * types of the numbers stored in
+                                     * input and output vectors are
+                                     * compatible and with similar
+                                     * accuracy.
+                                     */
+    template <class OutputVector, typename number>
     void set_dof_values (const Vector<number> &local_values,
-                        Vector<number>       &values) const;
+                        OutputVector         &values) const;
 
                                     /**
                                      *  Pointer to the #i#th line
@@ -531,10 +553,21 @@ class DoFObjectAccessor<1, dim> :  public DoFAccessor<dim>,
                                      * has the right size beforehand. This
                                      * function is only callable for active
                                      * cells.
-                                     */
-    template <typename number>
-    void get_dof_values (const Vector<number> &values,
-                        Vector<number>       &local_values) const;
+                                     *
+                                     * The input vector may be either
+                                     * a #Vector<float>#,
+                                     * #Vector<double>#, or a
+                                     * #BlockVector<...,double>#. It
+                                     * is in the responsibility of
+                                     * the caller to assure that the
+                                     * types of the numbers stored in
+                                     * input and output vectors are
+                                     * compatible and with similar
+                                     * accuracy.
+                                     */
+    template <class InputVector, typename number>
+    void get_dof_values (const InputVector &values,
+                        Vector<number>    &local_values) const;
 
                                     /**
                                      * This function is the counterpart to
@@ -557,10 +590,21 @@ class DoFObjectAccessor<1, dim> :  public DoFAccessor<dim>,
                                      *
                                      * It is assumed that both vectors already
                                      * have the right size beforehand.
-                                     */
-    template <typename number>
+                                     *
+                                     * The output vector may be either
+                                     * a #Vector<float>#,
+                                     * #Vector<double>#, or a
+                                     * #BlockVector<...,double>#. It
+                                     * is in the responsibility of
+                                     * the caller to assure that the
+                                     * types of the numbers stored in
+                                     * input and output vectors are
+                                     * compatible and with similar
+                                     * accuracy.
+                                     */
+    template <class OutputVector, typename number>
     void set_dof_values (const Vector<number> &local_values,
-                        Vector<number>       &values) const;
+                        OutputVector         &values) const;
 
                                     /**
                                      * Return the #i#th child as a DoF line
@@ -698,10 +742,21 @@ class DoFObjectAccessor<2, dim> :  public DoFAccessor<dim>,
                                      * has the right size beforehand. This
                                      * function is only callable for active
                                      * cells.
-                                     */
-    template <typename number>
-    void get_dof_values (const Vector<number> &values,
-                        Vector<number>       &local_values) const;
+                                     *
+                                     * The input vector may be either
+                                     * a #Vector<float>#,
+                                     * #Vector<double>#, or a
+                                     * #BlockVector<...,double>#. It
+                                     * is in the responsibility of
+                                     * the caller to assure that the
+                                     * types of the numbers stored in
+                                     * input and output vectors are
+                                     * compatible and with similar
+                                     * accuracy.
+                                     */
+    template <class InputVector, typename number>
+    void get_dof_values (const InputVector &values,
+                        Vector<number>    &local_values) const;
 
                                     /**
                                      * This function is the counterpart to
@@ -724,10 +779,21 @@ class DoFObjectAccessor<2, dim> :  public DoFAccessor<dim>,
                                      *
                                      * It is assumed that both vectors already
                                      * have the right size beforehand.
-                                     */
-    template <typename number>
+                                     *
+                                     * The output vector may be either
+                                     * a #Vector<float>#,
+                                     * #Vector<double>#, or a
+                                     * #BlockVector<...,double>#. It
+                                     * is in the responsibility of
+                                     * the caller to assure that the
+                                     * types of the numbers stored in
+                                     * input and output vectors are
+                                     * compatible and with similar
+                                     * accuracy.
+                                     */
+    template <class OutputVector, typename number>
     void set_dof_values (const Vector<number> &local_values,
-                        Vector<number>       &values) const;
+                        OutputVector         &values) const;
 
                                     /**
                                      *  Return a pointer to the #i#th line
@@ -873,10 +939,21 @@ class DoFObjectAccessor<3, dim> :  public DoFAccessor<dim>,
                                      * has the right size beforehand. This
                                      * function is only callable for active
                                      * cells.
-                                     */
-    template <typename number>
-    void get_dof_values (const Vector<number> &values,
-                        Vector<number>       &local_values) const;
+                                     *
+                                     * The input vector may be either
+                                     * a #Vector<float>#,
+                                     * #Vector<double>#, or a
+                                     * #BlockVector<...,double>#. It
+                                     * is in the responsibility of
+                                     * the caller to assure that the
+                                     * types of the numbers stored in
+                                     * input and output vectors are
+                                     * compatible and with similar
+                                     * accuracy.
+                                     */
+    template <class InputVector, typename number>
+    void get_dof_values (const InputVector &values,
+                        Vector<number>    &local_values) const;
 
                                     /**
                                      * This function is the counterpart to
@@ -899,10 +976,21 @@ class DoFObjectAccessor<3, dim> :  public DoFAccessor<dim>,
                                      *
                                      * It is assumed that both vectors already
                                      * have the right size beforehand.
-                                     */
-    template <typename number>
+                                     *
+                                     * The output vector may be either
+                                     * a #Vector<float>#,
+                                     * #Vector<double>#, or a
+                                     * #BlockVector<...,double>#. It
+                                     * is in the responsibility of
+                                     * the caller to assure that the
+                                     * types of the numbers stored in
+                                     * input and output vectors are
+                                     * compatible and with similar
+                                     * accuracy.
+                                     */
+    template <class OutputVector, typename number>
     void set_dof_values (const Vector<number> &local_values,
-                        Vector<number>       &values) const;
+                        OutputVector         &values) const;
 
                                     /**
                                      *  Return a pointer to the #i#th line
@@ -1074,9 +1162,9 @@ class DoFCellAccessor :  public DoFObjectAccessor<dim, dim>
                                      * cells by the finite element
                                      * objects.
                                      */
-    template <typename number>
-    void get_interpolated_dof_values (const Vector<number> &values,
-                                     Vector<number>       &interpolated_values) const;
+    template <class InputVector, typename number>
+    void get_interpolated_dof_values (const InputVector &values,
+                                     Vector<number>    &interpolated_values) const;
 
                                     /**
                                      * This, again, is the counterpart to
@@ -1139,10 +1227,21 @@ class DoFCellAccessor :  public DoFObjectAccessor<dim, dim>
                                      * presently only provided for
                                      * cells by the finite element
                                      * objects.
-                                     */
-    template <typename number>
+                                     *
+                                     * The output vector may be either
+                                     * a #Vector<float>#,
+                                     * #Vector<double>#, or a
+                                     * #BlockVector<...,double>#. It
+                                     * is in the responsibility of
+                                     * the caller to assure that the
+                                     * types of the numbers stored in
+                                     * input and output vectors are
+                                     * compatible and with similar
+                                     * accuracy.
+                                     */
+    template <class OutputVector, typename number>
     void set_dof_values_by_interpolation (const Vector<number> &local_values,
-                                         Vector<number>       &values) const;
+                                         OutputVector         &values) const;
     
                                     /**
                                      *  Exception
index 2b000be90ce6062d65a1657063dbd2a201357eff..b3ec0a389b582205cccd01cc915680ce6bdbebdd 100644 (file)
@@ -303,20 +303,35 @@ class FEValuesBase
                                      *
                                      * The function assumes that the
                                      * #values# object already has the
-                                     * right size. 
+                                     * right size.
+                                     *
+                                     * The actual data type of the
+                                     * input vector may be either a
+                                     * #Vector<double>#,
+                                     * #Vector<float>#, or
+                                     * #BlockVector<double,...>#.
                                      */
-    void get_function_values (const Vector<double> &fe_function,
-                             vector<double>       &values) const;
+    template <class InputVector>
+    void get_function_values (const InputVector &fe_function,
+                             vector<double>    &values) const;
 
                                     /**
-                                     * Access to vector valued finite element functions.
+                                     * Access to vector valued finite
+                                     * element functions.
                                      *
                                      * This function does the same as
                                      * the other #get_function_values#,
                                      * but applied to multi-component
                                      * elements.
-                                     */
-    void get_function_values (const Vector<double>    &fe_function,
+                                     *
+                                     * The actual data type of the
+                                     * input vector may be either a
+                                     * #Vector<double>#,
+                                     * #Vector<float>#, or
+                                     * #BlockVector<double,...>#.
+                                     */
+    template <class InputVector>
+    void get_function_values (const InputVector       &fe_function,
                              vector<Vector<double> > &values) const;
 
                                     /**
@@ -357,8 +372,15 @@ class FEValuesBase
                                      * The function assumes that the
                                      * #gradients# object already has the
                                      * right size.
-                                     */
-    void get_function_grads (const Vector<double>   &fe_function,
+                                     *
+                                     * The actual data type of the
+                                     * input vector may be either a
+                                     * #Vector<double>#,
+                                     * #Vector<float>#, or
+                                     * #BlockVector<double,...>#.
+                                     */
+    template <class InputVector>
+    void get_function_grads (const InputVector      &fe_function,
                             vector<Tensor<1,dim> > &gradients) const;
 
                                     /**
@@ -379,8 +401,15 @@ class FEValuesBase
                                      * the other #get_function_values#,
                                      * but applied to multi-component
                                      * elements.
-                                     */
-    void get_function_grads (const Vector<double>   &fe_function,
+                                     *
+                                     * The actual data type of the
+                                     * input vector may be either a
+                                     * #Vector<double>#,
+                                     * #Vector<float>#, or
+                                     * #BlockVector<double,...>#.
+                                     */
+    template <class InputVector>
+    void get_function_grads (const InputVector               &fe_function,
                             vector<vector<Tensor<1,dim> > > &gradients) const;
 
                                     /**
@@ -426,8 +455,15 @@ class FEValuesBase
                                      * The function assumes that the
                                      * #second_derivatives# object already has
                                      * the right size.
-                                     */
-    void get_function_2nd_derivatives (const Vector<double>   &fe_function,
+                                     *
+                                     * The actual data type of the
+                                     * input vector may be either a
+                                     * #Vector<double>#,
+                                     * #Vector<float>#, or
+                                     * #BlockVector<double,...>#.
+                                     */
+    template <class InputVector>
+    void get_function_2nd_derivatives (const InputVector      &fe_function,
                                       vector<Tensor<2,dim> > &second_derivatives) const;
 
                                     /**
index 7989c1cc92653623adab6cafad846fe0d93ca3c7..c6c9f088ab116fa8e553d7a304b0a84a3c12c529 100644 (file)
@@ -21,6 +21,7 @@
 #include <fe/fe.h>
 
 #include <lac/vector.h>
+#include <lac/block_vector.h>
 #include <lac/sparse_matrix.h>
 
 #include <vector>
@@ -117,10 +118,10 @@ distribute_local_to_global (const FullMatrix<double> &local_source,
 
 
 template <int dim>
-template <typename number>
+template <class InputVector, typename number>
 void
-DoFObjectAccessor<1,dim>::get_dof_values (const Vector<number> &values,
-                                         Vector<number>       &local_values) const
+DoFObjectAccessor<1,dim>::get_dof_values (const InputVector &values,
+                                         Vector<number>    &local_values) const
 {
   Assert (dim==1, ExcInternalError());
   
@@ -147,10 +148,11 @@ DoFObjectAccessor<1,dim>::get_dof_values (const Vector<number> &values,
 
 
 template <int dim>
-template <typename number>
+template <class OutputVector, typename number>
 void
 DoFObjectAccessor<1,dim>::set_dof_values (const Vector<number> &local_values,
-                                         Vector<number>       &values) const {
+                                         OutputVector         &values) const
+{
   Assert (dim==1, ExcInternalError());
   
   Assert (dof_handler != 0, DoFAccessor<1>::ExcInvalidObject());
@@ -269,10 +271,10 @@ distribute_local_to_global (const FullMatrix<double> &local_source,
 
 
 template <int dim>
-template <typename number>
+template <class InputVector, typename number>
 void
-DoFObjectAccessor<2,dim>::get_dof_values (const Vector<number> &values,
-                                         Vector<number>       &local_values) const
+DoFObjectAccessor<2,dim>::get_dof_values (const InputVector &values,
+                                         Vector<number>    &local_values) const
 {
   Assert (dim==2, ExcInternalError());
   
@@ -303,10 +305,10 @@ DoFObjectAccessor<2,dim>::get_dof_values (const Vector<number> &values,
 
 
 template <int dim>
-template <typename number>
+template <class OutputVector, typename number>
 void
 DoFObjectAccessor<2,dim>::set_dof_values (const Vector<number> &local_values,
-                                         Vector<number>       &values) const
+                                         OutputVector         &values) const
 {
   Assert (dim==2, ExcInternalError());
   
@@ -432,10 +434,10 @@ distribute_local_to_global (const FullMatrix<double> &local_source,
 
 
 template <int dim>
-template <typename number>
+template <class InputVector, typename number>
 void
-DoFObjectAccessor<3,dim>::get_dof_values (const Vector<number> &values,
-                                         Vector<number>       &local_values) const
+DoFObjectAccessor<3,dim>::get_dof_values (const InputVector &values,
+                                         Vector<number>    &local_values) const
 {
   Assert (dim==3, ExcInternalError());
 
@@ -470,10 +472,10 @@ DoFObjectAccessor<3,dim>::get_dof_values (const Vector<number> &values,
 
 
 template <int dim>
-template <typename number>
+template <class OutputVector, typename number>
 void
 DoFObjectAccessor<3,dim>::set_dof_values (const Vector<number> &local_values,
-                                         Vector<number>       &values) const
+                                         OutputVector         &values) const
 {
   Assert (dim==3, ExcInternalError());
 
@@ -548,10 +550,10 @@ DoFCellAccessor<3>::face (const unsigned int i) const
 
 
 template <int dim>
-template <typename number>
+template <class InputVector, typename number>
 void
-DoFCellAccessor<dim>::get_interpolated_dof_values (const Vector<number> &values,
-                                                  Vector<number>       &interpolated_values) const
+DoFCellAccessor<dim>::get_interpolated_dof_values (const InputVector &values,
+                                                  Vector<number>    &interpolated_values) const
 {
   const unsigned int dofs_per_cell = dof_handler->get_fe().dofs_per_cell;
   
@@ -610,10 +612,10 @@ DoFCellAccessor<dim>::get_interpolated_dof_values (const Vector<number> &values,
 
 
 template <int dim>
-template <typename number>
+template <class OutputVector, typename number>
 void
 DoFCellAccessor<dim>::set_dof_values_by_interpolation (const Vector<number> &local_values,
-                                                      Vector<number>       &values) const
+                                                      OutputVector         &values) const
 {
   const unsigned int dofs_per_cell = dof_handler->get_fe().dofs_per_cell;
   
@@ -645,63 +647,109 @@ DoFCellAccessor<dim>::set_dof_values_by_interpolation (const Vector<number> &loc
 };
 
 
+
+// --------------------------------------------------------------------------
 // explicit instantiations
 
 
 // for double
 template
 void
-DoFObjectAccessor<1,deal_II_dimension>::
-get_dof_values (const Vector<double> &,
-               Vector<double>       &) const;
+DoFObjectAccessor<1,deal_II_dimension>::get_dof_values (const Vector<double> &,
+                                                       Vector<double>       &) const;
 
 template
 void
-DoFObjectAccessor<1,deal_II_dimension>::
-set_dof_values (const Vector<double> &,
-               Vector<double>       &) const;
+DoFObjectAccessor<1,deal_II_dimension>::set_dof_values (const Vector<double> &,
+                                                       Vector<double>       &) const;
 
 
 // for float
 template
 void
-DoFObjectAccessor<1,deal_II_dimension>::
-get_dof_values (const Vector<float> &,
-               Vector<float>       &) const;
+DoFObjectAccessor<1,deal_II_dimension>::get_dof_values (const Vector<float> &,
+                                                       Vector<float>       &) const;
+
+template
+void
+DoFObjectAccessor<1,deal_II_dimension>::set_dof_values (const Vector<float> &,
+                                                       Vector<float>       &) const;
+
+
+// for block vector 2
+template
+void
+DoFObjectAccessor<1,deal_II_dimension>::get_dof_values (const BlockVector<2,double> &,
+                                                       Vector<float>       &) const;
+
+template
+void
+DoFObjectAccessor<1,deal_II_dimension>::set_dof_values (const Vector<double>  &,
+                                                       BlockVector<2,double> &) const;
+
+
+// for block vector 3
+template
+void
+DoFObjectAccessor<1,deal_II_dimension>::get_dof_values (const BlockVector<3,double> &,
+                                                       Vector<float>       &) const;
 
 template
 void
-DoFObjectAccessor<1,deal_II_dimension>::
-set_dof_values (const Vector<float> &,
-               Vector<float>       &) const;
+DoFObjectAccessor<1,deal_II_dimension>::set_dof_values (const Vector<double>  &,
+                                                       BlockVector<3,double> &) const;
+
+
+
 
 
 #if deal_II_dimension >= 2
 // for double
 template
 void
-DoFObjectAccessor<2,deal_II_dimension>::
-get_dof_values (const Vector<double> &,
-               Vector<double>       &) const;
+DoFObjectAccessor<2,deal_II_dimension>::get_dof_values (const Vector<double> &,
+                                                       Vector<double>       &) const;
 
 template
 void
-DoFObjectAccessor<2,deal_II_dimension>::
-set_dof_values (const Vector<double> &,
-               Vector<double>       &) const;
+DoFObjectAccessor<2,deal_II_dimension>::set_dof_values (const Vector<double> &,
+                                                       Vector<double>       &) const;
 
 // for float
 template
 void
-DoFObjectAccessor<2,deal_II_dimension>::
-get_dof_values (const Vector<float> &,
-               Vector<float>       &) const;
+DoFObjectAccessor<2,deal_II_dimension>::get_dof_values (const Vector<float> &,
+                                                       Vector<float>       &) const;
+
+template
+void
+DoFObjectAccessor<2,deal_II_dimension>::set_dof_values (const Vector<float> &,
+                                                       Vector<float>       &) const;
 
+
+// for block vector 2
 template
 void
-DoFObjectAccessor<2,deal_II_dimension>::
-set_dof_values (const Vector<float> &,
-               Vector<float>       &) const;
+DoFObjectAccessor<2,deal_II_dimension>::get_dof_values (const BlockVector<2,double> &,
+                                                       Vector<float>       &) const;
+
+template
+void
+DoFObjectAccessor<2,deal_II_dimension>::set_dof_values (const Vector<double>  &,
+                                                       BlockVector<2,double> &) const;
+
+
+// for block vector 3
+template
+void
+DoFObjectAccessor<2,deal_II_dimension>::get_dof_values (const BlockVector<3,double> &,
+                                                       Vector<float>       &) const;
+
+template
+void
+DoFObjectAccessor<2,deal_II_dimension>::set_dof_values (const Vector<double>  &,
+                                                       BlockVector<3,double> &) const;
+
 
 #endif
 
@@ -710,32 +758,54 @@ set_dof_values (const Vector<float> &,
 // for double
 template
 void
-DoFObjectAccessor<3,deal_II_dimension>::
-get_dof_values (const Vector<double> &,
-               Vector<double>       &) const;
+DoFObjectAccessor<3,deal_II_dimension>::get_dof_values (const Vector<double> &,
+                                                       Vector<double>       &) const;
 
 template
 void
-DoFObjectAccessor<3,deal_II_dimension>::
-set_dof_values (const Vector<double> &,
-               Vector<double>       &) const;
+DoFObjectAccessor<3,deal_II_dimension>::set_dof_values (const Vector<double> &,
+                                                       Vector<double>       &) const;
 
 // for float
 template
 void
-DoFObjectAccessor<3,deal_II_dimension>::
-get_dof_values (const Vector<float> &,
-               Vector<float>       &) const;
+DoFObjectAccessor<3,deal_II_dimension>::get_dof_values (const Vector<float> &,
+                                                       Vector<float>       &) const;
+
+template
+void
+DoFObjectAccessor<3,deal_II_dimension>::set_dof_values (const Vector<float> &,
+                                                       Vector<float>       &) const;
+
+
 
+// for block vector 2
+template
+void
+DoFObjectAccessor<3,deal_II_dimension>::get_dof_values (const BlockVector<2,double> &,
+                                                       Vector<float>       &) const;
+template
+void
+DoFObjectAccessor<3,deal_II_dimension>::set_dof_values (const Vector<double>  &,
+                                                       BlockVector<2,double> &) const;
+
+
+// for block vector 3
+template
+void
+DoFObjectAccessor<3,deal_II_dimension>::get_dof_values (const BlockVector<3,double> &,
+                                                       Vector<float>       &) const;
 template
 void
-DoFObjectAccessor<3,deal_II_dimension>::
-set_dof_values (const Vector<float> &,
-               Vector<float>       &) const;
+DoFObjectAccessor<3,deal_II_dimension>::set_dof_values (const Vector<double>  &,
+                                                       BlockVector<3,double> &) const;
+
+
 
 #endif
 
 
+
 template
 void
 DoFCellAccessor<deal_II_dimension>::
@@ -762,6 +832,19 @@ set_dof_values_by_interpolation(const Vector<float> &,
                                Vector<float>       &) const;
 
 
+template
+void
+DoFCellAccessor<deal_II_dimension>::
+get_interpolated_dof_values (const BlockVector<2,double> &,
+                            Vector<double>       &) const;
+
+template
+void
+DoFCellAccessor<deal_II_dimension>::
+set_dof_values_by_interpolation(const Vector<double> &,
+                               BlockVector<3,double> &) const;
+
+
 #if deal_II_dimension == 1
 template class DoFObjectAccessor<1, 1>;
 #endif
index b0ff1455be2e9c702cd7817f4105a849c75a71b1..ee631100d3a636426382ba068272260db6790646 100644 (file)
@@ -20,6 +20,7 @@
 #include <grid/tria_boundary.h>
 #include <dofs/dof_accessor.h>
 #include <lac/vector.h>
+#include <lac/block_vector.h>
 
 
 /*------------------------------- FEValuesBase ---------------------------*/
@@ -55,9 +56,11 @@ FEValuesBase<dim>::FEValuesBase (const unsigned int n_q_points,
 {};
 
 
+
 template <int dim>
 double FEValuesBase<dim>::shape_value (const unsigned int i,
-                                      const unsigned int j) const {
+                                      const unsigned int j) const
+{
   Assert (update_flags & update_values, ExcAccessToUninitializedField());
   Assert (selected_dataset<shape_values.size(),
          ExcIndexRange (selected_dataset, 0, shape_values.size()));
@@ -70,9 +73,11 @@ double FEValuesBase<dim>::shape_value (const unsigned int i,
 };
 
 
+
 template <int dim>
-void FEValuesBase<dim>::get_function_values (const Vector<double> &fe_function,
-                                            vector<double>       &values) const
+template <class InputVector>
+void FEValuesBase<dim>::get_function_values (const InputVector &fe_function,
+                                            vector<double>    &values) const
 {
   Assert (fe->n_components() == 1,
          ExcWrongNoOfComponents());
@@ -101,8 +106,10 @@ void FEValuesBase<dim>::get_function_values (const Vector<double> &fe_function,
 };
 
 
+
 template <int dim>
-void FEValuesBase<dim>::get_function_values (const Vector<double>     &fe_function,
+template <class InputVector>
+void FEValuesBase<dim>::get_function_values (const InputVector        &fe_function,
                                             vector< Vector<double> > &values) const
 {
   Assert (n_quadrature_points == values.size(),
@@ -134,6 +141,7 @@ void FEValuesBase<dim>::get_function_values (const Vector<double>     &fe_functi
 };
 
 
+
 template <int dim>
 const Tensor<1,dim> &
 FEValuesBase<dim>::shape_grad (const unsigned int i,
@@ -149,8 +157,10 @@ FEValuesBase<dim>::shape_grad (const unsigned int i,
 };
 
 
+
 template <int dim>
-void FEValuesBase<dim>::get_function_grads (const Vector<double>   &fe_function,
+template <class InputVector>
+void FEValuesBase<dim>::get_function_grads (const InputVector      &fe_function,
                                            vector<Tensor<1,dim> > &gradients) const
 {
   Assert (fe->n_components() == 1,
@@ -180,8 +190,11 @@ void FEValuesBase<dim>::get_function_grads (const Vector<double>   &fe_function,
       };
 };
 
+
+
 template <int dim>
-void FEValuesBase<dim>::get_function_grads (const Vector<double>            &fe_function,
+template <class InputVector>
+void FEValuesBase<dim>::get_function_grads (const InputVector               &fe_function,
                                            vector<vector<Tensor<1,dim> > > &gradients) const
 {
   Assert (n_quadrature_points == gradients.size(),
@@ -217,10 +230,12 @@ void FEValuesBase<dim>::get_function_grads (const Vector<double>            &fe_
 };
 
 
+
 template <int dim>
 const Tensor<2,dim> &
 FEValuesBase<dim>::shape_2nd_derivative (const unsigned int i,
-                                        const unsigned int j) const {
+                                        const unsigned int j) const
+{
   Assert (i<shape_2nd_derivatives.size(),
          ExcIndexRange (i, 0, shape_2nd_derivatives.size()));
   Assert (j<shape_2nd_derivatives[i].size(),
@@ -231,8 +246,10 @@ FEValuesBase<dim>::shape_2nd_derivative (const unsigned int i,
 };
 
 
+
 template <int dim>
-void FEValuesBase<dim>::get_function_2nd_derivatives (const Vector<double>   &fe_function,
+template <class InputVector>
+void FEValuesBase<dim>::get_function_2nd_derivatives (const InputVector      &fe_function,
                                                      vector<Tensor<2,dim> > &second_derivatives) const
 {
   Assert (fe->n_components() == 1,
@@ -263,6 +280,7 @@ void FEValuesBase<dim>::get_function_2nd_derivatives (const Vector<double>   &fe
 };
 
 
+
 template <int dim>
 const Point<dim> &
 FEValuesBase<dim>::quadrature_point (const unsigned int i) const
@@ -274,6 +292,7 @@ FEValuesBase<dim>::quadrature_point (const unsigned int i) const
 };
 
 
+
 template <int dim>
 const Point<dim> &
 FEValuesBase<dim>::support_point (const unsigned int i) const
@@ -285,6 +304,7 @@ FEValuesBase<dim>::support_point (const unsigned int i) const
 };
 
 
+
 template <int dim>
 double FEValuesBase<dim>::JxW (const unsigned int i) const
 {
@@ -343,8 +363,10 @@ FEValues<dim>::FEValues (const FiniteElement<dim> &fe,
 };
 
 
+
 template <int dim>
-void FEValues<dim>::reinit (const typename DoFHandler<dim>::cell_iterator &cell) {
+void FEValues<dim>::reinit (const typename DoFHandler<dim>::cell_iterator &cell)
+{
   present_cell = cell;
                                   // fill jacobi matrices and real
                                   // quadrature points
@@ -458,9 +480,11 @@ FEFaceValuesBase<dim>::FEFaceValuesBase (const unsigned int n_q_points,
 {};
 
 
+
 template <int dim>
 const Point<dim> &
-FEFaceValuesBase<dim>::normal_vector (const unsigned int i) const {
+FEFaceValuesBase<dim>::normal_vector (const unsigned int i) const
+{
   Assert (i<normal_vectors.size(), ExcIndexRange(i, 0, normal_vectors.size()));
   Assert (update_flags & update_normal_vectors,
          ExcAccessToUninitializedField());
@@ -522,7 +546,8 @@ FEFaceValues<dim>::FEFaceValues (const FiniteElement<dim> &fe,
 
 template <int dim>
 void FEFaceValues<dim>::reinit (const typename DoFHandler<dim>::cell_iterator &cell,
-                               const unsigned int                             face_no) {
+                               const unsigned int                             face_no)
+{
   present_cell  = cell;
   selected_dataset = face_no;
                                   // fill jacobi matrices and real
@@ -686,6 +711,7 @@ FESubfaceValues<dim>::FESubfaceValues (const FiniteElement<dim> &fe,
 };
 
 
+
 template <int dim>
 void FESubfaceValues<dim>::reinit (const typename DoFHandler<dim>::cell_iterator &cell,
                                   const unsigned int         face_no,
@@ -794,3 +820,79 @@ template class FEFaceValuesBase<deal_II_dimension>;
 template class FEFaceValues<deal_II_dimension>;
 template class FESubfaceValues<deal_II_dimension>;
 #endif
+
+
+//-----------------------------------------------------------------------------
+
+template
+void FEValuesBase<deal_II_dimension>::get_function_values (const Vector<double> &,
+                                            vector<double>       &) const;
+template
+void FEValuesBase<deal_II_dimension>::get_function_values (const Vector<float> &,
+                                            vector<double>      &) const;
+template
+void FEValuesBase<deal_II_dimension>::get_function_values (const BlockVector<2,double> &,
+                                            vector<double>      &) const;
+template
+void FEValuesBase<deal_II_dimension>::get_function_values (const BlockVector<3,double> &,
+                                            vector<double>      &) const;
+
+//-----------------------------------------------------------------------------
+
+template
+void FEValuesBase<deal_II_dimension>::get_function_values (const Vector<double> &,
+                                            vector<Vector<double> > &) const;
+template
+void FEValuesBase<deal_II_dimension>::get_function_values (const Vector<float> &,
+                                            vector<Vector<double> > &) const;
+template
+void FEValuesBase<deal_II_dimension>::get_function_values (const BlockVector<2,double> &,
+                                            vector<Vector<double> >     &) const;
+template
+void FEValuesBase<deal_II_dimension>::get_function_values (const BlockVector<3,double> &,
+                                            vector<Vector<double> >     &) const;
+
+//-----------------------------------------------------------------------------
+
+template
+void FEValuesBase<deal_II_dimension>::get_function_grads (const Vector<double> &,
+                                           vector<Tensor<1,deal_II_dimension> > &) const;
+template
+void FEValuesBase<deal_II_dimension>::get_function_grads (const Vector<float> &,
+                                            vector<Tensor<1,deal_II_dimension> > &) const;
+template
+void FEValuesBase<deal_II_dimension>::get_function_grads (const BlockVector<2,double> &,
+                                            vector<Tensor<1,deal_II_dimension> > &) const;
+template
+void FEValuesBase<deal_II_dimension>::get_function_grads (const BlockVector<3,double> &,
+                                            vector<Tensor<1,deal_II_dimension> > &) const;
+
+//-----------------------------------------------------------------------------
+
+template
+void FEValuesBase<deal_II_dimension>::get_function_grads (const Vector<double> &,
+                                            vector<vector<Tensor<1,deal_II_dimension> > > &) const;
+template
+void FEValuesBase<deal_II_dimension>::get_function_grads (const Vector<float> &,
+                                            vector<vector<Tensor<1,deal_II_dimension> > > &) const;
+template
+void FEValuesBase<deal_II_dimension>::get_function_grads (const BlockVector<2,double> &,
+                                            vector<vector<Tensor<1,deal_II_dimension> > > &) const;
+template
+void FEValuesBase<deal_II_dimension>::get_function_grads (const BlockVector<3,double> &,
+                                            vector<vector<Tensor<1,deal_II_dimension> > > &) const;
+
+//-----------------------------------------------------------------------------
+
+template
+void FEValuesBase<deal_II_dimension>::get_function_2nd_derivatives (const Vector<double> &,
+                                           vector<Tensor<2,deal_II_dimension> > &) const;
+template
+void FEValuesBase<deal_II_dimension>::get_function_2nd_derivatives (const Vector<float> &,
+                                            vector<Tensor<2,deal_II_dimension> > &) const;
+template
+void FEValuesBase<deal_II_dimension>::get_function_2nd_derivatives (const BlockVector<2,double> &,
+                                            vector<Tensor<2,deal_II_dimension> > &) const;
+template
+void FEValuesBase<deal_II_dimension>::get_function_2nd_derivatives (const BlockVector<3,double> &,
+                                            vector<Tensor<2,deal_II_dimension> > &) const;
index acc9dede48b2e4a4d2416349ba697a6328c889fc..b2a6ac9d181deb01d75926c6fa421e84016ceb51 100644 (file)
@@ -69,8 +69,8 @@ void FEQ1Mapping<1>::compute_jacobian_gradients (const DoFHandler<1>::cell_itera
 
 template <>
 void FEQ1Mapping<2>::compute_jacobian_matrices (const DoFHandler<2>::cell_iterator &cell,
-                                                   const vector<Point<2> > &unit_points,
-                                                   vector<Tensor<2,2> >    &jacobians) 
+                                               const vector<Point<2> > &unit_points,
+                                               vector<Tensor<2,2> >    &jacobians) 
 {
   const unsigned int dim = 2;
   

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.