Public Member Functions | Private Attributes | Friends | List of all members
oomph::LinearAlgebraDistribution Class Reference

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...
 
OomphCommunicatorcommunicator_pt () const
 const access to the communicator pointer More...
 
bool built () const
 if the communicator_pt is null then the distribution is not setup then false is returned, otherwise return true More...
 
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
 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) More...
 
OomphCommunicatorComm_pt
 the pointer to the MPI communicator object in this distribution More...
 

Friends

std::ostream & operator<< (std::ostream &stream, LinearAlgebraDistribution &dist)
 << operator More...
 

Detailed Description

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 63 of file linear_algebra_distribution.h.

Constructor & Destructor Documentation

◆ LinearAlgebraDistribution() [1/7]

oomph::LinearAlgebraDistribution::LinearAlgebraDistribution ( )
inline

Default Constructor - creates a Distribution that has not been setup.

Definition at line 68 of file linear_algebra_distribution.h.

◆ LinearAlgebraDistribution() [2/7]

oomph::LinearAlgebraDistribution::LinearAlgebraDistribution ( const OomphCommunicator comm,
const unsigned &  first_row_,
const unsigned &  n_row_local,
const unsigned &  n_row = 0 
)
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 73 of file linear_algebra_distribution.h.

References build().

◆ LinearAlgebraDistribution() [3/7]

oomph::LinearAlgebraDistribution::LinearAlgebraDistribution ( const OomphCommunicator comm,
const unsigned &  n_row,
const bool &  distributed_ = true 
)
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 85 of file linear_algebra_distribution.h.

References build().

◆ LinearAlgebraDistribution() [4/7]

oomph::LinearAlgebraDistribution::LinearAlgebraDistribution ( const OomphCommunicator *const  comm_pt,
const unsigned &  first_row_,
const unsigned &  n_row_local,
const unsigned &  n_row = 0 
)
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 96 of file linear_algebra_distribution.h.

References build().

◆ LinearAlgebraDistribution() [5/7]

oomph::LinearAlgebraDistribution::LinearAlgebraDistribution ( const OomphCommunicator *const  comm_pt,
const unsigned &  n_row,
const bool &  distributed_ = true 
)
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 108 of file linear_algebra_distribution.h.

References build().

◆ LinearAlgebraDistribution() [6/7]

oomph::LinearAlgebraDistribution::LinearAlgebraDistribution ( const LinearAlgebraDistribution old_dist)
inline

Copy Constructor.

Definition at line 117 of file linear_algebra_distribution.h.

References build().

◆ LinearAlgebraDistribution() [7/7]

oomph::LinearAlgebraDistribution::LinearAlgebraDistribution ( const LinearAlgebraDistribution old_dist_pt)
inline

pointer based copy constructor

Definition at line 124 of file linear_algebra_distribution.h.

References build().

◆ ~LinearAlgebraDistribution()

oomph::LinearAlgebraDistribution::~LinearAlgebraDistribution ( )
inline

Destructor.

Definition at line 131 of file linear_algebra_distribution.h.

References Comm_pt.

Member Function Documentation

◆ build() [1/4]

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 226 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().

◆ build() [2/4]

void oomph::LinearAlgebraDistribution::build ( const LinearAlgebraDistribution new_dist_pt)
inline

Copy the argument distribution. Also a helper method for the =assignment operator and copy constructor.

Definition at line 165 of file linear_algebra_distribution.h.

References build().

◆ build() [3/4]

void oomph::LinearAlgebraDistribution::build ( const OomphCommunicator *const  comm_pt,
const unsigned &  first_row,
const unsigned &  nrow_local,
const unsigned &  nrow = 0 
)

◆ build() [4/4]

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 176 of file linear_algebra_distribution.cc.

References Comm_pt, Distributed, First_row, Nrow, and Nrow_local.

◆ built()

bool oomph::LinearAlgebraDistribution::built ( ) const
inline

◆ clear()

void oomph::LinearAlgebraDistribution::clear ( )
inline

◆ communicator_pt()

OomphCommunicator* oomph::LinearAlgebraDistribution::communicator_pt ( ) const
inline

const access to the communicator pointer

Definition at line 335 of file linear_algebra_distribution.h.

References Comm_pt.

