Elmer FEM solver
Elmer is an open source finite element software for multiphysical problems
 All Classes Files Functions Variables Typedefs Macros Groups Pages
bandmatrix Module Reference

Public Member Functions

subroutine band_zeromatrix (A)
 
subroutine band_zerorow (A, n)
 
subroutine band_addtomatrixelement (A, i, j, value)
 
subroutine band_setmatrixelement (A, i, j, value)
 
real(kind=dp) function band_getmatrixelement (A, i, j)
 
subroutine band_gluelocalmatrix (A, N, Dofs, Indeces, LocalMatrix)
 
subroutine sband_setdirichlet (A, b, n, Value)
 
type(matrix_t) function, pointer band_creatematrix (N, Subband, Symmetric, AllocValues)
 
subroutine band_matrixvectormultiply (A, u, v)
 
subroutine band_matrixvectorprod (u, v, ipar)
 
subroutine band_diagprecondition (u, v, ipar)
 

Member Function/Subroutine Documentation

subroutine bandmatrix::band_addtomatrixelement ( type(matrix_t)  A,
integer  i,
integer  j,
real(kind=dp)  value 
)

Add a given element to the band matrix.

Parameters
aStructure holding matrix
irow number of the matrix element
jcolumn number of the matrix element
valueValue to be added

Referenced by solverutils::addtomatrixelement().

Here is the caller graph for this function:

type(matrix_t) function, pointer bandmatrix::band_creatematrix ( integer  N,
integer  Subband,
logical  Symmetric,
logical  AllocValues 
)

Create the structures required for a Band format matrix.

Parameters
nNumber of rows for the matrix
subbandMax(ABS(Col-Diag(Row))) of the matrix
symmetricSymmetric or not
allocvaluesShould the values arrays be allocated ?

References generalutils::allocatematrix(), and messages::fatal().

Referenced by elementutils::creatematrix().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine bandmatrix::band_diagprecondition ( real(kind=dp), dimension(*)  u,
real(kind=dp), dimension(*)  v,
integer, dimension(*)  ipar 
)

Diagonal preconditioning of a Band format matrix.

Parameters
uVector to be multiplied
vResult vector of the multiplication
iparstructure holding info from (HUTIter-iterative solver package)
real(kind=dp) function bandmatrix::band_getmatrixelement ( type(matrix_t)  A,
integer  i,
integer  j 
)

Get a given matrix element of a band matrix format.

Parameters
aStructure holding matrix
irow number of the matrix element
jcolumn number of the matrix element

Referenced by solverutils::getmatrixelement().

Here is the caller graph for this function:

subroutine bandmatrix::band_gluelocalmatrix ( type(matrix_t)  A,
integer  N,
integer  Dofs,
integer, dimension(:)  Indeces,
real(kind=dp), dimension(:,:)  LocalMatrix 
)

Add a set of values (.i.e. element stiffness matrix) to a Band format matrix.

Parameters
localmatrixA (N x Dofs) x ( N x Dofs) matrix holding the values to be added to the Band format matrix
aStructure holding matrix, values are affected in the process
nNumber of nodes in element
dofsNumber of degrees of freedom for one node
indecesMaps element node numbers to global (or partition) node numbers (to matrix rows and columns, if Dofs = 1)

References localmatrix().

Referenced by solverutils::updateglobalequations(), updateglobalpreconditioner(), and solverutils::updatemassmatrix().

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine bandmatrix::band_matrixvectormultiply ( type(matrix_t)  A,
real(kind=dp), dimension(*)  u,
real(kind=dp), dimension(*)  v 
)

Matrix vector product (v = Au) for a matrix given in Band format.

Parameters
uvector to be multiplied
vresult vector
aStructure holding the matrix

Referenced by solverutils::matrixvectormultiply().

Here is the caller graph for this function:

subroutine bandmatrix::band_matrixvectorprod ( real(kind=dp), dimension(*)  u,
real(kind=dp), dimension(*)  v,
integer, dimension(*)  ipar 
)

Matrix vector product (v = Au) for a matrix given in Band format. The matrix is accessed from a global variable GlobalMatrix.

Parameters
uVector to be multiplied
vResult vector of the multiplication
iparstructure holding info from (HUTIter-iterative solver package)
subroutine bandmatrix::band_setmatrixelement ( type(matrix_t)  A,
integer  i,
integer  j,
real(kind=dp)  value 
)

Set a given element in the band matrix.

Parameters
aStructure holding matrix
irow number of the matrix element
jcolumn number of the matrix element
valueValue to be set

Referenced by solverutils::setmatrixelement().

Here is the caller graph for this function:

subroutine bandmatrix::band_zeromatrix ( type(matrix_t)  A)

Zero a Band format matrix.

Parameters
aStructure holding matrix

Referenced by solverutils::initializetozero().

Here is the caller graph for this function:

subroutine bandmatrix::band_zerorow ( type(matrix_t)  A,
integer  n 
)

Zero given row from a Band format matrix.

Parameters
aStructure holding matrix
nRow number to be zerod

Referenced by solverutils::zerorow().

Here is the caller graph for this function:

subroutine bandmatrix::sband_setdirichlet ( type(matrix_t)  A,
real(kind=dp), dimension(:)  b,
integer, intent(in)  n,
real(kind=dp), intent(in)  Value 
)

Set value of unkown x_n to given value for symmetric band matrix. This is done by replacing the equation of the unknown by x_n = Value (i.e. zeroing the row of the unkown in the matrix, and setting diagonal to identity). Also the respective column is set to zero (except for the diagonal) to preserve symmetry, while also substituting the rhs by by rhs(i) = rhs(i) - A(i,n) * Value.

Parameters
aStructure holding matrix, values are affected in the process
bRHS vector
[in]valueValue for the unknown
[in]nOrdered number of the unkown (i.e. matrix row and column number)

Referenced by outletcompute(), solverutils::setdirichletpoint(), setdirichletpoints(), setelementvalues(), setlimitervalues(), and setpointvalues().

Here is the caller graph for this function:


The documentation for this module was generated from the following file: