Choral
Data Types | Functions/Subroutines
csr_mod Module Reference

DERIVED TYPE csr for sparse matrices More...

Data Types

interface  add
 addition of csr matrices (with the same pattern) More...
 
interface  add_submatrix
 add a CSR matrix M1 with pattern P1 to a CSR matrix M with pattern P such that P1 is included in P More...
 
interface  addentry
 
interface  clear
 destructor More...
 
interface  copy
 M2 = M1. More...
 
interface  csr
 Derived type for sparse matrices in CSR format. More...
 
interface  equal
 test csr matrices equality More...
 
interface  getrow
 
interface  invllt
 x solution of LLT x = b, L lower diagonal LT = transpose(L) More...
 
interface  matvecprod
 matrix vector product More...
 
interface  print
 print a short description More...
 
interface  reorder
 reorder rows and cols of a matrix following a permutation of the indexes More...
 
interface  scale
 scale the matrix entries More...
 
interface  setrow
 
interface  sortrows
 
interface  transpose
 transpose matrix of a sparse csr matrix More...
 
interface  valid
 is the csr matrix well defined ? More...
 
interface  write
 write the matrix to a file for a given format More...
 

Functions/Subroutines

subroutine csr_clear (m)
 Destructor for csr type. More...
 
type(csr) function csr_create (nnz)
  Constructor for the type csr More...
 
logical function csr_valid (m)
 Check if structure contents is correct. More...
 
subroutine csr_getrow (sz, v, p, n, g, ii)
 Extract row ii of the csr matrix g. More...
 
subroutine csr_setrow (m, val, col, ii)
 Set line ii. More...
 
logical function csr_equal (g1, g2)
 test equality More...
 
subroutine csr_addentry (m, rw, cl, val)
 Filling matrix m entries. More...
 
subroutine, public addentry_2 (m, rw, cl, val)
 Filling matrix m entries. More...
 
subroutine csr_write (m, fileName, fileFormat)
 Write matrix m to file 'fileName' with format 'fileFormat'. More...
 
subroutine matlab_write (m, fileName)
 Write matrix m to file 'fileName' matlab ASCII format. More...
 
subroutine csr_copy (s2, s1)
 s2 = s1, copy : s2 is nullified, reallocated and set to s1 More...
 
subroutine csr_matvecprod (v, m, u)
 Matrix Vector Product v <= m*u. More...
 
subroutine csr_add_2 (m, m1, a, m2, b)
 m <= a*m1 + b*m2 m1, m2 with identical patterns More...
 
subroutine csr_add_1 (m1, a, m2, b)
 m1 = a*m1 + b*m2 m1, m2 with identical patterns More...
 
subroutine csr_scale (m1, a)
 m1 <= a*m1 More...
 
subroutine csr_reorder (s, perm, permInv)
 reorder s according to the permutation perm, permInv = perm^{-1} More...
 
subroutine csr_sortrows (s)
 Sort column indexing per line (increasingly) More...
 
subroutine csr_invllt (x, LLT, b)
 x solution of LLT x = b, L lower diagonal LT = transpose(L) More...
 
subroutine invlt (x, LT, b)
 x solution of LT x = b, L lower diagonal LT = transpose(L) More...
 
subroutine invl (x, L, b)
 x solution of L x = b More...
 
subroutine, public getdiag (A)
 Computes the array Adiag of the indices in Acol for the diagonal coefficients. More...
 
subroutine csr_transpose (g2, g1)
 Transpose csr. More...
 
subroutine add_submatrix_1 (m, m2)
 m <= m + m2 m2(i,j) /=0 ==> m(i,j)/=0 More...
 
subroutine add_submatrix_2 (m, m2, a)
 m <= m + a*m2 m2(i,j) /=0 ==> m(i,j)/=0 More...
 
subroutine csr_print (g)
 Print a short description. More...
 

Detailed Description

DERIVED TYPE csr for sparse matrices

Storage of a sparse matrix with real entries with CSR format (Compressed Sparse Rows).

Author
Charles Pierre

Function/Subroutine Documentation

◆ add_submatrix_1()

