From: Justin O'Connor Date: Thu, 2 Dec 2021 20:35:02 +0000 (-0700) Subject: Make Trilinos vector non_negative function work in parallel. X-Git-Tag: v9.4.0-rc1~672^2~3 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=5a32e4e93defa975d6b4264134a585947774522e;p=dealii.git Make Trilinos vector non_negative function work in parallel. --- diff --git a/source/lac/trilinos_vector.cc b/source/lac/trilinos_vector.cc index 3e37356beb..7581f764fe 100644 --- a/source/lac/trilinos_vector.cc +++ b/source/lac/trilinos_vector.cc @@ -757,36 +757,25 @@ namespace TrilinosWrappers bool Vector::is_non_negative() const { - // if this vector is a parallel one, then - // we need to communicate to determine - // the answer to the current - // function. this still has to be - // implemented - AssertThrow(local_size() == size(), ExcNotImplemented()); // get a representation of the vector and // loop over all the elements - TrilinosScalar *start_ptr; - int leading_dimension; - int ierr = vector->ExtractView(&start_ptr, &leading_dimension); - AssertThrow(ierr == 0, ExcTrilinosError(ierr)); - - // TODO: This - // won't work in parallel like - // this. Find out a better way to - // this in that case. - const TrilinosScalar *ptr = start_ptr, *eptr = start_ptr + size(); - bool flag = true; + TrilinosScalar * start_ptr = (*vector)[0]; + const TrilinosScalar *ptr = start_ptr, *eptr = start_ptr + local_size(); + unsigned int flag = 0; while (ptr != eptr) { if (*ptr < 0.0) { - flag = false; + flag = 1; break; } ++ptr; } - return flag; + // in parallel, check that the vector + // is zero on _all_ processors. + const auto max_n_negative = Utilities::MPI::max(flag, get_mpi_communicator()); + return max_n_negative == 0; }