]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Temporary fix for Issue 191, allow sadd for Trilinos vectors
authordaniel.arndt <daniel.arndt@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 5 Mar 2014 11:48:30 +0000 (11:48 +0000)
committerdaniel.arndt <daniel.arndt@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 5 Mar 2014 11:48:30 +0000 (11:48 +0000)
git-svn-id: https://svn.dealii.org/trunk@32612 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/include/deal.II/lac/trilinos_vector_base.h
deal.II/source/lac/trilinos_vector_base.cc
tests/mpi/trilinos_sadd_01.with_trilinos=true.mpirun=2.output [moved from tests/mpi/trilinos_sadd_01.mpirun=2.output with 100% similarity]
tests/mpi/trilinos_sadd_02.cc [new file with mode: 0644]
tests/mpi/trilinos_sadd_02.with_trilinos=true.mpirun=1.output [new file with mode: 0644]
tests/mpi/trilinos_sadd_02.with_trilinos=true.mpirun=2.output [new file with mode: 0644]

index dfb77a229e69c80e2b1256133cdd9cfe8e248e27..84845480eb4f69d9ea508664a7357ab61c18a7e8 100644 (file)
@@ -1,7 +1,7 @@
 // ---------------------------------------------------------------------
 // $Id$
 //
-// Copyright (C) 2008 - 2014 by the deal.II authors
+// Copyright (C) 2008 - 2013 by the deal.II authors
 //
 // This file is part of the deal.II library.
 //
@@ -1383,7 +1383,7 @@ namespace TrilinosWrappers
             // use pre-allocated vector for non-local entries if it exists for
             // addition operation
             const TrilinosWrappers::types::int_type my_row = nonlocal_vector->Map().LID(static_cast<TrilinosWrappers::types::int_type>(row));
-            Assert(my_row != -1,
+            Assert(my_row != -1, 
                    ExcMessage("Attempted to write into off-processor vector entry "
                               "that has not be specified as being writable upon "
                                "initialization"));
@@ -1704,14 +1704,22 @@ namespace TrilinosWrappers
     // if we have ghost values, do not allow
     // writing to this vector at all.
     Assert (!has_ghost_elements(), ExcGhostsPresent());
-    Assert (local_size() == v.local_size(),
-            ExcDimensionMismatch(local_size(), v.local_size()));
+    Assert (size() == v.size(),
+            ExcDimensionMismatch (size(), v.size()));
 
     Assert (numbers::is_finite(s), ExcNumberNotFinite());
 
-    const int ierr = vector->Update(1., *(v.vector), s);
-
-    AssertThrow (ierr == 0, ExcTrilinosError(ierr));
+    if(local_size() == v.local_size())
+    {
+      const int ierr = vector->Update(1., *(v.vector), s);
+      AssertThrow (ierr == 0, ExcTrilinosError(ierr));
+    }
+    else
+    {
+      VectorBase tmp = v;
+      tmp *= s;
+      this->add(tmp, true);
+    }
   }
 
 
@@ -1725,15 +1733,23 @@ namespace TrilinosWrappers
     // if we have ghost values, do not allow
     // writing to this vector at all.
     Assert (!has_ghost_elements(), ExcGhostsPresent());
-    Assert (local_size() == v.local_size(),
-            ExcDimensionMismatch(local_size(), v.local_size()));
-
+    Assert (size() == v.size(),
+            ExcDimensionMismatch (size(), v.size()));
     Assert (numbers::is_finite(s), ExcNumberNotFinite());
     Assert (numbers::is_finite(a), ExcNumberNotFinite());
 
-    const int ierr = vector->Update(a, *(v.vector), s);
-
-    AssertThrow (ierr == 0, ExcTrilinosError(ierr));
+    if(local_size() == v.local_size())
+    {
+      const int ierr = vector->Update(a, *(v.vector), s);
+      AssertThrow (ierr == 0, ExcTrilinosError(ierr));
+    }
+    else
+    {
+      (*this)*=s;
+      VectorBase tmp = v;
+      tmp *= a;
+      this->add(tmp, true);
+    }
   }
 
 