subroutine csr_mod::add_submatrix_1 ( type(csr), intent(inout)  m,
type(csr), intent(in)  m2 
)
private

m <= m + m2 m2(i,j) /=0 ==> m(i,j)/=0

Definition at line 913 of file csr_mod.F90.

Here is the call graph for this function:

◆ add_submatrix_2()

subroutine csr_mod::add_submatrix_2 ( type(csr), intent(inout)  m,
type(csr), intent(in)  m2,
real(rp), intent(in)  a 
)
private

m <= m + a*m2 m2(i,j) /=0 ==> m(i,j)/=0

Definition at line 943 of file csr_mod.F90.

Here is the call graph for this function:

◆ addentry_2()

subroutine, public csr_mod::addentry_2 ( type (csr), intent(inout)  m,
integer, intent(in)  rw,
integer, intent(in)  cl,
real(rp), intent(in)  val 
)

Filling matrix m entries.

m(rw,cl) is already set !

Definition at line 418 of file csr_mod.F90.

◆ csr_add_1()

subroutine csr_mod::csr_add_1 ( type(csr), intent(inout)  m1,
real(rp), intent(in)  a,
type(csr), intent(in)  m2,
real(rp), intent(in)  b 
)
private

m1 = a*m1 + b*m2 m1, m2 with identical patterns

Definition at line 598 of file csr_mod.F90.

Here is the call graph for this function:

◆ csr_add_2()

subroutine csr_mod::csr_add_2 ( type(csr), intent(inout)  m,
type(csr), intent(in)  m1,
real(rp), intent(in)  a,
type(csr), intent(in)  m2,
real(rp), intent(in)  b 
)
private

m <= a*m1 + b*m2 m1, m2 with identical patterns

Definition at line 563 of file csr_mod.F90.

Here is the call graph for this function:

◆ csr_addentry()

subroutine csr_mod::csr_addentry ( type (csr), intent(inout)  m,
integer, intent(in)  rw,
integer, intent(in)  cl,
real(rp), intent(in)  val 
)
private

Filling matrix m entries.

If m(rw,cl) already set : m(rw,cl) = m(rw,cl) + val

If m(rw,cl) not set : m(rw,cl) = val

If the line rw is full : error

Definition at line 385 of file csr_mod.F90.

Here is the call graph for this function:

◆ csr_clear()

subroutine csr_mod::csr_clear ( type(csr), intent(inout)  m)
private

Destructor for csr type.

Definition at line 180 of file csr_mod.F90.

◆ csr_copy()

subroutine csr_mod::csr_copy ( type(csr), intent(inout)  s2,
type(csr), intent(in)  s1 
)
private

s2 = s1, copy : s2 is nullified, reallocated and set to s1

Definition at line 487 of file csr_mod.F90.

Here is the call graph for this function:

◆ csr_create()

type(csr) function csr_mod::csr_create ( integer, dimension(:), intent(in)  nnz)
private

Constructor for the type csr

  • OUTPUT:
    • m = csr matrix
  • INPUT:
    • nnz: integer array, nnz(i) = number of non zero entries on line i

Definition at line 204 of file csr_mod.F90.

Here is the call graph for this function:

◆ csr_equal()

logical function csr_mod::csr_equal ( type(csr), intent(in)  g1,
type(csr), intent(in)  g2 
)
private

test equality

g1 and g2 will be modified : rows will be sorted

Definition at line 338 of file csr_mod.F90.

Here is the call graph for this function:

◆ csr_getrow()

subroutine csr_mod::csr_getrow ( integer, intent(out)  sz,
real(rp), dimension(n), intent(out)  v,
integer, dimension(n), intent(out)  p,
integer, intent(in)  n,
type(csr), intent(in)  g,
integer, intent(in)  ii 
)
private

Extract row ii of the csr matrix g.

OUTPUT : sz = row ii size v(1:sz) = entries of row ii p(1:sz) = associated column indexes

INPUT : n = dimension of p and of v WARINIG : n must be larger than sz g = CSR matrix ii= row index

Definition at line 286 of file csr_mod.F90.

Here is the call graph for this function:

◆ csr_invllt()

