35 #include "qrm_common.h"
90 integer,
allocatable :: scol(:), srow(:)
91 integer,
allocatable :: mcperm(:), mrperm(:)
92 integer :: ncsing, nrsing
97 integer,
allocatable :: nel(:), isrow(:), tmp(:)
98 integer :: j, row, col, jj, ii,
i, curr, ncol_p, nrow_p, nz_p, nz_save, idx
101 character(len=*),
parameter :: name=
'qrm_detect_singletons'
116 call move_alloc(mcperm, nel)
117 __qrm_check_ret(name,
'qrm_arealloc',9999)
132 nel(j) = graph%jptr(j+1)-graph%jptr(j)
133 if (nel(j) .eq. 0)
then
138 else if (nel(j) .eq. 1)
then
142 isrow(ncsing) = graph%irn(graph%jptr(j))
146 if((ncsing .eq. 0) .or. (ncsing .eq. graph%n))
then
149 __qrm_check_ret(name,
'qrm_dealloc',9999)
150 call move_alloc(nel, mcperm)
156 __qrm_check_ret(name,
'qrm_spmat_convert',9999)
161 if(curr .gt. ncsing)
exit
164 if ((row .eq. 0))
then
168 if(graph_t%iptr(row) .lt.0)
then
175 graph_t%iptr(row) =
flip(graph_t%iptr(row))
177 do jj=
unflip(graph_t%iptr(row)),
unflip(graph_t%iptr(row+1)) -1
181 if (nel(j) .eq. 1)
then
188 do ii=graph%jptr(j), graph%jptr(j+1)-1
190 if(graph_t%iptr(
i) .gt. 0)
exit
201 if(isrow(
i) .ne. 0)
then
203 srow(curr) = isrow(
i)
213 call move_alloc(nel, mcperm)
225 if(graph_t%iptr(
i) .lt. 0) cycle
235 if (mcperm(j) .lt. 0) cycle
239 do ii= graph%jptr(j), graph%jptr(j+1)-1
242 if (graph_t%iptr(
i) .lt. 0) cycle
243 graph%irn(nz_p) = mrperm(
i)
246 graph%jptr(ncol_p) = nz_save
249 graph%jptr(ncol_p+1) = nz_save
252 if((ncsing + ncol_p) .ne. graph%n)
then
253 __qrm_prnt_dbg(
'("Inconsistency in singleton detection")')
259 if(mrperm(
i) .gt. 0)
then
265 mrperm(nrow_p+1:graph_t%m) = 0
268 nrsing = graph%m-nrow_p
275 __qrm_check_ret(name,
'spmat_destroy/dealloc',9999)
278 tmp(1:ncsing) = scol(1:ncsing)
280 call move_alloc(tmp, scol)
283 tmp(1:nrsing) = srow(1:nrsing)
285 call move_alloc(tmp, srow)
293 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.
subroutine dqrm_detect_singletons(graph, scol, srow, mrperm, mcperm, nrsing, ncsing)
This subroutine detects singletons in a matrix.
This module contains the interfaces of all non-typed routines.
subroutine qrm_err_act_restore(err_act)
Restores the value of the qrm_err_act variable.
integer function unflip(i)
This module contains the definition of the basic sparse matrix type and of the associated methods...
subroutine qrm_err_check()
This subroutine checks the errors stack. If something is found all the entries in the stack are poppe...
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_spmat_destroy(qrm_spmat, all)
This subroutine destroyes a qrm_spmat instance.
subroutine dqrm_spmat_convert(in_mat, out_mat, fmt, values)
This subroutine converts an input matrix into a different storage format. Optionally the values may b...
This type defines the data structure used to store a matrix.
Generic interface for the qrm_arealloc_i qrm_arealloc_s qrm_arealloc_d qrm_arealloc_c qrm_arealloc_z...
This module implements the memory handling routines. Pretty mucch allocations and deallocations...
subroutine qrm_err_act_save(err_act)
Saves a copy of the qrm_err_act variable.