# Distributed vectors and matrices

Distributed vectors and matrices are implemented in Seldon. Unitary tests are present in files src/Program/Unit/Algebra/distributed_vector.cc, src/Program/Unit/Algebra/distributed_matrix_test.cc and src/Program/Unit/Algebra/distributed_solver.cc for respectively distributed vectors, distributed matrices and distributed solvers. In Montjoie, distributed block-diagonal matrices are implemented. Below, we show a basic example with these matrices

```// for the constructor, you use the local number of rows
DistributedMatrix<double, General, BlockDiagRow> A(n, n);

// then you need to initialise the arrays used for parallel functions (e.g. matrix-vector product)
A.Init(nglob, &glob_number, &OverlapRow, &original_proc,
n, 1, &SharingProc, &SharingRows, MPI::COMM_WORLD);

// example with two processors
// processor 0 owns rows 0, 2, 5, 6, 8
// processor 1 owns rows 1, 2, 3, 4, 6, 7
// block 0 : (0, 1, 4)
// block 1 : (2, 3, 6)
// block 2 : (5, 7, 8)
Vector<IVect> pattern(3), proc(3);
for (int k = 0; k < 3; k++)
{
pattern(k).Reallocate(3);
proc(k).Reallocate(3);
}

if (MPI::COMM_WORLD.Get_rank() == 0)
{
// block 0 : (0, 1, 4)
// rows 1 and 4 are owned by processor 1, we put global numbers
pattern(0)(0) = 0; pattern(0)(1) = 1; pattern(0)(2) = 4;
proc(0)(0) = 0; proc(0)(1) = 1; proc(0)(2) = 1;

// block 1 : (2, 3, 6)
// local numbers of rows 2, 6 : 1, 3
pattern(1)(0) = 1; pattern(1)(1) = 3; pattern(1)(2) = 3;
proc(1)(0) = 0; proc(1)(1) = 1; proc(1)(2) = 0;

// block 2 : (5, 7, 8)
// local numbers of rows 5, 8 : 2, 4
pattern(2)(0) = 2; pattern(2)(1) = 7; pattern(2)(2) = 4;
proc(2)(0) = 0; proc(2)(1) = 1; proc(2)(2) = 0;
}
else
{
// block 0 : (0, 1, 4)
// local numbers of rows 1, 4 : 0 3
pattern(0)(0) = 0; pattern(0)(1) = 0; pattern(0)(2) = 3;
proc(0)(0) = 0; proc(0)(1) = 1; proc(0)(2) = 1;

// block 1 : (2, 3, 6)
// local numbers of rows 2, 3, 6 : 1, 2, 4
pattern(1)(0) = 1; pattern(1)(1) = 2; pattern(1)(2) = 4;
proc(1)(0) = 1; proc(1)(1) = 1; proc(1)(2) = 1;

// block 2 : (5, 7, 8)
// local number of row 7 : 5
pattern(2)(0) = 5; pattern(2)(1) = 5; pattern(2)(2) = 8;
proc(2)(0) = 0; proc(2)(1) = 1; proc(2)(2) = 0;
}

A.SetPattern(pattern, proc);

// once the pattern has been defined, you can add non-zero entries to the matrix

// if the column is distant => AddDistantInteraction

// another solution is to convert a distributed sparse matrix into a block-diagonal one
// in that case SetPattern is not needed, you just call ConvertToBlockDiagonal

// once the matrix is filled, you can compute its inverse
GetInverse(A);
```

## Methods for block-diagonal distributed matrices :

 Clear clears the matrix Init sets pointers to the arrays containing global row numbers and overlapped rows GetCommunicator returns the MPI communicator associated with the matrix GetGlobalM returns the global number of rows GetNodlScalar returns the number of rows for a scalar unknown GetNbScalarUnknowns returns the number of scalar unknowns GetProcNumber returns a pointer to the array containing the processor number associated with each row. GetGlobalRowNumber returns local to global numbering GetOverlapRowNumber returns the array containing the numbers of rows already handled by an another processor GetOverlapProcNumber returns the array containing the processor numbers of rows already handled by an another processor GetProcessorSharingRows returns the list of processors that share rows with the current processor GetSharingRowNumbers returns the list of rows shared with each processor SetPattern provides the structure of the block-diagonal matrix AddDistantInteraction adds a value to A(i, j) where i is local and j global AddRowDistantInteraction adds a value to A(i, j) where i is global and j local IndexGlobal returns the global row number of row j of block i WriteText writes the distributed matrix on several files (one per processor) in ascii format IsReadyForMltAdd returns true if the structure is ready to perform a matrix vector without preparation Assemble assembles the block-diagonal matrix

## Functions for distributed block-diagonal matrices :

 Mlt performs a matrix-vector product MltAdd performs a matrix-vector product Add adds a distributed matrix to another one ConvertToSparse converts a block-diagonal matrix to a sparse matrix ConvertToBlockDiagonal converts a sparse matrix to a block-diagonal matrix Copy Copies/Converts a distributed matrix into another one GetInverse replaces a block-diagonal matrix by its inverse GetCholesky replaces a block-diagonal matrix by its Cholesky factorisation SolveCholesky solves L x = b or L^T x = b for a block-diagonal matrix MltCholesky performs product with Cholesky factor L (or its transpose) EraseRow erases rows in the distributed matrix EraseCol erases columns in the distributed matrix