subroutine csr_mod::csr_invllt ( real(rp), dimension(:), intent(out)  x,
type(csr), intent(in)  LLT,
real(rp), dimension(:), intent(in)  b 
)
private

x solution of LLT x = b, L lower diagonal LT = transpose(L)

Definition at line 734 of file csr_mod.F90.

Here is the call graph for this function:

◆ csr_matvecprod()

subroutine csr_mod::csr_matvecprod ( real(rp), dimension(:), intent(out)  v,
type(csr), intent(in)  m,
real(rp), dimension(:), intent(in)  u 
)
private

Matrix Vector Product v <= m*u.

Definition at line 522 of file csr_mod.F90.

Here is the call graph for this function:

◆ csr_print()

subroutine csr_mod::csr_print ( type(csr), intent(in)  g)
private

Print a short description.

Definition at line 974 of file csr_mod.F90.

◆ csr_reorder()

subroutine csr_mod::csr_reorder ( type(csr), intent(inout)  s,
integer, dimension(:), intent(in)  perm,
integer, dimension(:), intent(in)  permInv 
)
private

reorder s according to the permutation perm, permInv = perm^{-1}

Definition at line 652 of file csr_mod.F90.

Here is the call graph for this function:

◆ csr_scale()

subroutine csr_mod::csr_scale ( type(csr), intent(inout)  m1,
real(rp), intent(in)  a 
)
private

m1 <= a*m1

Definition at line 629 of file csr_mod.F90.

Here is the call graph for this function:

◆ csr_setrow()

subroutine csr_mod::csr_setrow ( type(csr), intent(inout)  m,
real(rp), dimension(:), intent(in)  val,
integer, dimension(:), intent(in)  col,
integer, intent(in)  ii 
)
private

Set line ii.

Definition at line 319 of file csr_mod.F90.

◆ csr_sortrows()

subroutine csr_mod::csr_sortrows ( type(csr), intent(inout)  s)
private

Sort column indexing per line (increasingly)

Definition at line 699 of file csr_mod.F90.

Here is the call graph for this function:

◆ csr_transpose()

subroutine csr_mod::csr_transpose ( type(csr), intent(inout)  g2,
type(csr), intent(in)  g1 
)
private

Transpose csr.

Definition at line 870 of file csr_mod.F90.

Here is the call graph for this function:

◆ csr_valid()

logical function csr_mod::csr_valid ( type (csr), intent(in)  m)
private

Check if structure contents is correct.

Definition at line 242 of file csr_mod.F90.

◆ csr_write()

subroutine csr_mod::csr_write ( type(csr), intent(in)  m,
character(len=*), intent(in)  fileName,
character(len=*), intent(in)  fileFormat 
)
private

Write matrix m to file 'fileName' with format 'fileFormat'.

Definition at line 440 of file csr_mod.F90.

Here is the call graph for this function:

◆ getdiag()

subroutine, public csr_mod::getdiag ( type(csr), intent(inout)  A)

Computes the array Adiag of the indices in Acol for the diagonal coefficients.

If line ii has no diagonal entry, then Adiag(ii) = -1

Definition at line 835 of file csr_mod.F90.

Here is the call graph for this function:

◆ invl()

subroutine csr_mod::invl ( real(rp), dimension(:), intent(out)  x,
type(csr), intent(in)  L,
real(rp), dimension(:), intent(in)  b 
)
private

x solution of L x = b

Definition at line 796 of file csr_mod.F90.

Here is the call graph for this function:

◆ invlt()

subroutine csr_mod::invlt ( real(rp), dimension(:), intent(out)  x,
type(csr), intent(in)  LT,
real(rp), dimension(:), intent(in)  b 
)
private

x solution of LT x = b, L lower diagonal LT = transpose(L)

Definition at line 761 of file csr_mod.F90.

Here is the call graph for this function:

◆ matlab_write()

subroutine csr_mod::matlab_write ( type(csr), intent(in)  m,
character(len=*), intent(in)  fileName 
)
private

Write matrix m to file 'fileName' matlab ASCII format.

Definition at line 457 of file csr_mod.F90.

Here is the call graph for this function: