]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Templatize the create_mass_matrix functions also on the template argument of the...
authorhartmann <hartmann@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 22 Oct 2003 15:52:34 +0000 (15:52 +0000)
committerhartmann <hartmann@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 22 Oct 2003 15:52:34 +0000 (15:52 +0000)
git-svn-id: https://svn.dealii.org/trunk@8125 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/numerics/matrices.h
deal.II/deal.II/source/numerics/matrices.cc
deal.II/doc/news/c-4-0.html

index 2dd2d7cd061fba177b34f27440283ba274cecd94..582a19565536aae83b30a35a3b2ff2ba3867ed9f 100644 (file)
@@ -203,11 +203,11 @@ class MatrixCreator
                                      * See the general doc of this class
                                      * for more information.
                                      */
-    template <int dim>
+    template <int dim, typename number>
     static void create_mass_matrix (const Mapping<dim>       &mapping,
                                    const DoFHandler<dim>    &dof,
                                    const Quadrature<dim>    &q,
-                                   SparseMatrix<double>     &matrix,
+                                   SparseMatrix<number>     &matrix,
                                    const Function<dim> * const a = 0);
 
                                     /**
@@ -215,10 +215,10 @@ class MatrixCreator
                                      * function, see above, with
                                      * @p{mapping=MappingQ1<dim>()}.
                                      */
-    template <int dim>
+    template <int dim, typename number>
     static void create_mass_matrix (const DoFHandler<dim>    &dof,
                                    const Quadrature<dim>    &q,
-                                   SparseMatrix<double>     &matrix,
+                                   SparseMatrix<number>     &matrix,
                                    const Function<dim> * const a = 0);
 
                                     /**
@@ -234,11 +234,11 @@ class MatrixCreator
                                      * See the general doc of this
                                      * class for more information.
                                      */
-    template <int dim>
+    template <int dim, typename number>
     static void create_mass_matrix (const Mapping<dim>       &mapping,
                                    const DoFHandler<dim>    &dof,
                                    const Quadrature<dim>    &q,
-                                   SparseMatrix<double>     &matrix,
+                                   SparseMatrix<number>     &matrix,
                                    const Function<dim>      &rhs,
                                    Vector<double>           &rhs_vector,
                                    const Function<dim> * const a = 0);
@@ -248,10 +248,10 @@ class MatrixCreator
                                      * function, see above, with
                                      * @p{mapping=MappingQ1<dim>()}.
                                      */
-    template <int dim>
+    template <int dim, typename number>
     static void create_mass_matrix (const DoFHandler<dim>    &dof,
                                    const Quadrature<dim>    &q,
-                                   SparseMatrix<double>     &matrix,
+                                   SparseMatrix<number>     &matrix,
                                    const Function<dim>      &rhs,
                                    Vector<double>           &rhs_vector,
                                    const Function<dim> * const a = 0);
@@ -445,12 +445,12 @@ class MatrixCreator
                                      * used to synchronise access to
                                      * the matrix.
                                      */
-    template <int dim>
+    template <int dim, typename number>
     static
     void create_mass_matrix_1 (const Mapping<dim>       &mapping,
                               const DoFHandler<dim>    &dof,
                               const Quadrature<dim>    &q,
-                              SparseMatrix<double>     &matrix,
+                              SparseMatrix<number>     &matrix,
                               const Function<dim> * const a,
                               const IteratorRange<dim>  range,
                               Threads::ThreadMutex     &mutex);
@@ -465,12 +465,12 @@ class MatrixCreator
                                      * used to synchronise access to
                                      * the matrix.
                                      */
-    template <int dim>
+    template <int dim, typename number>
     static
     void create_mass_matrix_2 (const Mapping<dim>       &mapping,
                               const DoFHandler<dim>    &dof,
                               const Quadrature<dim>    &q,
-                              SparseMatrix<double>     &matrix,
+                              SparseMatrix<number>     &matrix,
                               const Function<dim>      &rhs,
                               Vector<double>           &rhs_vector,
                               const Function<dim> * const a,
index faa104f129bafbc1a1ff1e4e3ab35e473c3102bd..3e146b28ff00189292f88554f4756cc0707c5f73 100644 (file)
@@ -58,11 +58,11 @@ MatrixCreator::IteratorRange<dim>::IteratorRange (const iterator_pair &ip)
 
 
 
-template <int dim>
+template <int dim, typename number>
 void MatrixCreator::create_mass_matrix (const Mapping<dim>       &mapping,
                                        const DoFHandler<dim>    &dof,
                                        const Quadrature<dim>    &q,
-                                       SparseMatrix<double>     &matrix,
+                                       SparseMatrix<number>     &matrix,
                                        const Function<dim> * const coefficient)
 {
   Assert (matrix.m() == dof.n_dofs(),
@@ -88,11 +88,11 @@ void MatrixCreator::create_mass_matrix (const Mapping<dim>       &mapping,
   typedef void (*create_mass_matrix_1_t) (const Mapping<dim>       &mapping,
                                          const DoFHandler<dim>    &dof,
                                          const Quadrature<dim>    &q,
-                                         SparseMatrix<double>     &matrix,
+                                         SparseMatrix<number>     &matrix,
                                          const Function<dim> * const coefficient,
                                          const IteratorRange<dim>  range,
                                          Threads::ThreadMutex     &mutex);
-  create_mass_matrix_1_t p = &MatrixCreator::template create_mass_matrix_1<dim>;
+  create_mass_matrix_1_t p = &MatrixCreator::template create_mass_matrix_1<dim,number>;
   for (unsigned int thread=0; thread<n_threads; ++thread)
     threads += Threads::spawn (p)(mapping, dof, q, matrix, coefficient,
                                   thread_ranges[thread], mutex);
@@ -101,11 +101,11 @@ void MatrixCreator::create_mass_matrix (const Mapping<dim>       &mapping,
 
 
 
-template <int dim>
+template <int dim, typename number>
 void MatrixCreator::create_mass_matrix_1 (const Mapping<dim>       &mapping,
                                          const DoFHandler<dim>    &dof,
                                          const Quadrature<dim>    &q,
-                                         SparseMatrix<double>     &matrix,
+                                         SparseMatrix<number>     &matrix,
                                          const Function<dim> * const coefficient,
                                          const IteratorRange<dim>  range,
                                          Threads::ThreadMutex     &mutex)
@@ -125,7 +125,7 @@ void MatrixCreator::create_mass_matrix_1 (const Mapping<dim>       &mapping,
         coefficient->n_components==1 ||
         coefficient->n_components==n_components, ExcComponentMismatch());
   
-  FullMatrix<double>  cell_matrix (dofs_per_cell, dofs_per_cell);
+  FullMatrix<number>  cell_matrix (dofs_per_cell, dofs_per_cell);
   std::vector<double> coefficient_values (n_q_points);
   std::vector<Vector<double> > coefficient_vector_values (n_q_points,
                                                          Vector<double> (n_components));
@@ -222,10 +222,10 @@ void MatrixCreator::create_mass_matrix_1 (const Mapping<dim>       &mapping,
 }
 
 
-template <int dim>
+template <int dim, typename number>
 void MatrixCreator::create_mass_matrix (const DoFHandler<dim>    &dof,
                                        const Quadrature<dim>    &q,
-                                       SparseMatrix<double>     &matrix,
+                                       SparseMatrix<number>     &matrix,
                                        const Function<dim> * const coefficient)
 {
   Assert (DEAL_II_COMPAT_MAPPING, ExcCompatibility("mapping"));
@@ -234,11 +234,11 @@ void MatrixCreator::create_mass_matrix (const DoFHandler<dim>    &dof,
 }
 
 
-template <int dim>
+template <int dim, typename number>
 void MatrixCreator::create_mass_matrix (const Mapping<dim>       &mapping,
                                        const DoFHandler<dim>    &dof,
                                        const Quadrature<dim>    &q,
-                                       SparseMatrix<double>     &matrix,
+                                       SparseMatrix<number>     &matrix,
                                        const Function<dim>      &rhs,
                                        Vector<double>           &rhs_vector,
                                        const Function<dim> * const coefficient)
@@ -266,13 +266,13 @@ void MatrixCreator::create_mass_matrix (const Mapping<dim>       &mapping,
   typedef void (*create_mass_matrix_2_t) (const Mapping<dim>       &mapping,
                                          const DoFHandler<dim>    &dof,
                                          const Quadrature<dim>    &q,
-                                         SparseMatrix<double>     &matrix,
+                                         SparseMatrix<number>     &matrix,
                                          const Function<dim>      &rhs,
                                          Vector<double>           &rhs_vector,
                                          const Function<dim> * const coefficient,
                                          const IteratorRange<dim>  range,
                                          Threads::ThreadMutex     &mutex);
-  create_mass_matrix_2_t p = &MatrixCreator::template create_mass_matrix_2<dim>;
+  create_mass_matrix_2_t p = &MatrixCreator::template create_mass_matrix_2<dim,number>;
   for (unsigned int thread=0; thread<n_threads; ++thread)
     threads += Threads::spawn (p)(mapping, dof, q, matrix, rhs,
                                   rhs_vector, coefficient,
@@ -282,12 +282,12 @@ void MatrixCreator::create_mass_matrix (const Mapping<dim>       &mapping,
 
 
 
-template <int dim>
+template <int dim, typename number>
 void
 MatrixCreator::create_mass_matrix_2 (const Mapping<dim>       &mapping,
                                     const DoFHandler<dim>    &dof,
                                     const Quadrature<dim>    &q,
-                                    SparseMatrix<double>     &matrix,
+                                    SparseMatrix<number>     &matrix,
                                     const Function<dim>      &rhs,
                                     Vector<double>           &rhs_vector,
                                     const Function<dim> * const coefficient,
@@ -311,7 +311,7 @@ MatrixCreator::create_mass_matrix_2 (const Mapping<dim>       &mapping,
         coefficient->n_components==1 ||
         coefficient->n_components==n_components, ExcComponentMismatch());
   
-  FullMatrix<double>  cell_matrix (dofs_per_cell, dofs_per_cell);
+  FullMatrix<number>  cell_matrix (dofs_per_cell, dofs_per_cell);
   Vector<double>      local_rhs (dofs_per_cell);
   std::vector<double> rhs_values (fe_values.n_quadrature_points);
   std::vector<double> coefficient_values (n_q_points);
@@ -418,10 +418,10 @@ MatrixCreator::create_mass_matrix_2 (const Mapping<dim>       &mapping,
 }
 
 
-template <int dim>
+template <int dim, typename number>
 void MatrixCreator::create_mass_matrix (const DoFHandler<dim>    &dof,
                                        const Quadrature<dim>    &q,
-                                       SparseMatrix<double>     &matrix,
+                                       SparseMatrix<number>     &matrix,
                                        const Function<dim>      &rhs,
                                        Vector<double>           &rhs_vector,
                                        const Function<dim> * const coefficient)
@@ -1297,6 +1297,39 @@ void MatrixCreator::create_mass_matrix<deal_II_dimension>
  Vector<double>           &rhs_vector,
  const Function<deal_II_dimension> * const coefficient);
 
+
+template
+void MatrixCreator::create_mass_matrix<deal_II_dimension>
+(const Mapping<deal_II_dimension>       &mapping,
+ const DoFHandler<deal_II_dimension>    &dof,
+ const Quadrature<deal_II_dimension>    &q,
+ SparseMatrix<float>     &matrix,
+ const Function<deal_II_dimension> * const coefficient);
+template
+void MatrixCreator::create_mass_matrix<deal_II_dimension>
+(const DoFHandler<deal_II_dimension>    &dof,
+ const Quadrature<deal_II_dimension>    &q,
+ SparseMatrix<float>     &matrix,
+ const Function<deal_II_dimension> * const coefficient);
+template
+void MatrixCreator::create_mass_matrix<deal_II_dimension>
+(const Mapping<deal_II_dimension>       &mapping,
+ const DoFHandler<deal_II_dimension>    &dof,
+ const Quadrature<deal_II_dimension>    &q,
+ SparseMatrix<float>     &matrix,
+ const Function<deal_II_dimension>      &rhs,
+ Vector<double>           &rhs_vector,
+ const Function<deal_II_dimension> * const coefficient);
+template
+void MatrixCreator::create_mass_matrix<deal_II_dimension>
+(const DoFHandler<deal_II_dimension>    &dof,
+ const Quadrature<deal_II_dimension>    &q,
+ SparseMatrix<float>     &matrix,
+ const Function<deal_II_dimension>      &rhs,
+ Vector<double>           &rhs_vector,
+ const Function<deal_II_dimension> * const coefficient);
+
+
 #if deal_II_dimension != 1
 template
 void
index 982686e07b818383e2a861f0d4cf180697a7b286..dba522d04602735a491f5832882d6a83a2200f54 100644 (file)
@@ -256,6 +256,18 @@ contributor's names are abbreviated by WB (Wolfgang Bangerth), GK
 <h3>deal.II</h3>
 
 <ol>
+  <li> <p>
+       Improved: The <code
+       class="member">MatrixCreator::create_mass_matrix</code>
+       functions are now templatized also on the template argument of
+       the <code class="class">SparseMatrix</code> class. This allows
+       invoking this function for <code
+       class="class">SparseMatrix&lt;double&gt;</code> and <code
+       class="class">SparseMatrix&lt;float&gt;</code> objects.
+       <br>
+       (RH 2003/10/22)
+       </p>
+
   <li> <p>
        New: There is now also a function <code
        class="member">MGDoFCellAccessor::neighbor_child_on_subface</code>

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.