### GetCommunicator

#### Syntax

```  MPI::Comm& GetCommunicator();
```

This method returns the MPI communicator used to distribute the matrix.

### GetGlobalM

#### Syntax

``` int GetGlobalM() const;
```

This function returns the global number of rows.

#### Example :

```// default constructor
DistributedMatrix<double, General, BlockDiagRow> A;

// after filling nglob, global_row, overlap_row and overlap_proc
// calling Init in order to provide these datas
A.Init(nglob, &glob_number, &OverlapRow, &original_proc,
Nvol, nb_u, &ProcShared, &SharedRows, MPI::COMM_WORLD);

nglob = A.GetGlobalM();
```

### GetNodlScalar

#### Syntax

``` int GetNodlScalar() const;
```

This method returns the local number of rows for each scalar components. Several components can be specified when Init is called.

#### Example :

```// default constructor
DistributedMatrix<double, General, BlockDiagRow> A;

// after filling nglob, global_row, overlap_row and overlap_proc
// calling Init in order to provide these datas
A.Init(nglob, &glob_number, &OverlapRow, &original_proc,
Nvol, nb_u, &ProcShared, &SharedRows, MPI::COMM_WORLD);

Nvol = A.GetNodlScalar();
```

### GetNbScalarUnknowns

#### Syntax

``` int GetNbScalarUnknowns() const;
```

This method returns the local number of scalar components. Several components can be specified when Init is called.

#### Example :

```// default constructor
DistributedMatrix<double, General, BlockDiagRow> A;

// after filling nglob, global_row, overlap_row and overlap_proc
// calling Init in order to provide these datas
A.Init(nglob, &glob_number, &OverlapRow, &original_proc,
Nvol, nb_u, &ProcShared, &SharedRows, MPI::COMM_WORLD);

nb_u = A.GetNbScalarUnknowns();
```

#### Syntax

``` void AddDistantInteraction(int i, int jglob, int proc, const T& val);
```

This member function adds val for the local row i, and the global row jglob, proc being the processor that treats the global row jglob.

#### Example :

```// default constructor
DistributedMatrix<double, General, BlockDiagRow> A;

// after filling nglob, global_row, overlap_row and overlap_proc
// calling Init in order to provide these datas
A.Init(nglob, &global_row, &overlap_row, &overlap_proc,
Nvol, nb_u, &ProcShared, &SharedRows, MPI::COMM_WORLD);

// for local interaction, you have to use AddInteraction

// when the column is distant (ie located on another processor), you have to use AddDistantInteraction
// jglob : global column number
// proc : distant processor
// val : value to add
```

#### Syntax

``` void AddRowDistantInteraction(int iglob, int j, int proc, const T& val);
```

This member function adds val for the global row iglob, and the local column j, proc being the processor that treats the global row iglob.

#### Example :

```// default constructor
DistributedMatrix<double, General, BlockDiagRow> A;

// after filling nglob, global_row, overlap_row and overlap_proc
// calling Init in order to provide these datas
A.Init(nglob, &global_row, &overlap_row, &overlap_proc,
Nvol, nb_u, &ProcShared, &SharedRows, MPI::COMM_WORLD);

// for a local interaction, you have to use AddInteraction

// when the column is distant (ie located on another processor), you have to use AddDistantInteraction
// i : local row number
// jglob : global column number
// proc : distant processor
// val : value to add
// for the global matrix, it means that A_{global_row(i), jglob} is incremented with val

// when the row is distant, you have to use AddRowDistantInteraction
// iglob : global row number
// j : local column number
// proc : distant processor
// val : value to add
// for the global matrix, it means that A_{iglob, global_row(j)} is incremented with val

```

#### Syntax

``` bool IsReadyForMltAdd() const
```

This member function returns true if the structure is ready to perform a matrix-vector product. This function is not really useful since if the user performs a matrix-vector product with a structure that is not ready, the program will create all the arrays needed to perform the matrix-vector product. The function IsReadyForMltAdd can be used for informative purposes, for example to avoid taking into account the cost of constructing all these arrays (they induce MPI communications) if the structure was not ready.

#### Example :

```// default constructor
DistributedMatrix<double, General, BlockDiagRow> A;

// after filling nglob, global_row, overlap_row and overlap_proc
// calling Init in order to provide these datas
A.Init(nglob, &global_row, &overlap_row, &overlap_proc,
Nvol, nb_u, &ProcShared, &SharedRows, MPI::COMM_WORLD);

// for a local interaction, you have to use AddInteraction

// when the column is distant (ie located on another processor), you have to use AddDistantInteraction
// i : local row number
// jglob : global column number
// proc : distant processor
// val : value to add
// for the global matrix, it means that A_{global_row(i), jglob} is incremented with val

// when the row is distant, you have to use AddRowDistantInteraction
// iglob : global row number
// j : local column number
// proc : distant processor
// val : value to add
// for the global matrix, it means that A_{iglob, global_row(j)} is incremented with val

// the structure will be said "ready" for a matrix-vector operation
// if a first matrix-vector product has been performed
```

### Clear

#### Syntax

``` void Clear();
```

This function clears the matrix.

#### Example :

```// default constructor
DistributedMatrix<double, General, BlockDiagRow> A;

// the matrix is erased
A.Clear();
```

