]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Inst. interpolate function for REAL_SERIAL_VECTOR. Add TODO.
authoryoung <young@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 26 Mar 2014 13:40:34 +0000 (13:40 +0000)
committeryoung <young@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 26 Mar 2014 13:40:34 +0000 (13:40 +0000)
git-svn-id: https://svn.dealii.org/branches/branch_petscscalar_complex@32698 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/include/deal.II/numerics/vector_tools.templates.h
deal.II/source/numerics/vector_tools_interpolate.inst.in

index e6c9e9ed3d771e514642981bb6f56c577f83e8d0..e41bd1be6eac16f42e9eafdc0544fe10b9390a6c 100644 (file)
@@ -280,17 +280,20 @@ namespace VectorTools
 
 
 
-
+  // @whattodo TODO: To return functionality for complex numbers, this
+  // needs to be templated with Number and instantiated for
+  // REAL_SERIAL_VECTORS+double and
+  // COMPLEX_SERIAL_VECTORS+complex<double>
   template <int dim, class InVector, class OutVector, int spacedim>
   void
   interpolate (const DoFHandler<dim,spacedim>  &dof_1,
                const DoFHandler<dim,spacedim>  &dof_2,
-               const FullMatrix<double>        &transfer,
+               const FullMatrix<double>        &transfer, // <- FullMatrix<std::complex>
                const InVector                  &data_1,
                OutVector                       &data_2)
   {
     Vector<double> cell_data_1(dof_1.get_fe().dofs_per_cell); // <- Vector<std::complex>
-    Vector<double> cell_data_2(dof_2.get_fe().dofs_per_cell);
+    Vector<double> cell_data_2(dof_2.get_fe().dofs_per_cell); // <- Vector<std::complex>
 
     std::vector<short unsigned int> touch_count (dof_2.n_dofs(), 0); //TODO: check on datatype... kinda strange (UK)
     std::vector<types::global_dof_index>       local_dof_indices (dof_2.get_fe().dofs_per_cell);
@@ -301,13 +304,7 @@ namespace VectorTools
 
     for (; h != endh; ++h, ++l)
       {
-       // @whattodo Some of the problems in petsc_vector_base.h are
-       // instantiated here. Should look into DoFCellAccessor too. It
-       // looks like *maybe* vector_tools needs to be templated (run
-       // grep -R "double" on that file). Additionally, *this* file
-       // needs to be templated with "number".
-       /* h->get_dof_values(data_1, cell_data_1);  */
-       Assert ((false), ExcMessage ("This function is corrupt: @whattodo"));
+       h->get_dof_values(data_1, cell_data_1);  
 
         transfer.vmult(cell_data_2, cell_data_1);
 
@@ -318,8 +315,7 @@ namespace VectorTools
           {
             data_2(local_dof_indices[j]) += cell_data_2(j);
 
-            // count, how often we have
-            // added to this dof
+            // count, how often we have added to this dof
             Assert (touch_count[local_dof_indices[j]] < numbers::internal_face_boundary_id,
                     ExcInternalError());
             ++touch_count[local_dof_indices[j]];
index a21c69a76b34f5339a41d07a1fec4cbdb2fb52ac..9b466454111bf22566edfe6c396b68f6df214c86 100644 (file)
@@ -15,7 +15,7 @@
 // ---------------------------------------------------------------------
 
 
-for (VEC : SERIAL_VECTORS ; deal_II_dimension : DIMENSIONS; deal_II_space_dimension :  SPACE_DIMENSIONS)
+for (VEC : REAL_SERIAL_VECTORS ; deal_II_dimension : DIMENSIONS; deal_II_space_dimension :  SPACE_DIMENSIONS)
   {
 #if deal_II_dimension <= deal_II_space_dimension
     namespace VectorTools \{

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.