36 #include "qrm_common.h" 68 subroutine dqrm_rowperm(graph, cperm, rperm, nvar, stair)
76 integer :: cperm(:), rperm(:), stair(:), nvar(:)
79 integer,
allocatable :: mark(:)
80 integer :: ii, jj, i, j, node, nv, pnt, rpnt
83 character(len=*),
parameter :: name=
'qrm_rowperm' 85 call qrm_err_act_save(err_act)
88 __qrm_check_ret(name,
'qrm_aalloc',9999)
96 if (pnt .gt. graph%n)
exit 100 do jj = pnt, pnt+nv-1
102 do ii=graph%jptr(j), graph%jptr(j+1)-1
104 if(mark(i) .eq. 0)
then 105 stair(node) = stair(node)+1
113 if (pnt .gt. graph%n)
exit 118 if(mark(i) .eq. 0)
then 126 __qrm_check_ret(name,
'qrm_adealloc',9999)
128 call qrm_err_act_restore(err_act)
132 call qrm_err_act_restore(err_act)
133 if(err_act .eq. qrm_abort_)
then Generic interface for the qrm_adealloc_i, qrm_adealloc_2i, qrm_adealloc_s, qrm_adealloc_2s, qrm_adealloc_3s, qrm_adealloc_d, qrm_adealloc_2d, qrm_adealloc_3d, qrm_adealloc_c, qrm_adealloc_2c, qrm_adealloc_3c, qrm_adealloc_z, qrm_adealloc_2z, qrm_adealloc_3z, routines.
This module contains the interfaces of all non-typed routines.
This module contains the definition of the basic sparse matrix type and of the associated methods...
Generic interface for the qrm_aalloc_i, qrm_aalloc_2i, qrm_aalloc_s, qrm_aalloc_2s, qrm_aalloc_3s, qrm_aalloc_d, qrm_aalloc_2d, qrm_aalloc_3d, qrm_aalloc_c, qrm_aalloc_2c, qrm_aalloc_3c, qrm_aalloc_z, qrm_aalloc_2z, qrm_aalloc_3z, routines.
subroutine dqrm_rowperm(graph, cperm, rperm, nvar, stair)
This routine computes a row permutation that puts the input matrix into a staircase format...
This type defines the data structure used to store a matrix.
This module implements the memory handling routines. Pretty mucch allocations and deallocations...