### GetGlobalRowNumber

#### Syntax

``` IVect& GetGlobalRowNumber() const;
```

This function returns the global row numbers (local to global numbering).

#### Example :

```// default constructor
DistributedMatrix<double, General, BlockDiagRow> A;

// after filling nglob, global_row, overlap_row and overlap_proc
// calling Init in order to provide these datas
A.Init(nglob, &glob_number, &OverlapRow, &original_proc,
n, 1, &SharingProc, &SharingRows, MPI::COMM_WORLD);

IVect& global_row = A.GetGlobalRowNumber();
```

### GetOverlapRowNumber

#### Syntax

``` IVect& GetOverlapRowNumber() const;
```

This function returns the overlapped row numbers, i.e. rows that are already treated by another processor.

#### Example :

```// default constructor
DistributedMatrix<double, General, BlockDiagRow> A;

// after filling nglob, global_row, overlap_row and overlap_proc
// calling Init in order to provide these datas
A.Init(nglob, &glob_number, &OverlapRow, &original_proc,
n, 1, &SharingProc, &SharingRows, MPI::COMM_WORLD);

IVect& overlap_row = A.GetOverlapRowNumber();
```

### GetOverlapProcNumber

#### Syntax

``` IVect& GetOverlapProcNumber() const;
```

This function returns the processor numbers associated with overlapped row numbers, i.e. rows that are already treated by another processor.

#### Example :

```// default constructor
DistributedMatrix<double, General, BlockDiagRow> A;

// after filling nglob, global_row, overlap_row and overlap_proc
// calling Init in order to provide these datas
A.Init(nglob, &glob_number, &OverlapRow, &original_proc,
n, 1, &SharingProc, &SharingRows, MPI::COMM_WORLD);

IVect& original_proc = A.GetOverlaProcNumber();
```

### GetProcessorSharingRows

#### Syntax

``` IVect& GetProcessorSharingRows() const;
```

This function returns the processor numbers that share rows with the current processor.

#### Example :

```// default constructor
DistributedMatrix<double, General, BlockDiagRow> A;

// after filling nglob, global_row, overlap_row and overlap_proc
// calling Init in order to provide these datas
A.Init(nglob, &glob_number, &OverlapRow, &original_proc,
n, 1, &SharingProc, &SharingRows, MPI::COMM_WORLD);

IVect& SharingProc = A.GetProcessorSharingRows();
```

### GetSharingRowNumbers

#### Syntax

``` Vector<IVect>& GetSharingRowNumbers() const;
```

This function returns the row numbers that are shared with other processors.

#### Example :

```// default constructor
DistributedMatrix<double, General, BlockDiagRow> A;

// after filling nglob, global_row, overlap_row and overlap_proc
// calling Init in order to provide these datas
A.Init(nglob, &glob_number, &OverlapRow, &original_proc,
n, 1, &SharingProc, &SharingRows, MPI::COMM_WORLD);

IVect& SharingProc = A.GetSharingRowNumbers();
```

### WriteText

#### Syntax

```  void WriteText(string) const;
void WriteText(ostream&) const;
```

This member function writes the matrix on files or output streams in text format. Each processor writes its part of the matrix on a different file (ending with _P0.dat, _P1.dat, etc).

#### Example :

```// default constructor
DistributedMatrix<double, General, BlockDiagRow> A;

// you can write the matrix in files
// the matrix here will be written in files mat_P0.dat, mat_P1.dat, etc
A.WriteText("mat.dat");
```

#### Syntax

```  void MltAdd(const T& alpha, const DistributedMatrix<T>& A, const Vector<T>& x, const T& beta, Vector<T>& y, bool assemble);
void MltAdd(const T& alpha, const SeldonTranspose&, const DistributedMatrix<T>& A, const Vector<T>& x, const T& beta, Vector<T>& y, bool assemble);
```

This function can perform matrix-vector product (with the original matrix or with its transpose). The matrix-matrix product is currently not implemented. If the communicator contains only one processor, the sequential function MltAdd (with the class Matrix) will be called, otherwise an error will be displayed during the execution.

#### Example :

```// default constructor
DistributedMatrix<double, General, BlockDiagRow> A;

// global row numbers are provided with Init (glob_number)
A.Init(nglob, &glob_number, &OverlapRow, &original_proc,
Nvol, nb_u, &ProcShared, &SharedRows, MPI::COMM_WORLD);

// once the matrix is constructed, matrix-vector products can be performed
// each vector is stored as a usual Seldon vector (class Vector)
// each processor is assumed to store the values of the vector for the
// global rows given during the construction of the matrix (here, the array glob_number)
Vector<double> x(n), y(n), z;

// MltAdd will assume that the values of x for rows that are shared between processors are the same
// which is here not the case when calling FillRand
x.FillRand();
y.FillRand();
z = y;

// a solution is to call AssembleVector to ensure this property
AssembleVector(x, MPI::SUM, ProcShared, SharedRows, MPI::COMM_WORLD, Nvol, nb_u, 23);

double alpha = 2.5, beta = 0.4;

// matrix-vector product y = beta y + alpha A x
// the default value of the last argument (assemble) is true
// the result y is assembled and will contain the same result as a sequential matrix-vector product,
// the values being distributed in the different processors

// assemble = true, induces a final call to the function AssembleVector (to add values of shared rows)
// if the user wants to perform a matrix-vector product without performing this final assembling step, he puts false in the last optional argument
MltAdd(alpha, A, x, beta, z, false);

// and perform the assembling phase later
AssembleVector(z, MPI::SUM, ProcShared, SharedRows, MPI::COMM_WORLD, Nvol, nb_u, 23);

// MltAdd can be used to perform a matrix-vector with the transpose matrix or the transpose conjugate
MltAdd(alpha, SeldonTrans, A, x, beta, y);
```