Referenced by build(), oomph::Problem::calculate_continuation_derivatives_fd_helper(), oomph::Problem::calculate_continuation_derivatives_helper(), oomph::CRDoubleMatrixHelpers::concatenate(), oomph::DoubleVectorHelpers::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::SuperLUSolver::factorise(), oomph::MumpsSolver::factorise(), oomph::SuperLUSolver::factorise_distributed(), oomph::Problem::get_jacobian(), oomph::CRDoubleMatrix::global_matrix(), oomph::HypreInterface::hypre_matrix_setup(), oomph::HypreInterface::hypre_solver_setup(), oomph::BlockPreconditioner< MATRIX >::internal_get_block(), oomph::MGSolver< DIM >::interpolation_matrix_set(), oomph::HelmholtzMGPreconditioner< 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(), operator==(), oomph::DoubleMultiVector::output(), oomph::DoubleVector::output(), oomph::DoubleMultiVector::redistribute(), oomph::DoubleVector::redistribute(), oomph::CRDoubleMatrix::redistribute(), oomph::MGSolver< DIM >::self_test(), oomph::MatrixVectorProduct::setup(), oomph::Preconditioner::setup(), oomph::MGSolver< DIM >::setup_mg_structures(), oomph::HelmholtzMGPreconditioner< DIM >::setup_mg_structures(), oomph::MGSolver< DIM >::setup_smoothers(), oomph::HelmholtzMGPreconditioner< DIM >::setup_smoothers(), oomph::DenseLU::solve(), oomph::TrilinosAztecOOSolver::solve(), oomph::MumpsSolver::solve(), oomph::AugmentedProblemGMRES::solve_helper(), oomph::DoubleVectorHelpers::split(), and oomph::DoubleVectorHelpers::split_without_communication().

◆ distributed()

bool oomph::LinearAlgebraDistribution::distributed ( ) const
inline

access function to the distributed - indicates whether the distribution is serial or distributed

Definition at line 329 of file linear_algebra_distribution.h.

References Distributed.

Referenced by oomph::SuperLUSolver::backsub_serial(), oomph::SuperLUSolver::backsub_transpose_serial(), build(), oomph::Problem::calculate_continuation_derivatives_fd_helper(), oomph::Problem::calculate_continuation_derivatives_helper(), oomph::Smoother::check_validity_of_solve_helper_inputs(), 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::DistributableLinearAlgebraObject::distributed(), oomph::DoubleVectorHaloScheme::DoubleVectorHaloScheme(), oomph::Problem::get_eigenproblem_matrices(), oomph::Problem::get_jacobian(), oomph::BlockPreconditioner< MATRIX >::internal_get_block(), oomph::Problem::newton_solve_continuation(), operator==(), oomph::DoubleMultiVector::redistribute(), oomph::DoubleVector::redistribute(), oomph::CRDoubleMatrix::redistribute(), oomph::MatrixVectorProduct::setup(), oomph::DenseLU::solve(), oomph::SuperLUSolver::solve(), oomph::MumpsSolver::solve(), oomph::CG< MATRIX >::solve_helper(), oomph::BiCGStab< MATRIX >::solve_helper(), oomph::GMRES< MATRIX >::solve_helper(), oomph::AugmentedProblemGMRES::solve_helper(), oomph::HelmholtzGMRESMG< MATRIX >::solve_helper(), oomph::HelmholtzFGMRESMG< MATRIX >::solve_helper(), and oomph::SuperLUSolver::solve_transpose().

◆ first_row() [1/2]

unsigned oomph::LinearAlgebraDistribution::first_row ( ) const
inline

access function for the first row on this processor. If not distributed then this is just zero.

Definition at line 261 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(), build(), oomph::CRDoubleMatrixHelpers::concatenate(), oomph::DoubleVectorHelpers::concatenate(), 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::DistributableLinearAlgebraObject::first_row(), oomph::Problem::global_dof_pt(), oomph::HypreInterface::hypre_solve(), oomph::BlockPreconditioner< MATRIX >::index_in_block(), oomph::BlockPreconditioner< MATRIX >::internal_dof_number(), oomph::BlockPreconditioner< MATRIX >::internal_index_in_dof(), operator==(), oomph::Problem::parallel_sparse_assemble(), oomph::PitchForkHandler::PitchForkHandler(), rank_of_global_row(), oomph::DoubleMultiVector::redistribute(), oomph::DoubleVector::redistribute(), oomph::CRDoubleMatrix::redistribute(), and oomph::CRDoubleMatrix::sparse_indexed_output_with_offset().

◆ first_row() [2/2]

unsigned oomph::LinearAlgebraDistribution::first_row ( const unsigned &  p) const
inline

access function for the first row on the p-th processor

Definition at line 289 of file linear_algebra_distribution.h.

References Comm_pt, Distributed, and First_row.

◆ first_row_vector()

Vector<unsigned> oomph::LinearAlgebraDistribution::first_row_vector ( ) const
inline

return the first_row Vector

Definition at line 404 of file linear_algebra_distribution.h.

References First_row.

Referenced by build().

◆ global_to_local_row_map()

unsigned oomph::LinearAlgebraDistribution::global_to_local_row_map ( const unsigned &  global_i) const
inline

return the local index corresponding to the global index

Definition at line 365 of file linear_algebra_distribution.h.

References Nrow, and nrow_local().

◆ nrow()

