62 real(RP),
dimension(:),
allocatable :: invd
68 integer ,
dimension(:),
allocatable :: perm
71 integer ,
dimension(:),
allocatable :: perminv
122 type(
csr),
intent(inout) :: K
123 integer ,
intent(in) :: type
126 if(.NOT.
valid(k))
call quit(
"precond_mod: precond_create:& 127 & matrix 'K' not valid")
145 call quit(
'precond_mod: precond_create:& 146 & preconditioning type incorrect')
151 &
'precond_mod : create = ',&
163 call reorder(k, pc%perm, pc%permInv)
172 subroutine jacobi(invD, s)
173 real(RP) ,
dimension(:),
allocatable :: invD
174 type(
csr) ,
intent(inout) :: s
179 &
'precond_mod : jacobi' 184 if (minval(s%diag) <= 0 )
call quit( &
185 &
'precond_mod: jacobi: negative diagonal entries')
191 invd(ii) = 1._rp / s%a(jj)
201 type(
csr) ,
intent(in) :: s
202 integer ,
dimension(:),
allocatable :: perm, permInv
204 integer,
dimension(s%nl) :: deg, P, LS, LS0, degLS, new_i
205 logical,
dimension(s%nl) :: marked
206 integer :: ii, jj, jDeb, jFin, kk
207 integer :: nL, nL0, node
210 &
'precond_mod : reverse_cuthill_macKey' 220 deg(ii) = s%row(ii+1)-s%row(ii)
228 if (deg(ii)==jj)
then 233 marked(node) = .true.
247 jdeb = s%row( ls0(ii) )
248 jfin = s%row( ls0(ii)+1 ) - 1
254 if (.not.marked(node))
then 258 marked(node) = .true.
270 degls(ii) = deg(ls(ii))
277 p(kk) = ls(new_i(ii))
283 perm(ii) = p(s%nl+1-ii)
287 perminv(perm(ii)) = ii
296 subroutine icc0(LLT, s)
297 type(
csr),
intent(out) :: LLT
298 type(
csr),
intent(inout) :: s
300 integer :: i,k,j,m,mm,band, kDeb,kFin
301 real(RP) :: int1,int2,sum, t1, t2
304 &
' precond_mod : icc0 =' 314 band = max( band, s%col(s%row(i+1)-1)-i )
317 & int(band,4),
' : band width' 326 do m = s%diag(i),s%row(i+1)-1
329 do k = max(i-band-1,1),i-1
332 do mm = s%diag(k),s%row(k+1)-1
333 if (s%col(mm)==i) int1 = llt%a(mm)
334 if (s%col(mm)==j) int2 = llt%a(mm)
336 sum = sum - int1*int2
341 llt%a(m) = sum / llt%a(s%diag(i))
350 kfin = llt%row(i+1)-1
354 do m = llt%row(j), llt%diag(j)-1
355 if (llt%col(m)==i)
then subroutine precond_clear(pc)
subroutine, public icc0(LLT, s)
Incomplete Choleski decomposition with 0 fill-in for s.
subroutine, public reverse_cuthill_mackey(perm, permInv, s)
Reverse CutHill Mac Kee algorithm.
deallocate memory for real(RP) arrays
The type precond defines preconditioning for linear systems.
type(precond) function precond_create(K, type)
Constructor for the type precond
subroutine, public jacobi(invD, s)
Jacobi preconditioner.
REAL NUMBERS PRECISION IN CHORAL: selects simple/double/quad
Derived type for sparse matrices in CSR format.
ALGEBRAIC OPERATIONS ON SETS
integer, parameter pc_jacobi
Jacobi diagonal preconditioning.
allocate memory for real(RP) arrays
subroutine, public getdiag(A)
Computes the array Adiag of the indices in Acol for the diagonal coefficients.
reorder rows and cols of a matrix following a permutation of the indexes
integer, parameter pc_icc0
Incomplete Cholesky, order 0.
integer choral_verb
Verbosity level.
DEFINITION OF GLOBAL VARIABLES FOR THE LIBRARY CHORAL
DERIVED TYPE precond: preconditioning for linear systems
subroutine, public shellsort_dec(t, new_i, n)
Sort integer array : shell sort.
DERIVED TYPE csr for sparse matrices
integer, parameter pc_0
void preconditioning