#### Syntax

```  void Add(const T& alpha, const DistributedMatrix<T>& A, DistributedMatrix<T>& B);
```

This function adds a distributed matrix multiplied by a scalar to another distributed matrix. The two matrices are assumed to have the same global row numbers.

#### Example :

```// default constructor
DistributedMatrix<double, General, BlockDiagRow> A;

// global row numbers are provided with Init (glob_number)
A.Init(nglob, &glob_number, &OverlapRow, &original_proc,
Nvol, nb_u, &ProcShared, &SharedRows, MPI::COMM_WORLD);

// a way to have the same global row numbers
// is to use Init with the matrix A
DistributedMatrix<double, General, BlockDiagRow> B;

B.Init(A);

double alpha = 2.1;
// we compute B = B + alpha A
```

### Copy

#### Syntax

```  void Copy(const DistributedMatrix<T>& A, DistributedMatrix<T>&);
```

This function converts a distributed sparse matrix into another one. The format can be different (RowSparse, ColSparse, ArrayRowSparse, etc), the local part of the matrix is converted by calling the appropriate Copy function (located in the file Matrix_Conversions.cxx), the distributed part is always stored with the same format.

### EraseRow

#### Syntax

``` void EraseRow(IVect& num, DistributedMatrix& A);
```

This function erases some rows (local row numbers are provided in num) of the matrix A.

#### Example :

```// default constructor
DistributedMatrix<double, General, BlockDiagRow> A;

// after filling nglob, global_row, overlap_row and overlap_proc
// calling Init in order to provide these datas
A.Init(nglob, &global_row, &overlap_row, &overlap_proc, MPI::COMM_WORLD);

// then you can erase some rows
IVect num(2); num(0) = 2; num(1) = 5;
EraseRow(num, A);
```

### EraseCol

#### Syntax

``` void EraseCol(IVect& num, DistributedMatrix& A);
```

This function erases some columns (local column numbers are provided in num) of the matrix A.

#### Example :

```// default constructor
DistributedMatrix<double, General, BlockDiagRow> A;

// after filling nglob, global_row, overlap_row and overlap_proc
// calling Init in order to provide these datas
A.Init(nglob, &global_row, &overlap_row, &overlap_proc, MPI::COMM_WORLD);

// then you can erase some columns
IVect num(2); num(0) = 4; num(1) = 9;
EraseCol(num, A);
```

### Init

#### Syntax

``` void Init(int nglob, IVect* global_row, IVect* overlap_row, IVect* overlap_proc,
int Nvol, int nb_u, IVect* SharingProc, Vector<IVect>* SharingRows, MPI::Comm& comm);

void Init(const DistributedMatrix<T>&);

```

The function Init must be called after the constructor in order to provide the global row numbers (array global_row), the rows already treated by another processor (array overlap_row), the processors that treat originally these rows (array overlap_proc) and the MPI communicator. Nvol is the number of rows for a scalar component, nb_u the number of scalar components, SharingProc stores the processors that share rows with current processor, and SharingRows the local row numbers that are shared with processor SharingProc(i). These numbers are assumed to be sorted such that shared row numbers with processor j in processor i correspond to row numbers shared with processor i in processor j.

#### Example :

```// on each processor, you specify rows that are already treated by another processor
// for example if the processor 0 handles rows [0, 3, 5, 6] and
// processor 1 the rows [1 2, 4, 5, 7], the row 5 is already treated by
// processor 0. Of course, if each row is treated by a processor and only one
// this array should be left empty
IVect OverlapRow;
int nglob = 8;
int n = 4;
if (MPI::COMM_WORLD.Get_rank() == 1)
{
n = 5;
OverlapRow.Reallocate(1);
// be careful because OverlapRow stores local numbers
// here the global row 5 has a local number equal to 3 on processor 1
OverlapRow(0) = 3;
}

// for distributed matrices
// we also need to know the global row numbers
// and for each overlapped row, the original processor
IVect glob_number, original_proc;

// for shared rows, all the row numbers shared with processors
IVect SharingProc(1); Vector<IVect> SharingRows(1);
if (MPI::COMM_WORLD.Get_rank() == 0)
{
glob_number.Reallocate(4);
glob_number(0) = 0;
glob_number(1) = 3;
glob_number(2) = 5;
glob_number(3) = 6;

// global row 5 is shared with processor 1 (local number = 2)
SharingProc(0) = 1;
SharingRows(0).PushBack(2);
}
else
{
glob_number.Reallocate(5);
glob_number(0) = 1;
glob_number(1) = 2;
glob_number(2) = 4;
glob_number(3) = 5;
glob_number(4) = 7;

// row 5 is already handled by processor 0
original_proc.Reallocate(1);
original_proc(0) = 0;

// here SharingRows is similar to OverlapRow because all shared dofs of this processor are already treated by processor 0
// global row 5 is shared with processor 0 (local number = 3)
SharingProc(0) = 0;
SharingRows(0).PushBack(3);
}

// for the constructor, you use the local number of rows
DistributedMatrix<double, General, BlockDiagRow> A(n, n);

// then you need to initialise the arrays used for parallel functions (e.g. matrix-vector product)
A.Init(nglob, &glob_number, &OverlapRow, &original_proc,
n, 1, &SharingProc, &SharingRows, MPI::COMM_WORLD);
```

