]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Use the generic instantiation scheme instead of manually listing the number types.
authorkronbichler <kronbichler@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 11 Mar 2013 15:15:41 +0000 (15:15 +0000)
committerkronbichler <kronbichler@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 11 Mar 2013 15:15:41 +0000 (15:15 +0000)
git-svn-id: https://svn.dealii.org/trunk@28853 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/source/lac/CMakeLists.txt
deal.II/source/lac/lapack_full_matrix.cc
deal.II/source/lac/lapack_full_matrix.inst.in [new file with mode: 0644]
deal.II/source/lac/trilinos_sparse_matrix.cc
deal.II/source/lac/trilinos_sparse_matrix.inst.in [new file with mode: 0644]

index 7244f645ffd50a126cb522a3d19e6a4141a9a69c..b460e66faa9ea1e4793b8a7655279baeb8c0fb3b 100644 (file)
@@ -81,12 +81,14 @@ SET(inst_in_files
   chunk_sparse_matrix.inst.in
   constraint_matrix.inst.in
   full_matrix.inst.in
+  lapack_full_matrix.inst.in
   operator.inst.in
   parallel_vector.inst.in
   precondition_block.inst.in
   relaxation_block.inst.in
   sparse_matrix_ez.inst.in
   sparse_matrix.inst.in
+  trilinos_sparse_matrix.inst.in
   trilinos_vector_base.inst.in
   vector.inst.in
   vector_memory.inst.in
index 5582ef2b5a83c2e4b062c781f954e7b44ce0e90b..f61fcbd94af9e69f85b778ae33a9b4fb8549f240 100644 (file)
@@ -891,27 +891,7 @@ PreconditionLU<number>::Tvmult(BlockVector<number> &dst,
 
 
 
-template class LAPACKFullMatrix<double>;
-template LAPACKFullMatrix<double> &
-LAPACKFullMatrix<double>::operator = (const FullMatrix<double> &M);
-template LAPACKFullMatrix<double> &
-LAPACKFullMatrix<double>::operator = (const FullMatrix<float> &M);
-
-template class LAPACKFullMatrix<float>;
-template LAPACKFullMatrix<float> &
-LAPACKFullMatrix<float>::operator = (const FullMatrix<double> &M);
-template LAPACKFullMatrix<float> &
-LAPACKFullMatrix<float>::operator = (const FullMatrix<float> &M);
-
-template class LAPACKFullMatrix<long double>;
-template LAPACKFullMatrix<long double> &
-LAPACKFullMatrix<long double>::operator = (const FullMatrix<long double> &M);
-template LAPACKFullMatrix<long double> &
-LAPACKFullMatrix<long double>::operator = (const FullMatrix<double> &M);
-template LAPACKFullMatrix<long double> &
-LAPACKFullMatrix<long double>::operator = (const FullMatrix<float> &M);
-
-template class PreconditionLU<double>;
-template class PreconditionLU<float>;
+#include "lapack_full_matrix.inst"
+
 
 DEAL_II_NAMESPACE_CLOSE
diff --git a/deal.II/source/lac/lapack_full_matrix.inst.in b/deal.II/source/lac/lapack_full_matrix.inst.in
new file mode 100644 (file)
index 0000000..ddde87b
--- /dev/null
@@ -0,0 +1,25 @@
+//---------------------------------------------------------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2013 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.
+//
+//---------------------------------------------------------------------------
+
+
+for (S : REAL_SCALARS)
+  {
+    template class LAPACKFullMatrix<S>;
+    template class PreconditionLU<S>;
+  }
+
+for (S1, S2 : REAL_SCALARS)
+  {
+    template LAPACKFullMatrix<S1> &
+    LAPACKFullMatrix<S1>::operator = (const FullMatrix<S2> &M);
+  }
index f01ada3dcca7b9bab3e106eaf627e5bbec800d58..29030f6e47a829181a452c3b99585eb601d964d8 100644 (file)
@@ -1388,12 +1388,17 @@ namespace TrilinosWrappers
             sizeof(int)*local_size() +
             static_memory);
   }
+}
+
 
 
+// explicit instantiations
+#include "trilinos_sparse_matrix.inst"
 
 
-  // explicit instantiations
-  //
+// TODO: put these instantiations into generic file
+namespace TrilinosWrappers
+{
   template void
   SparseMatrix::reinit (const dealii::SparsityPattern &);
   template void
@@ -1442,64 +1447,6 @@ namespace TrilinosWrappers
                         const CompressedSetSparsityPattern &,
                         const bool);
 
-  template void
-  SparseMatrix::reinit (const dealii::SparseMatrix<float> &,
-                        const double,
-                        const bool,
-                        const dealii::SparsityPattern *);
-  template void
-  SparseMatrix::reinit (const dealii::SparseMatrix<double> &,
-                        const double,
-                        const bool,
-                        const dealii::SparsityPattern *);
-  template void
-  SparseMatrix::reinit (const dealii::SparseMatrix<long double> &,
-                        const double,
-                        const bool,
-                        const dealii::SparsityPattern *);
-
-  template void
-  SparseMatrix::reinit (const Epetra_Map &,
-                        const dealii::SparseMatrix<float> &,
-                        const double,
-                        const bool,
-                        const dealii::SparsityPattern *);
-  template void
-  SparseMatrix::reinit (const Epetra_Map &,
-                        const dealii::SparseMatrix<double> &,
-                        const double,
-                        const bool,
-                        const dealii::SparsityPattern *);
-  template void
-  SparseMatrix::reinit (const Epetra_Map &,
-                        const dealii::SparseMatrix<long double> &,
-                        const double,
-                        const bool,
-                        const dealii::SparsityPattern *);
-
-  template void
-  SparseMatrix::reinit (const Epetra_Map &,
-                        const Epetra_Map &,
-                        const dealii::SparseMatrix<float> &,
-                        const double,
-                        const bool,
-                        const dealii::SparsityPattern *);
-  template void
-  SparseMatrix::reinit (const Epetra_Map &,
-                        const Epetra_Map &,
-                        const dealii::SparseMatrix<double> &,
-                        const double,
-                        const bool,
-                        const dealii::SparsityPattern *);
-  template void
-  SparseMatrix::reinit (const Epetra_Map &,
-                        const Epetra_Map &,
-                        const dealii::SparseMatrix<long double> &,
-                        const double,
-                        const bool,
-                        const dealii::SparsityPattern *);
-
-
 }
 
 DEAL_II_NAMESPACE_CLOSE
diff --git a/deal.II/source/lac/trilinos_sparse_matrix.inst.in b/deal.II/source/lac/trilinos_sparse_matrix.inst.in
new file mode 100644 (file)
index 0000000..32fc182
--- /dev/null
@@ -0,0 +1,38 @@
+//---------------------------------------------------------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2013 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.
+//
+//---------------------------------------------------------------------------
+
+
+for (S : REAL_SCALARS)
+  {
+    namespace TrilinosWrappers
+    \{
+      template void
+      SparseMatrix::reinit (const dealii::SparseMatrix<S> &,
+                            const double,
+                            const bool,
+                            const dealii::SparsityPattern *);
+      template void
+      SparseMatrix::reinit (const Epetra_Map &,
+                            const dealii::SparseMatrix<S> &,
+                            const double,
+                            const bool,
+                            const dealii::SparsityPattern *);
+      template void
+      SparseMatrix::reinit (const Epetra_Map &,
+                            const Epetra_Map &,
+                            const dealii::SparseMatrix<S> &,
+                            const double,
+                            const bool,
+                            const dealii::SparsityPattern *);
+    \}
+  }

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.