// ---------------------------------------------------------------------
// $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.
//
// 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"));
// 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);
+ }
}
// 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);
+ }
}
// 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);
+ }
}
// 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);
+ }
}
// ---------------------------------------------------------------------
// $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.
//
*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));
}
}
--- /dev/null
+// ---------------------------------------------------------------------
+// $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();
+
+}