Another example does not construct these arrays but calls the method Init of a sparse distributed matrix (not block-diagonal) to construct them automatically.

```// local number of rows
int n = 4;
if (MPI::COMM_WORLD.Get_rank() == 1)
n = 5;

// each processor constructs its global numbers
// rows [0, 3, 5, 6] for processor 0
// rows [1 2, 4, 5, 7] for processor 1
// you notice that global row 5 is here shared by the two processors
IVect glob_number(n);
if (MPI::COMM_WORLD.Get_rank() == 0)
{
glob_number(0) = 0;
glob_number(1) = 3;
glob_number(2) = 5;
glob_number(3) = 6;
}
else
{
glob_number(0) = 1;
glob_number(1) = 2;
glob_number(2) = 4;
glob_number(3) = 5;
glob_number(4) = 7;
}

// for the constructor, you use the local number of rows
DistributedMatrix<double, General, ArrayRowSparse> Asp(n, n);
DistributedMatrix<double, General, BlockDiagRow> A(n, n);

// arrays overlap_row, original_proc, SharingProc, SharingRows
// are output arrays of Init here, they are constructed with the given global numbers
IVect overlap_row, original_proc;
IVect SharingProc;
Vector<IVect> SharingRows;
// these arrays are constructed by Init of a distributed sparse matrix (Asp here)
Asp.Init(glob_number, overlap_row, original_proc, SharingProc, SharingRows, MPI::COMM_WORLD);

// once these arrays have been constructed, you can use them for other distributed matrices that have the same distribution.
DistributedMatrix<double, General, BlockDiagRow> B(n, n);

B.Init(nglob, &glob_number, &OverlapRow, &original_proc,
n, 1, &SharingProc, &SharingRows, MPI::COMM_WORLD);

// another solution is to call Init with the original distributed matrix
A.Init(Asp);
```

### GetProcNumber

#### Syntax

``` int* GetProcNumber() const
```

This method returns a pointer storing the processor number for each "row". Here, we take into account all the rows that are coupled with the rows owned by the current processor.

#### Example :

```// for the constructor, you use the local number of rows
DistributedMatrix<double, General, BlockDiagRow> A(n, n);

// then you need to initialise the arrays used for parallel functions (e.g. matrix-vector product)
A.Init(nglob, &glob_number, &OverlapRow, &original_proc,
n, 1, &SharingProc, &SharingRows, MPI::COMM_WORLD);

A.SetPattern(pattern, proc);

int* proc_row = A.GetProcNumber();
```

### SetPattern

#### Syntax

```  void SetPattern(const Vector<IVect>& NumDof_Blocks);
void SetPattern(const Vector<IVect>& NumDof_Blocks, const Vector<IVect>& ProcDofs);
```

This method specifies the profile of the block-diagonal matrix. In the first syntax, the blocks are not shared between processors, and the user provides the local row numbers for each block. In the second syntax, the blocks may be shared between processors, the second argument gives the processor number for each row of the block. If the processor is distant (different from the current processor), the row number is the global row number, otherwise it is the local row number.

#### Example :

```// for the constructor, you use the local number of rows
DistributedMatrix<double, General, BlockDiagRow> A(n, n);

// then you need to initialise the arrays used for parallel functions (e.g. matrix-vector product)
A.Init(nglob, &glob_number, &OverlapRow, &original_proc,
n, 1, &SharingProc, &SharingRows, MPI::COMM_WORLD);

// example with two processors
// processor 0 owns rows 0, 2, 5, 6, 8
// processor 1 owns rows 1, 2, 3, 4, 6, 7
// block 0 : (0, 1, 4)
// block 1 : (2, 3, 6)
// block 2 : (5, 7, 8)
Vector<IVect> pattern(3), proc(3);
for (int k = 0; k < 3; k++)
{
pattern(k).Reallocate(3);
proc(k).Reallocate(3);
}

if (MPI::COMM_WORLD.Get_rank() == 0)
{
// block 0 : (0, 1, 4)
// rows 1 and 4 are owned by processor 1, we put global numbers
pattern(0)(0) = 0; pattern(0)(1) = 1; pattern(0)(2) = 4;
proc(0)(0) = 0; proc(0)(1) = 1; proc(0)(2) = 1;

// block 1 : (2, 3, 6)
// local numbers of rows 2, 6 : 1, 3
pattern(1)(0) = 1; pattern(1)(1) = 3; pattern(1)(2) = 3;
proc(1)(0) = 0; proc(1)(1) = 1; proc(1)(2) = 0;

// block 2 : (5, 7, 8)
// local numbers of rows 5, 8 : 2, 4
pattern(2)(0) = 2; pattern(2)(1) = 7; pattern(2)(2) = 4;
proc(2)(0) = 0; proc(2)(1) = 1; proc(2)(2) = 0;
}
else
{
// block 0 : (0, 1, 4)
// local numbers of rows 1, 4 : 0 3
pattern(0)(0) = 0; pattern(0)(1) = 0; pattern(0)(2) = 3;
proc(0)(0) = 0; proc(0)(1) = 1; proc(0)(2) = 1;

// block 1 : (2, 3, 6)
// local numbers of rows 2, 3, 6 : 1, 2, 4
pattern(1)(0) = 1; pattern(1)(1) = 2; pattern(1)(2) = 4;
proc(1)(0) = 1; proc(1)(1) = 1; proc(1)(2) = 1;

// block 2 : (5, 7, 8)
// local number of row 7 : 5
pattern(2)(0) = 5; pattern(2)(1) = 5; pattern(2)(2) = 8;
proc(2)(0) = 0; proc(2)(1) = 1; proc(2)(2) = 0;
}

A.SetPattern(pattern, proc);
```

