*/
template <typename number>
number
- master_compute_threshold (const Vector<number> &criteria,
- const std::pair<double,double> global_min_and_max,
- const double target_error,
- MPI_Comm mpi_communicator)
+ compute_threshold (const Vector<number> &criteria,
+ const std::pair<double,double> global_min_and_max,
+ const double target_error,
+ MPI_Comm mpi_communicator)
{
double interesting_range[2] = { global_min_and_max.first,
global_min_and_max.second
};
adjust_interesting_range (interesting_range);
+ const unsigned int master_mpi_rank = 0;
unsigned int iteration = 0;
do
{
MPI_Bcast (&interesting_range[0], 2, MPI_DOUBLE,
- 0, mpi_communicator);
+ master_mpi_rank, mpi_communicator);
if (interesting_range[0] == interesting_range[1])
{
double final_threshold = std::min (interesting_range[0],
global_min_and_max.second);
MPI_Bcast (&final_threshold, 1, MPI_DOUBLE,
- 0, mpi_communicator);
+ master_mpi_rank, mpi_communicator);
return final_threshold;
}
const double test_threshold
= (interesting_range[0] > 0
?
- std::exp((std::log(interesting_range[0]) +
- std::log(interesting_range[1])) / 2)
+ std::sqrt(interesting_range[0] * interesting_range[1])
:
(interesting_range[0] + interesting_range[1]) / 2);
double total_error;
MPI_Reduce (&my_error, &total_error, 1, MPI_DOUBLE,
- MPI_SUM, 0, mpi_communicator);
+ MPI_SUM, master_mpi_rank, mpi_communicator);
// now adjust the range. if we have to many cells, we take
// the upper half of the previous range, otherwise the lower
// half. if we have hit the right number, then set the range
- // to the exact value
+ // to the exact value.
+ // slave nodes also update their own interesting_range, however
+ // their results are not significant since the values will be
+ // overwritten by MPI_Bcast from the master node in next loop.
if (total_error > target_error)
interesting_range[0] = test_threshold;
else if (total_error < target_error)
Assert (false, ExcInternalError());
return -1;
}
-
-
- /**
- * The corresponding function to
- * the one above, to be run on the
- * slaves.
- */
- template <typename number>
- number
- slave_compute_threshold (const Vector<number> &criteria,
- MPI_Comm mpi_communicator)
- {
- do
- {
- double interesting_range[2] = { -1, -1 };
- MPI_Bcast (&interesting_range[0], 2, MPI_DOUBLE,
- 0, mpi_communicator);
-
- if (interesting_range[0] == interesting_range[1])
- {
- // we got the final value. the master adjusts the
- // value one more time, so receive it from there
- double final_threshold = -1;
- MPI_Bcast (&final_threshold, 1, MPI_DOUBLE,
- 0, mpi_communicator);
-
- return final_threshold;
- }
-
- // count how many elements
- // there are that are bigger
- // than the following trial
- // threshold
- const double test_threshold
- = (interesting_range[0] > 0
- ?
- std::exp((std::log(interesting_range[0]) +
- std::log(interesting_range[1])) / 2)
- :
- (interesting_range[0] + interesting_range[1]) / 2);
-
- double my_error = 0;
- for (unsigned int i=0; i<criteria.size(); ++i)
- if (criteria(i) > test_threshold)
- my_error += criteria(i);
-
- MPI_Reduce (&my_error, 0, 1, MPI_DOUBLE,
- MPI_SUM, 0, mpi_communicator);
- }
- while (true);
-
- Assert (false, ExcInternalError());
- return -1;
- }
}
}
const double total_error
= compute_global_sum (locally_owned_indicators,
mpi_communicator);
-
- // from here on designate a
- // master and slaves
double top_threshold, bottom_threshold;
- if (Utilities::MPI::this_mpi_process (mpi_communicator) == 0)
- {
- // this is the master
- // processor
- top_threshold
- =
- RefineAndCoarsenFixedFraction::
- master_compute_threshold (locally_owned_indicators,
- global_min_and_max,
- top_fraction_of_error *
- total_error,
- mpi_communicator);
- // compute bottom
- // threshold only if
- // necessary. otherwise
- // use a threshold lower
- // than the smallest
- // value we have locally
- if (bottom_fraction_of_error > 0)
- bottom_threshold
- =
- RefineAndCoarsenFixedFraction::
- master_compute_threshold (locally_owned_indicators,
- global_min_and_max,
- (1-bottom_fraction_of_error) *
- total_error,
- mpi_communicator);
- else
- {
- bottom_threshold = *std::min_element (criteria.begin(),
- criteria.end());
- bottom_threshold -= std::fabs(bottom_threshold);
- }
- }
+ top_threshold =
+ RefineAndCoarsenFixedFraction::
+ compute_threshold (locally_owned_indicators,
+ global_min_and_max,
+ top_fraction_of_error *
+ total_error,
+ mpi_communicator);
+ // compute bottom
+ // threshold only if
+ // necessary. otherwise
+ // use a threshold lower
+ // than the smallest
+ // value we have locally
+ if (bottom_fraction_of_error > 0)
+ bottom_threshold =
+ RefineAndCoarsenFixedFraction::
+ compute_threshold (locally_owned_indicators,
+ global_min_and_max,
+ (1-bottom_fraction_of_error) *
+ total_error,
+ mpi_communicator);
else
{
- // this is a slave
- // processor
- top_threshold
- =
- RefineAndCoarsenFixedFraction::
- slave_compute_threshold (locally_owned_indicators,
- mpi_communicator);
-
- // compute bottom
- // threshold only if
- // necessary. otherwise
- // use a threshold lower
- // than the smallest
- // value we have locally
- if (bottom_fraction_of_error > 0)
- bottom_threshold
- =
- RefineAndCoarsenFixedFraction::
- slave_compute_threshold (locally_owned_indicators,
- mpi_communicator);
- else
- {
- bottom_threshold = *std::min_element (criteria.begin(),
- criteria.end());
- bottom_threshold -= std::fabs(bottom_threshold);
- }
+ bottom_threshold = *std::min_element (criteria.begin(),
+ criteria.end());
+ bottom_threshold -= std::fabs(bottom_threshold);
}
// now refine the mesh