@@ -1749,18 +1765,29 @@ namespace TrilinosWrappers
     // if we have ghost values, do not allow
     // writing to this vector at all.
     Assert (!has_ghost_elements(), ExcGhostsPresent());
-    Assert (local_size() == v.local_size(),
-            ExcDimensionMismatch(local_size(), v.local_size()));
-    Assert (local_size() == w.local_size(),
-            ExcDimensionMismatch(local_size(), w.local_size()));
-
+    Assert (size() == v.size(),
+            ExcDimensionMismatch (size(), v.size()));
+    Assert (size() == w.size(),
+            ExcDimensionMismatch (size(), w.size()));
     Assert (numbers::is_finite(s), ExcNumberNotFinite());
     Assert (numbers::is_finite(a), ExcNumberNotFinite());
     Assert (numbers::is_finite(b), ExcNumberNotFinite());
-
-    const int ierr = vector->Update(a, *(v.vector), b, *(w.vector), s);
-
-    AssertThrow (ierr == 0, ExcTrilinosError(ierr));
+    
+    if(local_size() == v.local_size() && local_size() == w.local_size())
+    {
+      const int ierr = vector->Update(a, *(v.vector), b, *(w.vector), s);
+      AssertThrow (ierr == 0, ExcTrilinosError(ierr));
+    }
+    else
+    {
+      (*this)*=s;
+      VectorBase tmp = v;
+      tmp *= a;
+      this->add(tmp, true);
+      tmp = w;
+      tmp *= b;
+      this->add(tmp, true);
+    }
   }
 
 
@@ -1778,27 +1805,44 @@ namespace TrilinosWrappers
     // if we have ghost values, do not allow
     // writing to this vector at all.
     Assert (!has_ghost_elements(), ExcGhostsPresent());
-    Assert (local_size() == v.local_size(),
-            ExcDimensionMismatch(local_size(), v.local_size()));
-    Assert (local_size() == w.local_size(),
-            ExcDimensionMismatch(local_size(), w.local_size()));
-    Assert (local_size() == x.local_size(),
-            ExcDimensionMismatch(local_size(), x.local_size()));
-
+    Assert (size() == v.size(),
+            ExcDimensionMismatch (size(), v.size()));
+    Assert (size() == w.size(),
+            ExcDimensionMismatch (size(), w.size()));
+    Assert (size() == x.size(),
+            ExcDimensionMismatch (size(), x.size()));
     Assert (numbers::is_finite(s), ExcNumberNotFinite());
     Assert (numbers::is_finite(a), ExcNumberNotFinite());
     Assert (numbers::is_finite(b), ExcNumberNotFinite());
     Assert (numbers::is_finite(c), ExcNumberNotFinite());
 
-    // Update member can only
-    // input two other vectors so
-    // do it in two steps
-    const int ierr = vector->Update(a, *(v.vector), b, *(w.vector), s);
-    AssertThrow (ierr == 0, ExcTrilinosError(ierr));
-
-    const int jerr = vector->Update(c, *(x.vector), 1.);
-    Assert (jerr == 0, ExcTrilinosError(jerr));
-    (void)jerr; // removes -Wunused-parameter warning in optimized mode
+    if(local_size() == v.local_size()
+       && local_size() == w.local_size()
+       && local_size() == x.local_size())
+    {
+      // Update member can only
+      // input two other vectors so
+      // do it in two steps
+      const int ierr = vector->Update(a, *(v.vector), b, *(w.vector), s);
+      AssertThrow (ierr == 0, ExcTrilinosError(ierr));
+      
+      const int jerr = vector->Update(c, *(x.vector), 1.);
+      Assert (jerr == 0, ExcTrilinosError(jerr));
+      (void)jerr; // removes -Wunused-parameter warning in optimized mode
+    }
+    else
+    {
+      (*this)*=s;
+      VectorBase tmp = v;
+      tmp *= a;
+      this->add(tmp, true);
+      tmp = w;
+      tmp *= b;
+      this->add(tmp, true);
+      tmp = x;
+      tmp *= c;
+      this->add(tmp, true);
+    }
   }
 
 
index 95a9b1387baa949084d9cc0666fe08f2e8cf7c69..da1e491f16d3361264b4aaf534ad09cf5cad3331 100644 (file)
@@ -1,7 +1,7 @@
 // ---------------------------------------------------------------------
 // $Id$
 //
-// Copyright (C) 2008 - 2013 by the deal.II authors
+// Copyright (C) 2008 - 2014 by the deal.II authors
 //
 // This file is part of the deal.II library.
 //
@@ -291,15 +291,24 @@ namespace TrilinosWrappers
       *this += v;
     else
       {
+        Assert (!has_ghost_elements(), ExcGhostsPresent());
         AssertThrow (size() == v.size(),
                      ExcDimensionMismatch (size(), v.size()));
 
-        Epetra_Import data_exchange (vector->Map(), v.vector->Map());
+        //HACK: For some reason the following doesn't work properly, so do it manually
+        //  Epetra_Import data_exchange (vector->Map(), v.vector->Map());
+        //  int ierr = vector->Import(*v.vector, data_exchange, Add);
+        //  AssertThrow (ierr == 0, ExcTrilinosError(ierr));
+        //  last_action = Add;
 
-        int ierr = vector->Import(*v.vector, data_exchange, Add);
+        Epetra_MultiVector dummy(vector->Map(), 1, false);
+        Epetra_Import data_exchange (dummy.Map(), v.vector->Map());
+
+        int ierr = dummy.Import(*v.vector, data_exchange, Insert);
         AssertThrow (ierr == 0, ExcTrilinosError(ierr));
 
-        last_action = Insert;
+        ierr = vector->Update (1.0, dummy, 1.0);
+        AssertThrow (ierr == 0, ExcTrilinosError(ierr));
       }
   }
 