### IndexGlobal

#### Syntax

```  int IndexGlobal(int i, int j)
```

This method returns the global row number of row j of block i.

#### Example :

```// for the constructor, you use the local number of rows
DistributedMatrix<double, General, BlockDiagRow> A(n, n);

// then you need to initialise the arrays used for parallel functions (e.g. matrix-vector product)
A.Init(nglob, &glob_number, &OverlapRow, &original_proc,
n, 1, &SharingProc, &SharingRows, MPI::COMM_WORLD);

// example with two processors
// processor 0 owns rows 0, 2, 5, 6, 8
// processor 1 owns rows 1, 2, 3, 4, 6, 7
// block 0 : (0, 1, 4)
// block 1 : (2, 3, 6)
// block 2 : (5, 7, 8)
Vector<IVect> pattern(3), proc(3);
for (int k = 0; k < 3; k++)
{
pattern(k).Reallocate(3);
proc(k).Reallocate(3);
}

if (MPI::COMM_WORLD.Get_rank() == 0)
{
// block 0 : (0, 1, 4)
// rows 1 and 4 are owned by processor 1, we put global numbers
pattern(0)(0) = 0; pattern(0)(1) = 1; pattern(0)(2) = 4;
proc(0)(0) = 0; proc(0)(1) = 1; proc(0)(2) = 1;

// block 1 : (2, 3, 6)
// local numbers of rows 2, 6 : 1, 3
pattern(1)(0) = 1; pattern(1)(1) = 3; pattern(1)(2) = 3;
proc(1)(0) = 0; proc(1)(1) = 1; proc(1)(2) = 0;

// block 2 : (5, 7, 8)
// local numbers of rows 5, 8 : 2, 4
pattern(2)(0) = 2; pattern(2)(1) = 7; pattern(2)(2) = 4;
proc(2)(0) = 0; proc(2)(1) = 1; proc(2)(2) = 0;
}
else
{
// block 0 : (0, 1, 4)
// local numbers of rows 1, 4 : 0 3
pattern(0)(0) = 0; pattern(0)(1) = 0; pattern(0)(2) = 3;
proc(0)(0) = 0; proc(0)(1) = 1; proc(0)(2) = 1;

// block 1 : (2, 3, 6)
// local numbers of rows 2, 3, 6 : 1, 2, 4
pattern(1)(0) = 1; pattern(1)(1) = 2; pattern(1)(2) = 4;
proc(1)(0) = 1; proc(1)(1) = 1; proc(1)(2) = 1;

// block 2 : (5, 7, 8)
// local number of row 7 : 5
pattern(2)(0) = 5; pattern(2)(1) = 5; pattern(2)(2) = 8;
proc(2)(0) = 0; proc(2)(1) = 1; proc(2)(2) = 0;
}

A.SetPattern(pattern, proc);

if (MPI::COMM_WORLD.Get_rank() == 0)
{
// it should return 6 (global row number of row 2 of block 1)
cout << A.IndexGlobal(1, 2) << endl;
}
```

### Assemble

#### Syntax

```  void Assemble()
```

This method assembles the block-diagonal matrix (values shared between processors are added).

#### Example :

