Describes the distribution of a distributable linear algebra type object. Typically this is a container (such as a DoubleVector) or an operator (e.g Preconditioner or LinearSolver). This object is used in both serial and parallel implementations. In the serial context (no MPI) this just contains an integer indicating the number of rows. In parallel either each processor holds a subset of the set of global rows. (each processor contains only a single continuous block of rows - parametised with variables denoting the first row and the number of local rows) or, all rows are be duplicated across all processors. In parallel this object also contains an OomphCommunicator object which primarily contains the MPI_Comm communicator associated with this object. More...
#include <linear_algebra_distribution.h>
Public Member Functions | |
LinearAlgebraDistribution () | |
Default Constructor - creates a Distribution that has not been setup. More... | |
LinearAlgebraDistribution (const OomphCommunicator &comm, const unsigned &first_row_, const unsigned &n_row_local, const unsigned &n_row=0) | |
Constructor. Takes the first_row, nrow_local (both for this processor) and nrow as arguments. If nrow is not provided or equal to 0 then it will be computed automatically. More... | |
LinearAlgebraDistribution (const OomphCommunicator &comm, const unsigned &n_row, const bool &distributed_=true) | |
Constructor. Takes the number of global rows and uniformly distributes them over the processors if distributed = true (default), if distributed = false then every row is duplicated on every processor. More... | |
LinearAlgebraDistribution (const OomphCommunicator *const comm_pt, const unsigned &first_row_, const unsigned &n_row_local, const unsigned &n_row=0) | |
Constructor. Takes the first_row, nrow_local (both for this processor) and nrow as arguments. If nrow is not provided or equal to 0 then it will be computed automatically. More... | |
LinearAlgebraDistribution (const OomphCommunicator *const comm_pt, const unsigned &n_row, const bool &distributed_=true) | |
Constructor. Takes the number of global rows and uniformly distributes them over the processors if distributed = true (default), if distributed = false then every row is duplicated on every processor. More... | |
LinearAlgebraDistribution (const LinearAlgebraDistribution &old_dist) | |
Copy Constructor. More... | |
LinearAlgebraDistribution (const LinearAlgebraDistribution *old_dist_pt) | |
pointer based copy constructor More... | |
~LinearAlgebraDistribution () | |
Destructor. More... | |
void | operator= (const LinearAlgebraDistribution &old_dist) |
Assignment Operator. More... | |
void | build (const OomphCommunicator *const comm_pt, const unsigned &first_row, const unsigned &nrow_local, const unsigned &nrow=0) |
Sets the distribution. Takes first_row, nrow_local and nrow as arguments. If nrow is not provided or equal to 0 then it is computed automatically. More... | |
void | build (const OomphCommunicator *const comm_pt, const unsigned &nrow, const bool &distributed=true) |
Build the LinearAlgebraDistribution. if distributed = true (default) then uniformly distribute nrow over all processors where processors 0 holds approximately the first nrow/n_proc, processor 1 holds the next nrow/n_proc and so on... or if distributed = false then every row is held on every processor. More... | |
void | build (const LinearAlgebraDistribution &new_dist) |
Copy the argument distribution. Also a helper method for the =assignment operator and copy constructor. More... | |
void | build (const LinearAlgebraDistribution *new_dist_pt) |
Copy the argument distribution. Also a helper method for the =assignment operator and copy constructor. More... | |
void | clear () |
clears the distribution More... | |
unsigned | nrow () const |
access function to the number of global rows. More... | |
unsigned | nrow_local () const |
access function for the num of local rows on this processor. If no MPI then Nrow is returned. More... | |
unsigned | nrow_local (const unsigned &p) const |
access function for the num of local rows on this processor. If no MPI the nrow is returned More... | |
unsigned | first_row () const |
access function for the first row on this processor. If not distributed then this is just zero. More... | |
unsigned | first_row (const unsigned &p) const |
access function for the first row on the p-th processor More... | |
bool | distributed () const |
access function to the distributed - indicates whether the distribution is serial or distributed More... | |
OomphCommunicator * | communicator_pt () const |
const access to the communicator pointer More... | |
bool | built () const |
bool | operator== (const LinearAlgebraDistribution &other_dist) const |
== Operator More... | |
bool | operator!= (const LinearAlgebraDistribution &other_dist) const |
!= operator More... | |
unsigned | global_to_local_row_map (const unsigned &global_i) const |
return the local index corresponding to the global index More... | |
unsigned | rank_of_global_row (const unsigned i) const |
return the processor rank of the global row number i More... | |
Vector< unsigned > | nrow_local_vector () const |
return the nrow_local Vector More... | |
Vector< unsigned > | first_row_vector () const |
return the first_row Vector More... | |
Private Attributes | |
unsigned | Nrow |
the number of global rows More... | |
Vector< unsigned > | Nrow_local |
the number of local rows on the processor More... | |
Vector< unsigned > | First_row |
the first row on this processor More... | |
bool | Distributed |
OomphCommunicator * | Comm_pt |
the pointer to the MPI communicator object in this distribution More... | |
Friends | |
std::ostream & | operator<< (std::ostream &stream, LinearAlgebraDistribution &dist) |
<< operator More... | |
Describes the distribution of a distributable linear algebra type object. Typically this is a container (such as a DoubleVector) or an operator (e.g Preconditioner or LinearSolver). This object is used in both serial and parallel implementations. In the serial context (no MPI) this just contains an integer indicating the number of rows. In parallel either each processor holds a subset of the set of global rows. (each processor contains only a single continuous block of rows - parametised with variables denoting the first row and the number of local rows) or, all rows are be duplicated across all processors. In parallel this object also contains an OomphCommunicator object which primarily contains the MPI_Comm communicator associated with this object.
Definition at line 68 of file linear_algebra_distribution.h.
|
inline |
Default Constructor - creates a Distribution that has not been setup.
Definition at line 75 of file linear_algebra_distribution.h.
Referenced by oomph::DistributableLinearAlgebraObject::DistributableLinearAlgebraObject().
|
inline |
Constructor. Takes the first_row, nrow_local (both for this processor) and nrow as arguments. If nrow is not provided or equal to 0 then it will be computed automatically.
Definition at line 80 of file linear_algebra_distribution.h.
References build().
|
inline |
Constructor. Takes the number of global rows and uniformly distributes them over the processors if distributed = true (default), if distributed = false then every row is duplicated on every processor.
Definition at line 92 of file linear_algebra_distribution.h.
References build().
|
inline |
Constructor. Takes the first_row, nrow_local (both for this processor) and nrow as arguments. If nrow is not provided or equal to 0 then it will be computed automatically.
Definition at line 103 of file linear_algebra_distribution.h.
References build().
|
inline |
Constructor. Takes the number of global rows and uniformly distributes them over the processors if distributed = true (default), if distributed = false then every row is duplicated on every processor.
Definition at line 115 of file linear_algebra_distribution.h.
References build().
|
inline |
|
inline |
pointer based copy constructor
Definition at line 131 of file linear_algebra_distribution.h.
References build().
|
inline |
void oomph::LinearAlgebraDistribution::build | ( | const OomphCommunicator *const | comm_pt, |
const unsigned & | first_row, | ||
const unsigned & | local_nrow, | ||
const unsigned & | global_nrow = 0 |
||
) |
Sets the distribution. Takes first_row, nrow_local and nrow as arguments. If nrow is not provided or equal to 0 then it is computed automatically.
Sets the distribution. Takes first_row, local_nrow and global_nrow as arguments. If global_nrow is not provided or equal to 0 then it is computed automatically
Definition at line 41 of file linear_algebra_distribution.cc.
References Comm_pt, Distributed, first_row(), First_row, Nrow, and Nrow_local.
Referenced by oomph::PeriodicOrbitAssemblyHandler< NNODE_1D >::adapt_temporal_mesh(), oomph::Problem::assign_eqn_numbers(), build(), oomph::LinearAlgebraDistributionHelpers::concatenate(), oomph::HypreHelpers::create_HYPRE_Matrix(), oomph::FoldHandler::FoldHandler(), oomph::BlockPreconditioner< CRDoubleMatrix >::get_concatenated_block(), oomph::CRDoubleMatrix::get_matrix_transpose(), oomph::HopfHandler::HopfHandler(), LinearAlgebraDistribution(), oomph::SumOfMatrices::multiply(), operator=(), oomph::BlockPitchForkLinearSolver::resolve(), oomph::MumpsSolver::solve(), oomph::FoldHandler::solve_augmented_block_system(), oomph::FoldHandler::solve_block_system(), oomph::HopfHandler::solve_complex_system(), oomph::FoldHandler::solve_full_system(), oomph::HopfHandler::solve_full_system(), oomph::HopfHandler::solve_standard_system(), oomph::Problem::synchronise_eqn_numbers(), oomph::FoldHandler::~FoldHandler(), and oomph::HopfHandler::~HopfHandler().
void oomph::LinearAlgebraDistribution::build | ( | const OomphCommunicator *const | comm_pt, |
const unsigned & | nrow, | ||
const bool & | distributed = true |
||
) |
Build the LinearAlgebraDistribution. if distributed = true (default) then uniformly distribute nrow over all processors where processors 0 holds approximately the first nrow/n_proc, processor 1 holds the next nrow/n_proc and so on... or if distributed = false then every row is held on every processor.
Uniformly distribute global_nrow over all processors where processors 0 holds approximately the first global_nrow/n_proc, processor 1 holds the next global_nrow/n_proc and so on...
Definition at line 170 of file linear_algebra_distribution.cc.
References Comm_pt, Distributed, First_row, Nrow, and Nrow_local.
void oomph::LinearAlgebraDistribution::build | ( | const LinearAlgebraDistribution & | new_dist | ) |
Copy the argument distribution. Also a helper method for the =assignment operator and copy constructor.
helper method for the =assignment operator and copy constructor
Definition at line 222 of file linear_algebra_distribution.cc.
References Comm_pt, communicator_pt(), distributed(), Distributed, First_row, first_row_vector(), nrow(), Nrow, Nrow_local, and nrow_local_vector().
|
inline |
Copy the argument distribution. Also a helper method for the =assignment operator and copy constructor.
Definition at line 172 of file linear_algebra_distribution.h.
References build().
|
inline |
if the communicator_pt is null then the distribution is not setup then false is returned, otherwise return true
Definition at line 350 of file linear_algebra_distribution.h.
References Comm_pt, and operator==().
Referenced by oomph::CRDoubleMatrix::add(), oomph::MumpsSolver::backsub(), oomph::SuperLUSolver::backsub_distributed(), oomph::DoubleVector::build(), oomph::DoubleMultiVector::build(), oomph::ComplexDampedJacobi< MATRIX >::complex_solve_helper(), oomph::LinearAlgebraDistributionHelpers::concatenate(), oomph::CRDoubleMatrixHelpers::concatenate(), oomph::CRDoubleMatrixHelpers::concatenate_without_communication(), oomph::TrilinosEpetraHelpers::create_distributed_epetra_vector(), oomph::TrilinosEpetraHelpers::multiply(), oomph::LagrangeEnforcedFlowPreconditioner::preconditioner_solve(), oomph::NavierStokesSchurComplementPreconditioner::preconditioner_solve(), oomph::BlockPitchForkLinearSolver::resolve(), oomph::DoubleMultiVector::shallow_build(), oomph::MumpsSolver::solve(), oomph::CG< MATRIX >::solve(), oomph::CG< MATRIX >::solve_helper(), and oomph::BiCGStab< MATRIX >::solve_helper().
|
inline |
clears the distribution
Definition at line 178 of file linear_algebra_distribution.h.
References Comm_pt, First_row, Nrow, and Nrow_local.
Referenced by oomph::BlockPreconditioner< CRDoubleMatrix >::block_setup(), oomph::BlockPreconditioner< CRDoubleMatrix >::clear_block_preconditioner_base(), oomph::DoubleVectorHelpers::concatenate_without_communication(), oomph::MGSolver< DIM >::setup_mg_structures(), and oomph::DoubleVectorHelpers::split_without_communication().
|
inline |
const access to the communicator pointer
Definition at line 343 of file linear_algebra_distribution.h.
References Comm_pt.
Referenced by oomph::SuperLUSolver::backsub_distributed(), build(), oomph::Problem::calculate_continuation_derivatives_fd_helper(), oomph::Problem::calculate_continuation_derivatives_helper(), oomph::DoubleVectorHelpers::concatenate(), oomph::CRDoubleMatrixHelpers::concatenate(), oomph::DoubleVectorHelpers::concatenate_without_communication(), oomph::CRDoubleMatrixHelpers::concatenate_without_communication(), oomph::TrilinosEpetraHelpers::copy_to_oomphlib_vector(), oomph::TrilinosEpetraHelpers::create_distributed_epetra_matrix(), oomph::TrilinosEpetraHelpers::create_distributed_epetra_matrix_for_aztecoo(), oomph::TrilinosEpetraHelpers::create_distributed_epetra_vector(), oomph::TrilinosEpetraHelpers::create_epetra_map(), oomph::HypreHelpers::create_HYPRE_Matrix(), oomph::HypreHelpers::create_HYPRE_Vector(), oomph::DoubleVectorHaloScheme::DoubleVectorHaloScheme(), oomph::MumpsSolver::factorise(), oomph::SuperLUSolver::factorise(), oomph::SuperLUSolver::factorise_distributed(), oomph::Problem::get_jacobian(), oomph::HypreInterface::hypre_matrix_setup(), oomph::BlockPreconditioner< CRDoubleMatrix >::internal_get_block(), oomph::HelmholtzMGPreconditioner< DIM >::interpolation_matrix_set(), oomph::MGSolver< DIM >::interpolation_matrix_set(), oomph::TrilinosEpetraHelpers::multiply(), oomph::DenseDoubleMatrix::multiply(), oomph::CCDoubleMatrix::multiply(), oomph::CRDoubleMatrix::multiply_transpose(), oomph::DenseDoubleMatrix::multiply_transpose(), oomph::CCDoubleMatrix::multiply_transpose(), Anasazi::MultiVecTraits< double, oomph::DoubleMultiVector >::MvTransMv(), oomph::Problem::newton_solve_continuation(), oomph::DoubleVector::output(), oomph::DoubleMultiVector::output(), oomph::DoubleVector::redistribute(), oomph::DoubleMultiVector::redistribute(), oomph::CRDoubleMatrix::redistribute(), oomph::MGSolver< DIM >::self_test(), oomph::Preconditioner::setup(), oomph::MatrixVectorProduct::setup(), oomph::HelmholtzMGPreconditioner< DIM >::setup_mg_structures(), oomph::MGSolver< DIM >::setup_mg_structures(), oomph::HelmholtzMGPreconditioner< DIM >::setup_smoothers(), oomph::MGSolver< DIM >::setup_smoothers(), oomph::MumpsSolver::solve(), oomph::DenseLU::solve(), oomph::TrilinosAztecOOSolver::solve(), oomph::DoubleVectorHelpers::split(), and oomph::DoubleVectorHelpers::split_without_communication().
|
inline |
access function to the distributed - indicates whether the distribution is serial or distributed
Definition at line 337 of file linear_algebra_distribution.h.
References Distributed.
Referenced by oomph::SuperLUSolver::backsub_serial(), build(), oomph::Problem::calculate_continuation_derivatives_fd_helper(), oomph::Problem::calculate_continuation_derivatives_helper(), oomph::HelmholtzSmoother::check_validity_of_solve_helper_inputs(), oomph::Smoother::check_validity_of_solve_helper_inputs(), oomph::ComplexGMRES< MATRIX >::complex_solve_helper(), oomph::LinearAlgebraDistributionHelpers::concatenate(), oomph::TrilinosEpetraHelpers::create_distributed_epetra_matrix(), oomph::TrilinosEpetraHelpers::create_distributed_epetra_vector(), oomph::TrilinosEpetraHelpers::create_epetra_map(), oomph::HypreHelpers::create_HYPRE_Matrix(), oomph::HypreHelpers::create_HYPRE_Vector(), oomph::DoubleVectorHaloScheme::DoubleVectorHaloScheme(), oomph::Problem::get_eigenproblem_matrices(), oomph::Problem::get_jacobian(), oomph::BlockPreconditioner< CRDoubleMatrix >::internal_get_block(), oomph::Problem::newton_solve_continuation(), oomph::operator<<(), operator=(), oomph::ILUZeroPreconditioner< CCDoubleMatrix >::preconditioner_solve(), oomph::ILUZeroPreconditioner< CRDoubleMatrix >::preconditioner_solve(), oomph::DoubleVector::redistribute(), oomph::DoubleMultiVector::redistribute(), oomph::CRDoubleMatrix::redistribute(), oomph::MatrixVectorProduct::setup(), oomph::MumpsSolver::solve(), oomph::DenseLU::solve(), oomph::SuperLUSolver::solve(), oomph::CG< MATRIX >::solve_helper(), oomph::BiCGStab< MATRIX >::solve_helper(), oomph::GMRES< MATRIX >::solve_helper(), oomph::HelmholtzGMRESMG< MATRIX >::solve_helper(), and oomph::HelmholtzFGMRESMG< MATRIX >::solve_helper().
|
inline |
access function for the first row on this processor. If not distributed then this is just zero.
Definition at line 268 of file linear_algebra_distribution.h.
References Comm_pt, Distributed, and First_row.
Referenced by oomph::NavierStokesSchurComplementPreconditioner::assemble_inv_press_and_veloc_mass_matrix_diagonal(), oomph::PressureBasedSolidLSCPreconditioner::assemble_mass_matrix_diagonal(), oomph::BlockPreconditioner< CRDoubleMatrix >::block_setup(), build(), oomph::DoubleVectorHelpers::concatenate(), oomph::CRDoubleMatrixHelpers::concatenate(), oomph::TrilinosEpetraHelpers::create_distributed_epetra_matrix(), oomph::TrilinosEpetraHelpers::create_distributed_epetra_matrix_for_aztecoo(), oomph::TrilinosEpetraHelpers::create_epetra_map(), oomph::HypreHelpers::create_HYPRE_Matrix(), oomph::HypreHelpers::create_HYPRE_Vector(), oomph::DoubleVectorHaloScheme::DoubleVectorHaloScheme(), oomph::Problem::global_dof_pt(), oomph::operator<<(), operator=(), oomph::Problem::parallel_sparse_assemble(), oomph::PitchForkHandler::PitchForkHandler(), rank_of_global_row(), oomph::DoubleVector::redistribute(), oomph::DoubleMultiVector::redistribute(), oomph::CRDoubleMatrix::redistribute(), and oomph::CRDoubleMatrix::sparse_indexed_output_with_offset().
|
inline |
access function for the first row on the p-th processor
Definition at line 297 of file linear_algebra_distribution.h.
References Comm_pt, Distributed, and First_row.
|
inline |
return the first_row Vector
Definition at line 412 of file linear_algebra_distribution.h.
References First_row.
Referenced by build().
|
inline |
return the local index corresponding to the global index
Definition at line 373 of file linear_algebra_distribution.h.
References Nrow, and nrow_local().
|
inline |
access function to the number of global rows.
Definition at line 193 of file linear_algebra_distribution.h.
References Nrow.
Referenced by oomph::PressureBasedSolidLSCPreconditioner::assemble_mass_matrix_diagonal(), oomph::MumpsSolver::backsub(), oomph::SuperLUSolver::backsub_distributed(), build(), oomph::SuperLUSolver::clean_up_memory(), oomph::TrilinosEpetraHelpers::create_distributed_epetra_matrix(), oomph::TrilinosEpetraHelpers::create_distributed_epetra_vector(), oomph::TrilinosEpetraHelpers::create_epetra_map(), oomph::SuperLUSolver::factorise_distributed(), oomph::Problem::get_eigenproblem_matrices(), oomph::Problem::get_jacobian(), oomph::Problem::get_residuals(), oomph::CRDoubleMatrix::multiply(), oomph::CRDoubleMatrix::multiply_transpose(), oomph::Problem::ndof(), oomph::operator<<(), operator=(), oomph::DoubleVector::redistribute(), oomph::DoubleMultiVector::redistribute(), oomph::CRDoubleMatrix::redistribute(), and oomph::Problem::synchronise_eqn_numbers().
|
inline |
access function for the num of local rows on this processor. If no MPI then Nrow is returned.
Definition at line 200 of file linear_algebra_distribution.h.
References Comm_pt, Distributed, Nrow, and Nrow_local.
Referenced by oomph::Problem::adapt(), oomph::Problem::adaptive_unsteady_newton_solve(), oomph::OomphLibPreconditionerEpetraOperator::ApplyInverse(), oomph::Problem::arc_length_step_solve_helper(), oomph::NavierStokesSchurComplementPreconditioner::assemble_inv_press_and_veloc_mass_matrix_diagonal(), oomph::PressureBasedSolidLSCPreconditioner::assemble_mass_matrix_diagonal(), oomph::BlockPreconditioner< CRDoubleMatrix >::block_setup(), oomph::Problem::calculate_continuation_derivatives_fd_helper(), oomph::Problem::calculate_continuation_derivatives_helper(), oomph::CRDoubleMatrixHelpers::concatenate(), oomph::CRDoubleMatrixHelpers::concatenate_without_communication(), oomph::CRDoubleMatrix::CRDoubleMatrix(), oomph::TrilinosEpetraHelpers::create_distributed_epetra_matrix(), oomph::TrilinosEpetraHelpers::create_distributed_epetra_matrix_for_aztecoo(), oomph::TrilinosEpetraHelpers::create_epetra_map(), oomph::HypreHelpers::create_HYPRE_Matrix(), oomph::HypreHelpers::create_HYPRE_Vector(), oomph::Problem::get_eigenproblem_matrices(), oomph::Problem::get_hessian_vector_products(), oomph::Problem::get_jacobian(), oomph::Problem::get_residuals(), oomph::Problem::global_dof_pt(), global_to_local_row_map(), oomph::Problem::newton_solve(), oomph::Problem::newton_solve_continuation(), oomph::operator<<(), operator=(), oomph::Problem::parallel_sparse_assemble(), oomph::PitchForkHandler::PitchForkHandler(), rank_of_global_row(), oomph::DoubleVector::redistribute(), oomph::DoubleMultiVector::redistribute(), oomph::CRDoubleMatrix::redistribute(), oomph::Problem::restore_dof_values(), oomph::PitchForkHandler::solve_block_system(), oomph::PitchForkHandler::solve_full_system(), oomph::Problem::store_current_dof_values(), and oomph::PitchForkHandler::~PitchForkHandler().
|
inline |
access function for the num of local rows on this processor. If no MPI the nrow is returned
Definition at line 229 of file linear_algebra_distribution.h.
References Comm_pt, Distributed, Nrow, and Nrow_local.
|
inline |
return the nrow_local Vector
Definition at line 406 of file linear_algebra_distribution.h.
References Nrow_local.
Referenced by build().
|
inline |
|
inline |
Assignment Operator.
Definition at line 144 of file linear_algebra_distribution.h.
References build(), distributed(), first_row(), nrow(), and nrow_local().
bool oomph::LinearAlgebraDistribution::operator== | ( | const LinearAlgebraDistribution & | other_dist | ) | const |
== Operator
operator==
Definition at line 266 of file linear_algebra_distribution.cc.
References Comm_pt, Distributed, First_row, i, Nrow, and Nrow_local.
Referenced by built().
|
inline |
return the processor rank of the global row number i
Definition at line 395 of file linear_algebra_distribution.h.
References first_row(), and nrow_local().
Referenced by oomph::DoubleVectorHelpers::concatenate(), oomph::CRDoubleMatrixHelpers::concatenate(), oomph::DoubleVectorHaloScheme::DoubleVectorHaloScheme(), and oomph::Problem::parallel_sparse_assemble().
|
friend |
<< operator
Definition at line 315 of file linear_algebra_distribution.cc.
Referenced by operator!=().
|
private |
the pointer to the MPI communicator object in this distribution
Definition at line 434 of file linear_algebra_distribution.h.
Referenced by build(), built(), clear(), communicator_pt(), first_row(), nrow_local(), operator==(), and ~LinearAlgebraDistribution().
|
private |
flag to indicate whether this distribution describes an object that is distributed over the processors of Comm_pt (true) or duplicated over the processors of Comm_pt (false)
Definition at line 431 of file linear_algebra_distribution.h.
Referenced by build(), distributed(), first_row(), nrow_local(), and operator==().
|
private |
the first row on this processor
Definition at line 426 of file linear_algebra_distribution.h.
Referenced by build(), clear(), first_row(), first_row_vector(), and operator==().
|
private |
the number of global rows
Definition at line 420 of file linear_algebra_distribution.h.
Referenced by build(), clear(), global_to_local_row_map(), nrow(), nrow_local(), and operator==().
|
private |
the number of local rows on the processor
Definition at line 423 of file linear_algebra_distribution.h.
Referenced by build(), clear(), nrow_local(), nrow_local_vector(), and operator==().