diff --git a/tests/mpi/trilinos_sadd_02.cc b/tests/mpi/trilinos_sadd_02.cc
new file mode 100644 (file)
index 0000000..a5723be
--- /dev/null
@@ -0,0 +1,111 @@
+// ---------------------------------------------------------------------
+// $Id$
+//
+// Copyright (C) 2011 - 2013 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+
+
+// check correct behaviour of sadd of Trilinos vectors
+// if they have different Epetra maps 
+
+#include "../tests.h"
+#include <deal.II/base/utilities.h>
+#include <deal.II/base/index_set.h>
+#include <deal.II/lac/trilinos_vector.h>
+#include <fstream>
+#include <iostream>
+#include <vector>
+#include <algorithm>
+
+void test ()
+{
+  TrilinosWrappers::MPI::Vector ghosted, distributed;
+  //All processes should own 10 entries
+  const int entries_per_process = 10;
+  const int n_proc = Utilities::MPI::n_mpi_processes(MPI_COMM_WORLD);
+  IndexSet locally_owned(entries_per_process*n_proc);
+  IndexSet locally_relevant(entries_per_process*n_proc);
+  const int my_id = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD);
+  const int begin_index = my_id*entries_per_process;
+  const int end_index = (my_id+1)*entries_per_process;
+  const int local_begin = std::max(0, begin_index-entries_per_process/2);
+  const int local_end = entries_per_process*n_proc;
+
+  locally_owned.add_range(begin_index, end_index);
+  locally_relevant.add_range
+    (local_begin, local_end);
+  distributed.reinit(locally_owned, MPI_COMM_WORLD);
+  ghosted.reinit (locally_owned, locally_relevant, MPI_COMM_WORLD);
+
+  distributed=1.;
+  ghosted=distributed;
+  distributed.sadd (2., ghosted);
+  ghosted = distributed;
+  deallog << "sadd (s, v)" <<std::endl;
+  if (my_id==0)
+    for (unsigned int i=local_begin; i<local_end; ++i)
+      deallog << i << ": " << ghosted(i) << std::endl;
+
+  distributed=1.;
+  ghosted=distributed;
+  distributed.sadd (2., 3., ghosted);
+  ghosted = distributed;
+  deallog << "sadd (s, a, v)" <<std::endl;
+  if (my_id==0)
+    for (unsigned int i=local_begin; i<local_end; ++i)
+      deallog << i << ": " << ghosted(i) << std::endl;
+
+  distributed=1.;
+  ghosted=distributed;
+  distributed.sadd (2., 3., ghosted, 4., ghosted);
+  ghosted = distributed;
+  deallog << "sadd (s, a, v, b, w)" <<std::endl;
+  if (my_id==0)
+    for (unsigned int i=local_begin; i<local_end; ++i)
+      deallog << i << ": " << ghosted(i) << std::endl;
+
+  distributed=1.;
+  ghosted=distributed;
+  distributed.sadd (2., 3., ghosted, 4., ghosted, 5., ghosted);
+  ghosted = distributed;
+  deallog << "sadd (s, a, v, b, w, c, x)" <<std::endl;
+  if (my_id==0)
+    for (unsigned int i=local_begin; i<local_end; ++i)
+      deallog << i << ": " << ghosted(i) << std::endl;
+}
+
+
+
+int main (int argc, char **argv)
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv);
+  
+  unsigned int myid = Utilities::MPI::this_mpi_process (MPI_COMM_WORLD);
+  deallog.push(Utilities::int_to_string(myid));
+  
+  if (myid == 0)
+  {
+    std::ofstream logfile("output");
+    deallog.attach(logfile);
+    deallog << std::setprecision(4);
+    deallog.depth_console(0);
+    deallog.threshold_double(1.e-10);
+    
+    test();
+  }
+  else
+    test();
+  
+}
diff --git a/tests/mpi/trilinos_sadd_02.with_trilinos=true.mpirun=1.output b/tests/mpi/trilinos_sadd_02.with_trilinos=true.mpirun=1.output
new file mode 100644 (file)
index 0000000..ddb6d88
--- /dev/null
@@ -0,0 +1,45 @@
+
+DEAL:0::sadd (s, v)
+DEAL:0::0: 3.000
+DEAL:0::1: 3.000
+DEAL:0::2: 3.000
+DEAL:0::3: 3.000
+DEAL:0::4: 3.000
+DEAL:0::5: 3.000
+DEAL:0::6: 3.000
+DEAL:0::7: 3.000
+DEAL:0::8: 3.000
+DEAL:0::9: 3.000
+DEAL:0::sadd (s, a, v)
+DEAL:0::0: 5.000
+DEAL:0::1: 5.000
+DEAL:0::2: 5.000
+DEAL:0::3: 5.000
+DEAL:0::4: 5.000
+DEAL:0::5: 5.000
+DEAL:0::6: 5.000
+DEAL:0::7: 5.000
+DEAL:0::8: 5.000
+DEAL:0::9: 5.000
+DEAL:0::sadd (s, a, v, b, w)
+DEAL:0::0: 9.000
+DEAL:0::1: 9.000
+DEAL:0::2: 9.000
+DEAL:0::3: 9.000
+DEAL:0::4: 9.000
+DEAL:0::5: 9.000
+DEAL:0::6: 9.000
+DEAL:0::7: 9.000
+DEAL:0::8: 9.000
+DEAL:0::9: 9.000
+DEAL:0::sadd (s, a, v, b, w, c, x)
+DEAL:0::0: 14.00
+DEAL:0::1: 14.00
+DEAL:0::2: 14.00
+DEAL:0::3: 14.00
+DEAL:0::4: 14.00
+DEAL:0::5: 14.00
+DEAL:0::6: 14.00
+DEAL:0::7: 14.00
+DEAL:0::8: 14.00
+DEAL:0::9: 14.00
diff --git a/tests/mpi/trilinos_sadd_02.with_trilinos=true.mpirun=2.output b/tests/mpi/trilinos_sadd_02.with_trilinos=true.mpirun=2.output
new file mode 100644 (file)
index 0000000..73a234d
--- /dev/null
@@ -0,0 +1,85 @@
+
+DEAL:0::sadd (s, v)
+DEAL:0::0: 3.000
+DEAL:0::1: 3.000
+DEAL:0::2: 3.000
+DEAL:0::3: 3.000
+DEAL:0::4: 3.000
+DEAL:0::5: 3.000
+DEAL:0::6: 3.000
+DEAL:0::7: 3.000
+DEAL:0::8: 3.000
+DEAL:0::9: 3.000
+DEAL:0::10: 3.000
+DEAL:0::11: 3.000
+DEAL:0::12: 3.000
+DEAL:0::13: 3.000
+DEAL:0::14: 3.000
+DEAL:0::15: 3.000
+DEAL:0::16: 3.000
+DEAL:0::17: 3.000
+DEAL:0::18: 3.000
+DEAL:0::19: 3.000
+DEAL:0::sadd (s, a, v)
+DEAL:0::0: 5.000
+DEAL:0::1: 5.000
+DEAL:0::2: 5.000
+DEAL:0::3: 5.000
+DEAL:0::4: 5.000
+DEAL:0::5: 5.000
+DEAL:0::6: 5.000
+DEAL:0::7: 5.000
+DEAL:0::8: 5.000
+DEAL:0::9: 5.000
+DEAL:0::10: 5.000
+DEAL:0::11: 5.000
+DEAL:0::12: 5.000
+DEAL:0::13: 5.000
+DEAL:0::14: 5.000
+DEAL:0::15: 5.000
+DEAL:0::16: 5.000
+DEAL:0::17: 5.000
+DEAL:0::18: 5.000
+DEAL:0::19: 5.000
+DEAL:0::sadd (s, a, v, b, w)
+DEAL:0::0: 9.000
+DEAL:0::1: 9.000
+DEAL:0::2: 9.000
+DEAL:0::3: 9.000
+DEAL:0::4: 9.000
+DEAL:0::5: 9.000
+DEAL:0::6: 9.000
+DEAL:0::7: 9.000
+DEAL:0::8: 9.000
+DEAL:0::9: 9.000
+DEAL:0::10: 9.000
+DEAL:0::11: 9.000
+DEAL:0::12: 9.000
+DEAL:0::13: 9.000
+DEAL:0::14: 9.000
+DEAL:0::15: 9.000
+DEAL:0::16: 9.000
+DEAL:0::17: 9.000
+DEAL:0::18: 9.000
+DEAL:0::19: 9.000
+DEAL:0::sadd (s, a, v, b, w, c, x)
+DEAL:0::0: 14.00
+DEAL:0::1: 14.00
+DEAL:0::2: 14.00
+DEAL:0::3: 14.00
+DEAL:0::4: 14.00
+DEAL:0::5: 14.00
+DEAL:0::6: 14.00
+DEAL:0::7: 14.00
+DEAL:0::8: 14.00
+DEAL:0::9: 14.00
+DEAL:0::10: 14.00
+DEAL:0::11: 14.00
+DEAL:0::12: 14.00
+DEAL:0::13: 14.00
+DEAL:0::14: 14.00
+DEAL:0::15: 14.00
+DEAL:0::16: 14.00
+DEAL:0::17: 14.00
+DEAL:0::18: 14.00
+DEAL:0::19: 14.00

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.