```// for the constructor, you use the local number of rows
DistributedMatrix<double, General, BlockDiagRow> A(n, n);

// then you need to initialise the arrays used for parallel functions (e.g. matrix-vector product)
A.Init(nglob, &glob_number, &OverlapRow, &original_proc,
n, 1, &SharingProc, &SharingRows, MPI::COMM_WORLD);

// example with two processors
// processor 0 owns rows 0, 2, 5, 6, 8
// processor 1 owns rows 1, 2, 3, 4, 6, 7
// block 0 : (0, 1, 4)
// block 1 : (2, 3, 6)
// block 2 : (5, 7, 8)
Vector<IVect> pattern(3), proc(3);
for (int k = 0; k < 3; k++)
{
pattern(k).Reallocate(3);
proc(k).Reallocate(3);
}

if (MPI::COMM_WORLD.Get_rank() == 0)
{
// block 0 : (0, 1, 4)
// rows 1 and 4 are owned by processor 1, we put global numbers
pattern(0)(0) = 0; pattern(0)(1) = 1; pattern(0)(2) = 4;
proc(0)(0) = 0; proc(0)(1) = 1; proc(0)(2) = 1;

// block 1 : (2, 3, 6)
// local numbers of rows 2, 6 : 1, 3
pattern(1)(0) = 1; pattern(1)(1) = 3; pattern(1)(2) = 3;
proc(1)(0) = 0; proc(1)(1) = 1; proc(1)(2) = 0;

// block 2 : (5, 7, 8)
// local numbers of rows 5, 8 : 2, 4
pattern(2)(0) = 2; pattern(2)(1) = 7; pattern(2)(2) = 4;
proc(2)(0) = 0; proc(2)(1) = 1; proc(2)(2) = 0;
}
else
{
// block 0 : (0, 1, 4)
// local numbers of rows 1, 4 : 0 3
pattern(0)(0) = 0; pattern(0)(1) = 0; pattern(0)(2) = 3;
proc(0)(0) = 0; proc(0)(1) = 1; proc(0)(2) = 1;

// block 1 : (2, 3, 6)
// local numbers of rows 2, 3, 6 : 1, 2, 4
pattern(1)(0) = 1; pattern(1)(1) = 2; pattern(1)(2) = 4;
proc(1)(0) = 1; proc(1)(1) = 1; proc(1)(2) = 1;

// block 2 : (5, 7, 8)
// local number of row 7 : 5
pattern(2)(0) = 5; pattern(2)(1) = 5; pattern(2)(2) = 8;
proc(2)(0) = 0; proc(2)(1) = 1; proc(2)(2) = 0;
}

A.SetPattern(pattern, proc);

// then the matrix can be assembled
A.Assemble();
```

### Mlt

#### Syntax

```  void Mlt(const DistributedMatrix<T>& A, const Vector<T>& x, Vector<T>& y, bool assemble);
void Mlt(const SeldonTranspose& trans, const DistributedMatrix<T>& A, const Vector<T>& x, Vector<T>& y, bool assemble);
void Mlt(const T& alpha, DistributedMatrix<T>& A);
void Mlt(const DistributedMatrix<T>& A, Vector<T>& x, bool assemble);
```

This function can perform matrix-vector product (with the original matrix or with its transpose) and can be used to multiply the matrix by a scalar. The matrix-matrix product is currently not implemented. If the communicator contains only one processor, the sequential function Mlt (with the class Matrix) will be called, otherwise an error will be displayed during the execution.

#### Example :

```// default constructor
DistributedMatrix<double, General, BlockDiagRow> A;

// global row numbers are provided with Init (glob_number)
A.Init(nglob, &glob_number, &OverlapRow, &original_proc,
Nvol, nb_u, &ProcShared, &SharedRows, MPI::COMM_WORLD);

// SetPattern is called to provide the profile of A
A.SetPattern(NumDofs, ProcDofs);

// for local interaction, you can use AddInteraction

// when the column is distant, you can use AddDistantInteraction

// when the row is distant, you can use AddRowDistantInteraction

// once the matrix is constructed, matrix-vector products can be performed
// each vector is stored as a usual Seldon vector (class Vector)
// each processor is assumed to store the values of the vector for the
// global rows given during the construction of the matrix (here, the array glob_number)
Vector<double> x(n), y(n), z;

// Mlt will assume that the values of x for rows that are shared between processors are the same
// which is here not the case when calling FillRand
x.FillRand();
y.FillRand();
z = y;

// a solution is to call AssembleVector to ensure this property
AssembleVector(x, MPI::SUM, ProcShared, SharedRows, MPI::COMM_WORLD, Nvol, nb_u, 23);

// classical matrix-vector product y = A x
// the default value of the last argument (assemble) is true
// the result y is assembled and will contain the same result as a sequential matrix-vector product,
// the values being distributed in the different processors
Mlt(A, x, y);

// assemble = true, induces a final call to the function AssembleVector (to add values of shared rows)
// if the user wants to perform a matrix-vector product without performing this final assembling step, he puts false in the last optional argument
Mlt(A, x, z, false);

// and perform the assembling phase later
AssembleVector(z, MPI::SUM, ProcShared, SharedRows, MPI::COMM_WORLD, Nvol, nb_u, 23);

// Mlt can be used to perform a matrix-vector with the transpose matrix or the transpose conjugate
Mlt(SeldonTrans, A, x, y);

// Finally Mlt can be used to multiply a distributed matrix by a scalar coefficient
Mlt(-2.1, A);

// the last syntax performs y = A x, on input x is provided and overwritten with y
y = x;
Mlt(A, y);
```

### ConvertToSparse

#### Syntax

```  void ConvertToSparse(const DistributedMatrix<T>& A, DistributedMatrix<T>& B);
```

This function converts a block-diagonal matrix to a sparse distributed matrix. It is equivalent to call the function Copy (which handles other conversions as well).

#### Example :

```// default constructor
DistributedMatrix<double, General, BlockDiagRow> A;

// global row numbers are provided with Init (glob_number)
A.Init(nglob, &glob_number, &OverlapRow, &original_proc,
Nvol, nb_u, &ProcShared, &SharedRows, MPI::COMM_WORLD);

// SetPattern is called to provide the profile of A
A.SetPattern(NumDofs, ProcDofs);

// for local interaction, you can use AddInteraction

// when the column is distant, you can use AddDistantInteraction

// when the row is distant, you can use AddRowDistantInteraction

// once the matrix is constructed, the matrix can be converted
// with Copy or ConvertToSparse
DistributedMatrix<double, General, ArrayRowSparse> B;
ConvertToSparse(A, B);
```

