35 #include "qrm_common.h"
76 character,
optional,
intent(in) :: transp
79 integer,
allocatable :: parent(:), rc(:), srow(:), scol(:)
80 integer,
allocatable :: mcperm(:), mrperm(:), stair(:), nvar(:)
81 integer :: ncsing, nrsing
82 integer,
pointer :: cperm_p(:)=>null(), rperm_p(:)=>null()
83 integer,
pointer :: cperm(:)=>null(), rperm(:)=>null(), tmp(:)
84 integer ::
i, j, ii, lastrc, cnt, ierr
85 real(kind(1.d0)) :: t1, t2
90 character(len=*),
parameter :: name=
'qrm_analyse'
94 __qrm_prnt_dbg(
'("Entering the analysis driver")')
97 __qrm_check_ret(name,
'qrm_check_spmat',9999)
100 __qrm_check_ret(name,
'qrm_adata_destroy',9999)
105 nullify(cperm_p,rperm_p,cperm,rperm,tmp)
107 if(present(transp))
then
115 qrm_mat%irn => qrm_mat%jcn
118 qrm_mat%m = qrm_mat%n
123 __qrm_check_ret(name,
'qrm_compute_graph',9999)
125 call
qrm_get(qrm_mat,
'qrm_sing',
i)
131 __qrm_check_ret(name,
'qrm_detect_singletons',9999)
133 __qrm_prnt_dbg(
'("row/column singletons ",i10,x,i10)')nrsing,ncsing
141 call qrm_aalloc(graph%adata%cperm, graph%n)
142 call qrm_aalloc(graph%adata%rperm, graph%m)
143 __qrm_check_ret(name,
'qrm_alloc',9999)
144 cperm => graph%adata%cperm
145 rperm => graph%adata%rperm
147 graph%adata%ncsing = ncsing
148 graph%adata%nrsing = nrsing
152 if ( (ncsing .eq. qrm_mat%n) .or. (nrsing .eq. qrm_mat%m) ) goto 9998
158 __qrm_check_ret(name,
'qrm_do_ordering',9998)
161 call qrm_aalloc(parent, graph%n)
162 __qrm_check_ret(name,
'qrm_aalloc',9998)
165 __qrm_check_ret(name,
'qrm_elim_tree',9998)
169 __qrm_check_ret(name,
'qrm_postorder',9998)
171 call qrm_aalloc(rc, graph%n)
172 __qrm_check_ret(name,
'qrm_aalloc',9998)
176 __qrm_check_ret(name,
'qrm_rowcount',9998)
178 call qrm_aalloc(nvar, graph%n)
179 __qrm_check_ret(name,
'qrm_aalloc',9998)
182 call
qrm_amalg_tree(graph%n, parent, rc, cperm, nvar, qrm_mat%icntl(3), qrm_mat%rcntl(1))
183 __qrm_check_ret(name,
'qrm_amalg_tree',9998)
186 call qrm_aalloc(stair, graph%n)
187 __qrm_check_ret(name,
'qrm_aalloc',9998)
191 __qrm_check_ret(name,
'qrm_rowperm',9998)
195 __qrm_check_ret(name,
'qrm_compress_data',9998)
199 __qrm_check_ret(name,
'qrm_symbolic',9998)
204 if (ncsing .gt. 0)
then
208 __qrm_check_ret(name,
'qrm_attach_singletons',9998)
213 qrm_mat%gstats(1:3) = graph%gstats(1:3)
226 qrm_mat%irn => qrm_mat%jcn
229 qrm_mat%m = qrm_mat%n
233 call qrm_adealloc(mcperm)
234 call qrm_adealloc(mrperm)
235 call qrm_adealloc(srow)
236 call qrm_adealloc(scol)
237 call qrm_adealloc(parent)
238 call qrm_adealloc(rc)
239 call qrm_adealloc(nvar)
240 call qrm_adealloc(stair)
242 __qrm_check_ret(name,
'qrm_dealloc',9999)
245 qrm_mat%adata%ok = .true.
251 if(err_act .eq. qrm_abort_)
then
Generif interface for the ::_qrm_spmat_destroy routine.
This module contains all the error management routines and data.
This module contains the interfaces of all non-typed routines.
This module contains the generic interfaces for all the analysis routines.
subroutine _qrm_do_ordering(graph, cperm, cperm_in)
This routine computes (through different methods) a column permutation of the input matrix in order t...
subroutine qrm_adata_move(adata_in, adata_out)
make s a move of an qrm_adata_type instance
subroutine qrm_err_act_restore(err_act)
Restores the value of the qrm_err_act variable.
subroutine _qrm_attach_singletons(qrm_mat, scol, srow, mrperm, mcperm, nrsing, ncsing)
This subroutine merges the results of the singletons detection into the results of the analysis phase...
subroutine _qrm_rowperm(graph, cperm, rperm, nvar, stair)
This routine computes a row permutation that puts the input matrix into a staircase format...
This module contains generic interfaces for a number of auxiliary tools.
subroutine _qrm_detect_singletons(graph, scol, srow, mrperm, mcperm, nrsing, ncsing)
This subroutine detects singletons in a matrix.
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_reorder_tree routine.
subroutine qrm_adata_destroy(adata)
Frees an qrm_adata_type instance.
subroutine _qrm_analyse(qrm_mat, transp)
This is the driver routine for the analysis phase.
subroutine _qrm_compute_graph(qrm_mat, graph)
Computes the adjacency graph of a matrix.
Generic interface for the ::qrm_compress_data routine.
Generic interface for the ::qrm_amalg_tree routine.
subroutine _qrm_rowcount(graph, parent, porder, rc)
This subroutine computes the rowcount of the R factor.
This type defines the data structure used to store a matrix.
Generic interface for the ::qrm_postorder routine.
subroutine _qrm_check_spmat(qrm_spmat, op)
Check the compatibility and correctness of icntl and rcntl parameters.
This module contains the definition of the basic sparse matrix type and of the associated methods...
Generif interface for the ::_qrm_pgeti, ::_qrm_pgetr and.
This module contains various string handling routines.
subroutine _qrm_elim_tree(graph, cperm, parent)
This subroutine builds the elimination tree for A'A.
subroutine _qrm_symbolic(graph)
This subroutine computes the symbolic QR factorization of a matrix.
subroutine qrm_err_act_save(err_act)
Saves a copy of the qrm_err_act variable.