unsigned oomph::LinearAlgebraDistribution::nrow ( ) const
inline

◆ nrow_local() [1/2]

unsigned oomph::LinearAlgebraDistribution::nrow_local ( ) const
inline

access function for the num of local rows on this processor. If no MPI then Nrow is returned.

Definition at line 193 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::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::create_new_linear_algebra_distribution(), oomph::Problem::get_hessian_vector_products(), oomph::Problem::global_dof_pt(), global_to_local_row_map(), oomph::HypreInterface::hypre_solve(), oomph::BlockPreconditioner< MATRIX >::index_in_block(), oomph::BlockPreconditioner< MATRIX >::internal_dof_number(), oomph::BlockPreconditioner< MATRIX >::internal_index_in_dof(), oomph::Problem::newton_solve(), oomph::Problem::newton_solve_continuation(), oomph::DistributableLinearAlgebraObject::nrow_local(), operator==(), oomph::Problem::parallel_sparse_assemble(), oomph::PitchForkHandler::PitchForkHandler(), rank_of_global_row(), oomph::DoubleMultiVector::redistribute(), oomph::DoubleVector::redistribute(), oomph::CRDoubleMatrix::redistribute(), oomph::BlockPitchForkLinearSolver::resolve(), oomph::Problem::restore_dof_values(), oomph::BlockPitchForkLinearSolver::solve(), oomph::PitchForkHandler::solve_block_system(), oomph::PitchForkHandler::solve_full_system(), oomph::Problem::store_current_dof_values(), and oomph::PitchForkHandler::~PitchForkHandler().

◆ nrow_local() [2/2]

unsigned oomph::LinearAlgebraDistribution::nrow_local ( const unsigned &  p) const
inline

access function for the num of local rows on this processor. If no MPI the nrow is returned

Definition at line 222 of file linear_algebra_distribution.h.

References Comm_pt, Distributed, Nrow, and Nrow_local.

◆ nrow_local_vector()

Vector<unsigned> oomph::LinearAlgebraDistribution::nrow_local_vector ( ) const
inline

return the nrow_local Vector

Definition at line 398 of file linear_algebra_distribution.h.

References Nrow_local.

Referenced by build().

◆ operator!=()

bool oomph::LinearAlgebraDistribution::operator!= ( const LinearAlgebraDistribution other_dist) const
inline

!= operator

Definition at line 355 of file linear_algebra_distribution.h.

◆ operator=()

void oomph::LinearAlgebraDistribution::operator= ( const LinearAlgebraDistribution old_dist)
inline

Assignment Operator.

Definition at line 137 of file linear_algebra_distribution.h.

References build().

◆ operator==()

bool oomph::LinearAlgebraDistribution::operator== ( const LinearAlgebraDistribution other_dist) const

== Operator

operator==

Definition at line 268 of file linear_algebra_distribution.cc.

References Comm_pt, communicator_pt(), distributed(), Distributed, first_row(), First_row, i, nrow(), Nrow, nrow_local(), and Nrow_local.

◆ rank_of_global_row()

unsigned oomph::LinearAlgebraDistribution::rank_of_global_row ( const unsigned  i) const
inline

Friends And Related Function Documentation

◆ operator<<

std::ostream& operator<< ( std::ostream &  stream,
LinearAlgebraDistribution dist 
)
friend

<< operator

Definition at line 318 of file linear_algebra_distribution.cc.

Member Data Documentation

◆ Comm_pt

OomphCommunicator* oomph::LinearAlgebraDistribution::Comm_pt
private

the pointer to the MPI communicator object in this distribution

Definition at line 425 of file linear_algebra_distribution.h.

Referenced by build(), built(), clear(), communicator_pt(), first_row(), nrow_local(), operator==(), and ~LinearAlgebraDistribution().

◆ Distributed

bool oomph::LinearAlgebraDistribution::Distributed
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 422 of file linear_algebra_distribution.h.

Referenced by build(), distributed(), first_row(), nrow_local(), and operator==().

◆ First_row

Vector<unsigned> oomph::LinearAlgebraDistribution::First_row
private

the first row on this processor

Definition at line 417 of file linear_algebra_distribution.h.

Referenced by build(), clear(), first_row(), first_row_vector(), and operator==().

◆ Nrow

unsigned oomph::LinearAlgebraDistribution::Nrow
private

the number of global rows

Definition at line 411 of file linear_algebra_distribution.h.

Referenced by build(), clear(), global_to_local_row_map(), nrow(), nrow_local(), and operator==().

◆ Nrow_local

Vector<unsigned> oomph::LinearAlgebraDistribution::Nrow_local
private

the number of local rows on the processor

Definition at line 414 of file linear_algebra_distribution.h.

Referenced by build(), clear(), nrow_local(), nrow_local_vector(), and operator==().


The documentation for this class was generated from the following files: