]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add test for sadd with partly empty Trilinos vectors 497/head
authorDaniel Arndt <d.arndt@math.uni-goettingen.de>
Thu, 5 Feb 2015 10:58:50 +0000 (11:58 +0100)
committerDaniel Arndt <d.arndt@math.uni-goettingen.de>
Thu, 5 Feb 2015 10:58:50 +0000 (11:58 +0100)
tests/mpi/trilinos_sadd_03.cc [new file with mode: 0644]
tests/mpi/trilinos_sadd_03.with_trilinos=true.mpirun=3.output [new file with mode: 0644]

diff --git a/tests/mpi/trilinos_sadd_03.cc b/tests/mpi/trilinos_sadd_03.cc
new file mode 100644 (file)
index 0000000..f1387f9
--- /dev/null
@@ -0,0 +1,163 @@
+// ---------------------------------------------------------------------
+//
+// 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 there are processes that doesn't own anything 
+
+#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, ghosted2, distributed, distributed2;
+  const unsigned int my_id = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD);
+  const unsigned int n_proc = Utilities::MPI::n_mpi_processes(MPI_COMM_WORLD);
+  //All processes except for the last one own elements
+  const int entries_per_process = (my_id==n_proc-1)?0:10;
+  IndexSet locally_owned(entries_per_process*(n_proc-1));
+  IndexSet locally_relevant(entries_per_process*(n_proc-1));
+  const int begin_index = my_id*entries_per_process;
+  const unsigned int end_index = (my_id+1)*entries_per_process;
+  const unsigned int local_begin = std::max(0, begin_index-entries_per_process/2);
+  const unsigned int local_end = entries_per_process*(n_proc-1);
+
+  if (my_id != n_proc-1)
+  {
+    locally_owned.add_range(begin_index, end_index);
+    locally_relevant.add_range(local_begin, local_end);
+  }
+  distributed.reinit(locally_owned, MPI_COMM_WORLD);
+  distributed2.reinit(locally_owned, MPI_COMM_WORLD);
+  ghosted.reinit (locally_owned, locally_relevant, MPI_COMM_WORLD);
+  ghosted2.reinit (locally_owned, locally_relevant, MPI_COMM_WORLD);
+
+  //sadd(s,v) with ghosted vector
+  distributed=1.;
+  ghosted=distributed;
+  ghosted *= 1.6;
+  distributed.sadd (2., ghosted);
+  ghosted = distributed;
+  deallog << "sadd (s, v) ghosted" <<std::endl;
+  if (my_id==0)
+    for (unsigned int i=local_begin; i<local_end; ++i)
+      deallog << i << ": " << ghosted(i) << std::endl;
+
+  //sadd(s,v) with distributed vector
+  distributed=1.;
+  distributed.sadd (1.5, distributed);
+  ghosted = distributed;
+  deallog << "sadd (s, v) distributed" <<std::endl;
+  if (my_id==0)
+    for (unsigned int i=local_begin; i<local_end; ++i)
+      deallog << i << ": " << ghosted(i) << std::endl;
+
+  //sadd(s,a,v) with ghosted vector
+  distributed=1.;
+  ghosted=distributed;
+  ghosted *= 1.6;
+  distributed.sadd (2., 3., ghosted);
+  ghosted = distributed;
+  deallog << "sadd (s, a, v) ghosted" <<std::endl;
+  if (my_id==0)
+    for (unsigned int i=local_begin; i<local_end; ++i)
+      deallog << i << ": " << ghosted(i) << std::endl;
+
+  //sadd(s,a,v) with distributed vector
+  distributed=1.;
+  distributed.sadd (2., 3., distributed);
+  ghosted = distributed;
+  deallog << "sadd (s, a, v) distributed" <<std::endl;
+  if (my_id==0)
+    for (unsigned int i=local_begin; i<local_end; ++i)
+      deallog << i << ": " << ghosted(i) << std::endl;
+
+  //sadd(s,a,v,b,w) with ghosted vector
+  distributed=1.;
+  ghosted=distributed;
+  ghosted *= 1.6;
+  ghosted2 = distributed;
+  ghosted2 *= 0.354;
+  distributed.sadd (2., 3., ghosted, 4., ghosted2);
+  ghosted = distributed;
+  deallog << "sadd (s, a, v, b, w) ghosted" <<std::endl;
+  if (my_id==0)
+    for (unsigned int i=local_begin; i<local_end; ++i)
+      deallog << i << ": " << ghosted(i) << std::endl;
+
+  //sadd(s,a,v,b,w) with distributed vector
+  distributed=1.;
+  distributed2=1.6;
+  distributed.sadd (2., 3., distributed, 4., distributed2);
+  ghosted = distributed;
+  deallog << "sadd (s, a, v, b, w) distributed" <<std::endl;
+  if (my_id==0)
+    for (unsigned int i=local_begin; i<local_end; ++i)
+      deallog << i << ": " << ghosted(i) << std::endl;
+
+  //sadd(s,a,v,b,w,c,x) with ghosted vector
+  distributed=1.;
+  ghosted=distributed;
+  distributed2 = 7.22;
+  ghosted2=distributed2;
+  distributed.sadd (2., 3., ghosted, 4., ghosted2, 5., distributed2);
+  ghosted = distributed;
+  deallog << "sadd (s, a, v, b, w, c, x) ghosted" <<std::endl;
+  if (my_id==0)
+    for (unsigned int i=local_begin; i<local_end; ++i)
+      deallog << i << ": " << ghosted(i) << std::endl;
+
+  //sadd(s,a,v,b,w,c,x) with distributed vector
+  distributed=1.;
+  distributed2 = 7.22;
+  distributed.sadd (2., 3., distributed, 4., distributed, 5., distributed2);
+  ghosted = distributed;
+  deallog << "sadd (s, a, v, b, w, c, x) distributed" <<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_03.with_trilinos=true.mpirun=3.output b/tests/mpi/trilinos_sadd_03.with_trilinos=true.mpirun=3.output
new file mode 100644 (file)
index 0000000..86ec565
--- /dev/null
@@ -0,0 +1,169 @@
+
+DEAL:0::sadd (s, v) ghosted
+DEAL:0::0: 3.600
+DEAL:0::1: 3.600
+DEAL:0::2: 3.600
+DEAL:0::3: 3.600
+DEAL:0::4: 3.600
+DEAL:0::5: 3.600
+DEAL:0::6: 3.600
+DEAL:0::7: 3.600
+DEAL:0::8: 3.600
+DEAL:0::9: 3.600
+DEAL:0::10: 3.600
+DEAL:0::11: 3.600
+DEAL:0::12: 3.600
+DEAL:0::13: 3.600
+DEAL:0::14: 3.600
+DEAL:0::15: 3.600
+DEAL:0::16: 3.600
+DEAL:0::17: 3.600
+DEAL:0::18: 3.600
+DEAL:0::19: 3.600
+DEAL:0::sadd (s, v) distributed
+DEAL:0::0: 2.500
+DEAL:0::1: 2.500
+DEAL:0::2: 2.500
+DEAL:0::3: 2.500
+DEAL:0::4: 2.500
+DEAL:0::5: 2.500
+DEAL:0::6: 2.500
+DEAL:0::7: 2.500
+DEAL:0::8: 2.500
+DEAL:0::9: 2.500
+DEAL:0::10: 2.500
+DEAL:0::11: 2.500
+DEAL:0::12: 2.500
+DEAL:0::13: 2.500
+DEAL:0::14: 2.500
+DEAL:0::15: 2.500
+DEAL:0::16: 2.500
+DEAL:0::17: 2.500
+DEAL:0::18: 2.500
+DEAL:0::19: 2.500
+DEAL:0::sadd (s, a, v) ghosted
+DEAL:0::0: 6.800
+DEAL:0::1: 6.800
+DEAL:0::2: 6.800
+DEAL:0::3: 6.800
+DEAL:0::4: 6.800
+DEAL:0::5: 6.800
+DEAL:0::6: 6.800
+DEAL:0::7: 6.800
+DEAL:0::8: 6.800
+DEAL:0::9: 6.800
+DEAL:0::10: 6.800
+DEAL:0::11: 6.800
+DEAL:0::12: 6.800
+DEAL:0::13: 6.800
+DEAL:0::14: 6.800
+DEAL:0::15: 6.800
+DEAL:0::16: 6.800
+DEAL:0::17: 6.800
+DEAL:0::18: 6.800
+DEAL:0::19: 6.800
+DEAL:0::sadd (s, a, v) distributed
+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) ghosted
+DEAL:0::0: 8.216
+DEAL:0::1: 8.216
+DEAL:0::2: 8.216
+DEAL:0::3: 8.216
+DEAL:0::4: 8.216
+DEAL:0::5: 8.216
+DEAL:0::6: 8.216
+DEAL:0::7: 8.216
+DEAL:0::8: 8.216
+DEAL:0::9: 8.216
+DEAL:0::10: 8.216
+DEAL:0::11: 8.216
+DEAL:0::12: 8.216
+DEAL:0::13: 8.216
+DEAL:0::14: 8.216
+DEAL:0::15: 8.216
+DEAL:0::16: 8.216
+DEAL:0::17: 8.216
+DEAL:0::18: 8.216
+DEAL:0::19: 8.216
+DEAL:0::sadd (s, a, v, b, w) distributed
+DEAL:0::0: 11.40
+DEAL:0::1: 11.40
+DEAL:0::2: 11.40
+DEAL:0::3: 11.40
+DEAL:0::4: 11.40
+DEAL:0::5: 11.40
+DEAL:0::6: 11.40
+DEAL:0::7: 11.40
+DEAL:0::8: 11.40
+DEAL:0::9: 11.40
+DEAL:0::10: 11.40
+DEAL:0::11: 11.40
+DEAL:0::12: 11.40
+DEAL:0::13: 11.40
+DEAL:0::14: 11.40
+DEAL:0::15: 11.40
+DEAL:0::16: 11.40
+DEAL:0::17: 11.40
+DEAL:0::18: 11.40
+DEAL:0::19: 11.40
+DEAL:0::sadd (s, a, v, b, w, c, x) ghosted
+DEAL:0::0: 69.98
+DEAL:0::1: 69.98
+DEAL:0::2: 69.98
+DEAL:0::3: 69.98
+DEAL:0::4: 69.98
+DEAL:0::5: 69.98
+DEAL:0::6: 69.98
+DEAL:0::7: 69.98
+DEAL:0::8: 69.98
+DEAL:0::9: 69.98
+DEAL:0::10: 69.98
+DEAL:0::11: 69.98
+DEAL:0::12: 69.98
+DEAL:0::13: 69.98
+DEAL:0::14: 69.98
+DEAL:0::15: 69.98
+DEAL:0::16: 69.98
+DEAL:0::17: 69.98
+DEAL:0::18: 69.98
+DEAL:0::19: 69.98
+DEAL:0::sadd (s, a, v, b, w, c, x) distributed
+DEAL:0::0: 45.10
+DEAL:0::1: 45.10
+DEAL:0::2: 45.10
+DEAL:0::3: 45.10
+DEAL:0::4: 45.10
+DEAL:0::5: 45.10
+DEAL:0::6: 45.10
+DEAL:0::7: 45.10
+DEAL:0::8: 45.10
+DEAL:0::9: 45.10
+DEAL:0::10: 45.10
+DEAL:0::11: 45.10
+DEAL:0::12: 45.10
+DEAL:0::13: 45.10
+DEAL:0::14: 45.10
+DEAL:0::15: 45.10
+DEAL:0::16: 45.10
+DEAL:0::17: 45.10
+DEAL:0::18: 45.10
+DEAL:0::19: 45.10

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.