#include <densemat.h>
Public Member Functions | |
void | add_entry (double e, long ri, long ci) |
void | addmat_dm (double c, densemat &dm) |
void | alloc (long m) |
void | assemble_dense_from_scr (long *scradr, long *scrci, double *scra, long neq) |
void | cg (double *x, double *y, long ni, double err, long &ani, double &ares, double zero, long iv) |
void | cg_prec (densemat &dm, double *x, double *y, long ni, double err, long &ani, double &ares, double zero, long iv, long tprec, double par) |
void | cg_prec_new (precond &pr, double *x, double *y, long ni, double err, long &ani, double &ares, double zero, long iv) |
void | changedecomp () |
void | copy (densemat *dm) |
void | copy_dm (densemat &dm) |
void | dealloc (void) |
long | decomp () |
densemat (void) | |
void | diag_scale (double *d) |
double | estim_spect_radius () |
void | gemp (double *x, double *y, long m, double limit, long pivot) |
void | gempkon (double *b, double *c, double *x, double *y, long m, double zero, long tc) |
double | give_entry (long ri, long ci) |
long | give_negm () |
void | glocalize (matrix &b, long *rcn, long *ccn) |
void | ill (double *x, double *y, double zero, double limit, long tc) |
void | initiate (long ndof, long mespr) |
void | ker (double *r, long &dim, long *se, long ense, double limit) |
void | ll (double *x, double *y, double zero, long tc) |
void | localize (matrix &b, long *cn) |
void | localized (double *b, long *cn, long m) |
void | lu (double *x, double *y, double zero, long tc) |
void | maxmb (long nm, double *ma, double *mb, double *mc) |
void | mult_localize (long nm, long *ncn1, long *ncn2, long *mcn) |
void | mxv_dm (double *b, double *c) |
void | nullmat () |
void | printdiag (FILE *out) |
void | printmat (FILE *out) |
void | scalmat_dm (double c) |
void | setfact () |
void | setnotfact () |
double * | status () |
~densemat (void) | |
Public Attributes | |
double * | a |
array containing global matrix | |
double * | a_plus |
array containing pseudoinverse (inverse) of global matrix | |
long * | abvns |
address of first entry of i-th base vector in bvns | |
double * | bvns |
base vectors of null space of global matrix | |
long | decompid |
decomposition (factorization) indicator | |
long * | iirm |
array with number of row indices of dependent rows in global matix | |
long | n |
number of rows (columns) | |
long | nbvns |
number of base vectors of null space of global matrix | |
long | negm |
number of entries in the global matrix |
class densemat
this class serves for storage of dense matrices
JK
Definition at line 18 of file densemat.h.
~densemat | ( | void | ) |
void add_entry | ( | double | e, | |
long | ri, | |||
long | ci | |||
) |
function adds required matrix entry
e | - matrix entry | |
ri,ci | - row and column indices |
JK, 24.7.2005
Definition at line 1381 of file densemat.cpp.
Referenced by gmatrix::add_entry(), and aggregator::coarse_matrix().
void addmat_dm | ( | double | c, | |
densemat & | dm | |||
) |
function adds premultiplied matrix stored in dm by coefficient c to actual matrix
c | - multiplicative coefficient | |
dm | - another dense matrix storage |
JK
Definition at line 300 of file densemat.cpp.
void alloc | ( | long | m | ) |
function allocates array containing matrix
m | - number of rows/columns |
JK
Definition at line 30 of file densemat.cpp.
References a, memset(), n, and negm.
Referenced by aggregator::coarse_matrix(), initiate(), transtop::mod_view_factors(), saddpoint::solve_system(), and gmatrix::solve_system().
void assemble_dense_from_scr | ( | long * | scradr, | |
long * | scrci, | |||
double * | scra, | |||
long | neq | |||
) |
function assembles matrix in dense format form scr format
scradr - array with adr in scr format scrci - array with ci in scr format scra - array with a in scr format
12.02.2009 JB
Definition at line 1456 of file densemat.cpp.
void cg | ( | double * | x, | |
double * | y, | |||
long | ni, | |||
double | err, | |||
long & | ani, | |||
double & | ares, | |||
double | zero, | |||
long | iv | |||
) |
function solves system of linear algebraic equations by conjugate gradient method, matrix is stored as dense notice: this function can be used only for symmetric positive definite matrices function was implemented for testing other conjugate gradient methods performed on more complicated schemes of matrix storage
x | - vector of unknowns | |
y | - vector of right hand side | |
ni | - maximum number of iterations | |
err | - required error (norm of residual vector) | |
ani | - number of performed iterations | |
ares | - attained error (norm of residual vector) | |
zero | - computer zero | |
iv | - initial values indicator |
iv=0 - initial vector is zero vector iv=1 - initial vector is taken from x array
JK, 17.12.2005
Definition at line 969 of file densemat.cpp.
References mxv_dm(), n, p, print_err(), and ss().
Referenced by gmatrix::solve_system().
void cg_prec | ( | densemat & | dm, | |
double * | x, | |||
double * | y, | |||
long | ni, | |||
double | err, | |||
long & | ani, | |||
double & | ares, | |||
double | zero, | |||
long | iv, | |||
long | tprec, | |||
double | par | |||
) |
x | - vector of unknowns | |
y | - vector of right hand side | |
ni | - maximum number of iterations | |
err | - required error (norm of residula vector) | |
ani | - number of computed iterations | |
ares | - reached residual norm | |
zero | - computer zero | |
iv | - initial values indicator | |
tprec | - type of preconditioner | |
par | - preconditioning parameter |
iv=0 - initial vector is assumed as zero vector iv=1 - initial vector is taken from x array
tprec=1 - diagonal (Jacobi) preconditioner tprec=10 - incomplete decomposition
par: relaxation parameter - for SSOR preconditioner weight coefficient - for incomplete decomposition
JK, 17.12.2005
Definition at line 1072 of file densemat.cpp.
References h(), ill(), mxv_dm(), n, p, print_err(), and ss().
void cg_prec_new | ( | precond & | pr, | |
double * | x, | |||
double * | y, | |||
long | ni, | |||
double | err, | |||
long & | ani, | |||
double & | ares, | |||
double | zero, | |||
long | iv | |||
) |
Definition at line 1195 of file densemat.cpp.
References h(), mxv_dm(), n, p, precond::preconditioning(), print_err(), and ss().
Referenced by gmatrix::solve_system().
void changedecomp | ( | ) |
function changes indicator of decomposition (factorization)
JK
Definition at line 115 of file densemat.cpp.
References decompid.
Referenced by gmatrix::changedecomp().
void copy | ( | densemat * | dm | ) |
Definition at line 44 of file densemat.cpp.
Referenced by gmatrix::solve_system().
void copy_dm | ( | densemat & | dm | ) |
function creates a copy of the matrix
dm | - dense matrix which will be set up by data |
JK, 15.3.2007
Definition at line 69 of file densemat.cpp.
void dealloc | ( | void | ) |
Definition at line 38 of file densemat.cpp.
References a.
Referenced by saddpoint::solve_system(), and gmatrix::solve_system().
long decomp | ( | ) |
function returns indicator of decomposition (factorization)
JK
Definition at line 105 of file densemat.cpp.
References decompid.
Referenced by gmatrix::back_substitution(), gmatrix::decomp(), gmatrix::decompose_matrix(), and gmatrix::solve_system().
void diag_scale | ( | double * | d | ) |
function scales matrix by its diagonal elements
d | - array, where scale factors from diagonal are stored |
JK, 28.5.2008
Definition at line 1404 of file densemat.cpp.
Referenced by gmatrix::diag_scale().
double estim_spect_radius | ( | ) |
function estimates spectral radius
the estimates are based on the Gershgorin circle theorem for details see G.H. Golub, C.F. Van Loan: Matrix computations. The Johns Hopkins University Press, 3rd edition, page 320, 1996
JK, 27.8.2008
Definition at line 1428 of file densemat.cpp.
Referenced by gmatrix::estim_spect_radius().
void gemp | ( | double * | x, | |
double * | y, | |||
long | m, | |||
double | zero, | |||
long | pivot | |||
) |
function solves equation system by Gauss elimination method matrix is stored as dense matrix x and y are also dense matrices, stored by rows right hand side vectors and vectors of unknowns are columns!
x | - array containing vectors of unknowns | |
y | - array containing right hand sides | |
m | - number of right hand sides | |
zero | - computer zero | |
pivot | - type of pivoting |
pivot=1 - pivot is searched only when there is zero on the diagonal pivot=2 - pivot is searched every time
JK, 23.7.2001
Definition at line 372 of file densemat.cpp.
References a, f, n, and print_err().
Referenced by saddpoint::solve_system(), and gmatrix::solve_system().
void gempkon | ( | double * | b, | |
double * | c, | |||
double * | x, | |||
double * | y, | |||
long | m, | |||
double | zero, | |||
long | tc | |||
) |
function reduces equation system by Gauss elimination method matrix is stored as dense matrix x and y are also dense matrices, stored by rows right hand side vectors and vectors of unknowns are columns!
it is condensation strategy
b | - array containing condensed matrix | |
x | - array containing left hand side | |
y | - array containing right hand side | |
m | - number of remaining unknowns | |
zero | - computer zero | |
tc | - type of computation |
tc=1 - function decomposes matrix and creates condensed matrix tc=2 - function computes backward substitution
JK, 2.6.2002
Definition at line 553 of file densemat.cpp.
Referenced by gmatrix::condense().
double give_entry | ( | long | ri, | |
long | ci | |||
) |
function returns required matrix entry
ri,ci | - row and column indices |
JK, 24.7.2005
Definition at line 1366 of file densemat.cpp.
Referenced by gmatrix::give_entry().
long give_negm | ( | ) |
function returns number of stored matrix entries
JK, 16.8.2007
Definition at line 1391 of file densemat.cpp.
References negm.
Referenced by gmatrix::give_negm().
void glocalize | ( | matrix & | b, | |
long * | rcn, | |||
long * | ccn | |||
) |
function localizes general rectangular (non-square) matrix into global matrix
b | - local matrix | |
rcn | - row code numbers | |
ccn | - column code numbers |
JK
Definition at line 201 of file densemat.cpp.
References a, matrix::m, n, and matrix::n.
Referenced by gmatrix::glocalize().
void ill | ( | double * | x, | |
double * | y, | |||
double | zero, | |||
double | limit, | |||
long | tc | |||
) |
function computes incomplete Cholesky decomposition of matrix and solves system of equations
notice: this function can be used only for symmetric positive definite matrices function was implemented for testing other Cholesky decompositions performed on more complicated schemes of matrix storage
x | - left hand side | |
y | - right hand side | |
zero | - computer zero | |
limit | - threshold for fill-in | |
tc | - type of computation |
tc = 1 - decomposition and back substitution tc = 2 - decomposition of matrix tc = 3 - back substitution
JK, 17.12.2005
Definition at line 796 of file densemat.cpp.
References a, decompid, n, and print_err().
Referenced by gmatrix::back_incomplete_fact(), cg_prec(), gmatrix::incomplete_fact(), and gmatrix::solve_system().
void initiate | ( | long | ndof, | |
long | mespr | |||
) |
function initiates densemat storage
ndof | - number of rows/columns of the matrix | |
mespr | - message printing indicator |
JK
Definition at line 255 of file densemat.cpp.
References alloc(), negm, nullmat(), and status().
Referenced by gmatrix::initiate().
void ker | ( | double * | r, | |
long & | dim, | |||
long * | se, | |||
long | ense, | |||
double | limit | |||
) |
function computes kernel of matrix
r | - array containing base vectors of kernel | |
dim | - dimension of kernel | |
se | - array containing numbers (indices) of singular equations | |
ense | - estimated number of singular equations | |
limit | - defines zero on diagonal |
21.4.2003, JK
Definition at line 872 of file densemat.cpp.
References a, n, nullv(), and print_err().
Referenced by seqfeti::kernel(), and gmatrix::kernel().
void ll | ( | double * | x, | |
double * | y, | |||
double | zero, | |||
long | tc | |||
) |
function computes Cholesky decomposition of matrix and solves system of equations
notice: this function can be used only for symmetric positive definite matrices function was implemented for testing other Cholesky decompositions performed on more complicated schemes of matrix storage
x | - left hand side | |
y | - right hand side | |
zero | - computer zero | |
tc | - type of computation |
tc = 1 - decomposition and back substitution tc = 2 - decomposition of matrix tc = 3 - back substitution
JK, 17.12.2005
Definition at line 714 of file densemat.cpp.
References a, decompid, n, and print_err().
Referenced by aggregator::prepare_boss(), and gmatrix::solve_system().
void localize | ( | matrix & | b, | |
long * | cn | |||
) |
function localizes local matrix b into global matrix
b | - local matrix | |
cn | - array containing code numbers |
JK
Definition at line 153 of file densemat.cpp.
References a, matrix::m, and n.
Referenced by gmatrix::localize().
void localized | ( | double * | b, | |
long * | cn, | |||
long | m | |||
) |
function localizes local matrix stored in array into global matrix
b | - array containing local matrix | |
cn | - array containing code numbers | |
m | - number of rows/columns of the matrix b |
JK
Definition at line 177 of file densemat.cpp.
Referenced by gmatrix::localized().
void lu | ( | double * | x, | |
double * | y, | |||
double | zero, | |||
long | tc | |||
) |
function decomposes matrix A into L.U form
x | - left hand side | |
y | - right hand side | |
zero | - computer zero | |
tc | - type of computation |
tc = 1 - decomposition and back substitution tc = 2 - decomposition of matrix tc = 3 - back substitution
JK, 6.1.2000
Definition at line 635 of file densemat.cpp.
References a, decompid, n, and print_err().
Referenced by gmatrix::back_substitution(), gmatrix::decompose_matrix(), transtop::heat_fluxes(), gmatrix::inverse_iteration(), transtop::mod_view_factors(), and gmatrix::solve_system().
void maxmb | ( | long | nm, | |
double * | ma, | |||
double * | mb, | |||
double * | mc | |||
) |
function multiplies matrix a stored in dense format with matrix b stored in dense format ma - matrix a mb - matrix b nm - number of colums of matrix a and number of rows of matrix b mc - result of computation
JB
Definition at line 336 of file densemat.cpp.
References n.
void mult_localize | ( | long | nm, | |
long * | ncn1, | |||
long * | ncn2, | |||
long * | mcn | |||
) |
function localizes contributions from Lagrange multipliers to the skyline storage
nm | - number of Lagrange multipliers | |
ncn1 | - nodal code numbers of the first node (on one side of interface) | |
ncn2 | - nodal code numbers of the second node (on the other side of interface) | |
mcn | - code numbers of Lagrange multipliers defined between the previous nodes |
JK, 18.8.2008
Definition at line 226 of file densemat.cpp.
Referenced by gmatrix::mult_localize().
void mxv_dm | ( | double * | b, | |
double * | c | |||
) |
function multiplies matrix stored as dense matrix by vector b, resulting vector is c
b | - array containing vector b | |
c | - array containing vector c |
JK
Definition at line 276 of file densemat.cpp.
Referenced by cg(), cg_prec(), cg_prec_new(), feti1::dirichletprec(), seqfeti::dirichletprec(), and transtop::heat_fluxes().
void nullmat | ( | ) |
function fills array by zeros
JK
Definition at line 135 of file densemat.cpp.
Referenced by initiate().
void printdiag | ( | FILE * | out | ) |
function prints diagonal entries of the matrix
out | - output stream |
JK
Definition at line 1349 of file densemat.cpp.
Referenced by gmatrix::printdiag().
void printmat | ( | FILE * | out | ) |
function prints matrix into output file
out | - output stream |
JK
Definition at line 1306 of file densemat.cpp.
Referenced by gmatrix::printmat().
void scalmat_dm | ( | double | c | ) |
function multiplies components of the matrix by coefficient c
c | - multiplicative coefficient |
JK
Definition at line 317 of file densemat.cpp.
void setfact | ( | ) |
Definition at line 121 of file densemat.cpp.
References decompid.
void setnotfact | ( | ) |
Definition at line 125 of file densemat.cpp.
References decompid.
double * status | ( | ) |
function returns status of array a
JK
Definition at line 95 of file densemat.cpp.
References a.
Referenced by initiate().
double* a |
array containing global matrix
Definition at line 75 of file densemat.h.
Referenced by add_entry(), addmat_dm(), alloc(), assemble_dense_from_scr(), copy(), copy_dm(), dealloc(), densemat(), diag_scale(), estim_spect_radius(), gemp(), gempkon(), give_entry(), glocalize(), ill(), ker(), ll(), localize(), localized(), lu(), transtop::mod_view_factors(), mult_localize(), mxv_dm(), nullmat(), printdiag(), printmat(), scalmat_dm(), comprow::select_submatrix(), saddpoint::solve_system(), gmatrix::solve_system(), status(), seqfeti::subdomain_matrices(), feti1::subdomain_matrix(), and ~densemat().
double* a_plus |
array containing pseudoinverse (inverse) of global matrix
Definition at line 80 of file densemat.h.
long* abvns |
address of first entry of i-th base vector in bvns
Definition at line 86 of file densemat.h.
double* bvns |
base vectors of null space of global matrix
Definition at line 84 of file densemat.h.
long decompid |
decomposition (factorization) indicator
Definition at line 77 of file densemat.h.
Referenced by changedecomp(), decomp(), densemat(), ill(), ll(), lu(), gmatrix::prepmat(), setfact(), and setnotfact().
long* iirm |
array with number of row indices of dependent rows in global matix
Definition at line 88 of file densemat.h.
long n |
number of rows (columns)
Definition at line 71 of file densemat.h.
Referenced by add_entry(), addmat_dm(), alloc(), assemble_dense_from_scr(), cg(), cg_prec(), cg_prec_new(), copy(), copy_dm(), densemat(), diag_scale(), estim_spect_radius(), gemp(), gempkon(), give_entry(), glocalize(), ill(), ker(), ll(), localize(), localized(), lu(), maxmb(), mult_localize(), mxv_dm(), nullmat(), printdiag(), printmat(), scalmat_dm(), comprow::select_submatrix(), seqfeti::subdomain_matrices(), and feti1::subdomain_matrix().
long nbvns |
number of base vectors of null space of global matrix
Definition at line 82 of file densemat.h.
long negm |
number of entries in the global matrix
Definition at line 73 of file densemat.h.
Referenced by alloc(), copy(), copy_dm(), give_negm(), initiate(), seqfeti::subdomain_matrices(), and feti1::subdomain_matrix().