### ConvertToBlockDiagonal

#### Syntax

```  void ConvertToBlockDiagonal(const DistributedMatrix<T>& A, DistributedMatrix<T>& B);
```

This function converts a sparse distributed matrix to a block-diagonal matrix. It is equivalent to call the function Copy (which handles other conversions as well).

#### Example :

```// default constructor
DistributedMatrix<double, General, ArrayRowSparse> A;

// global row numbers are provided with Init (glob_number)
A.Init(nglob, &glob_number, &OverlapRow, &original_proc,
Nvol, nb_u, &ProcShared, &SharedRows, MPI::COMM_WORLD);

// A is allocated with the loca number of rows
A.Reallocate(n, n);

// for local interaction, you can use AddInteraction

// when the column is distant, you can use AddDistantInteraction

// when the row is distant, you can use AddRowDistantInteraction

// once the matrix is constructed, the matrix can be converted
// with Copy or ConvertToBlockDiagonal
DistributedMatrix<double, General, BlockDiagRow> B;
ConvertToBlockDiagonal(A, B);
```

### GetInverse

#### Syntax

```  void GetInverse(DistributedMatrix<T>& B);
```

This function replaces a block-diagonal matrix by its inverse.

#### Example :

```// default constructor
DistributedMatrix<double, General, BlockDiagRow> A;

// global row numbers are provided with Init (glob_number)
A.Init(nglob, &glob_number, &OverlapRow, &original_proc,
Nvol, nb_u, &ProcShared, &SharedRows, MPI::COMM_WORLD);

// SetPattern is called to provide the profile of A
A.SetPattern(NumDofs, ProcDofs);

// for local interaction, you can use AddInteraction

// when the column is distant, you can use AddDistantInteraction

// when the row is distant, you can use AddRowDistantInteraction

// once the matrix is constructed, the inverse can be computed
GetInverse(A);
```

### GetCholesky

#### Syntax

```  void GetCholesky(DistributedMatrix<T>& B);
```

This function replaces a block-diagonal matrix by its Cholesky factorization. It is only available for BlockDiagRowSym storage, since the matrix is assumed to be symmetric definite positive.

#### Example :

```// default constructor
DistributedMatrix<double, Symmetric, BlockDiagRowSym> A;

// global row numbers are provided with Init (glob_number)
A.Init(nglob, &glob_number, &OverlapRow, &original_proc,
Nvol, nb_u, &ProcShared, &SharedRows, MPI::COMM_WORLD);

// SetPattern is called to provide the profile of A
A.SetPattern(NumDofs, ProcDofs);

// for local interaction, you can use AddInteraction

// when the column is distant, you can use AddDistantInteraction

// when the row is distant, you can use AddRowDistantInteraction

// once the matrix is constructed, the Cholesky factorization can be computed
GetCholesky(A);
```

### SolveCholesky

#### Syntax

```  void SolveCholesky(const SeldonTranspose&, const DistributedMatrix<T>& A, Vector<T>& x);
```

This function solves system L x = b or L^T x = b, where A = L L^T (Cholesky factorization). GetCholesky is assumed to have been called previously such that A stores the Cholesky factor L.

#### Example :

```// default constructor
DistributedMatrix<double, Symmetric, BlockDiagRowSym> A;

// global row numbers are provided with Init (glob_number)
A.Init(nglob, &glob_number, &OverlapRow, &original_proc,
Nvol, nb_u, &ProcShared, &SharedRows, MPI::COMM_WORLD);

// SetPattern is called to provide the profile of A
A.SetPattern(NumDofs, ProcDofs);

// for local interaction, you can use AddInteraction

// when the column is distant, you can use AddDistantInteraction

// when the row is distant, you can use AddRowDistantInteraction

// once the matrix is constructed, the Cholesky factorization can be computed
GetCholesky(A);

// and used to solve  L L^T x  = b
Vector<T> x, b(n);
b.FillRand();

x = b;
SolveCholesky(SeldonNoTrans, A, x);
SolveCholesky(SeldonTrans, A, x);
```

### MltCholesky

#### Syntax

```  void MltCholesky(const SeldonTranspose&, const DistributedMatrix<T>& A, Vector<T>& x);
```

This function computes matrix-vector product y = L x or y = L^T x, where A = L L^T (Cholesky factorization). GetCholesky is assumed to have been called previously such that A stores the Cholesky factor L.

#### Example :

```// default constructor
DistributedMatrix<double, Symmetric, BlockDiagRowSym> A;

// global row numbers are provided with Init (glob_number)
A.Init(nglob, &glob_number, &OverlapRow, &original_proc,
Nvol, nb_u, &ProcShared, &SharedRows, MPI::COMM_WORLD);

// SetPattern is called to provide the profile of A
A.SetPattern(NumDofs, ProcDofs);

// for local interaction, you can use AddInteraction

// when the column is distant, you can use AddDistantInteraction

// when the row is distant, you can use AddRowDistantInteraction

// once the matrix is constructed, the Cholesky factorization can be computed
GetCholesky(A);

// and used to compute y = L L^T x
Vector<T> y, x(n);
x.FillRand();

y = x;
MltCholesky(SeldonTrans, A, y);
MltCholesky(SeldonNoTrans, A, y);
```