QR_MUMPS
 All Classes Files Functions Variables Enumerations Enumerator Pages
qrm_analyse.F90
Go to the documentation of this file.
1 !! ##############################################################################################
2 !!
3 !! Copyright 2012 CNRS, INPT
4 !!
5 !! This file is part of qr_mumps.
6 !!
7 !! qr_mumps is free software: you can redistribute it and/or modify
8 !! it under the terms of the GNU Lesser General Public License as
9 !! published by the Free Software Foundation, either version 3 of
10 !! the License, or (at your option) any later version.
11 !!
12 !! qr_mumps is distributed in the hope that it will be useful,
13 !! but WITHOUT ANY WARRANTY; without even the implied warranty of
14 !! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 !! GNU Lesser General Public License for more details.
16 !!
17 !! You can find a copy of the GNU Lesser General Public License
18 !! in the qr_mumps/doc directory.
19 !!
20 !! ##############################################################################################
21 
22 
23 !! ##############################################################################################
33 
34 
35 #include "qrm_common.h"
36 
38 
64 subroutine _qrm_analyse(qrm_mat, transp)
65 
66  use _qrm_spmat_mod
67  use qrm_common_mod
68  use qrm_error_mod
69  use _qrm_utils_mod
70  use qrm_string_mod
71  use _qrm_analysis_mod, savesym => _qrm_analyse
72 
73  implicit none
74 
75  type(_qrm_spmat_type), target :: qrm_mat
76  character, optional, intent(in) :: transp
77 
78  type(_qrm_spmat_type), target :: graph
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
86  logical :: sing
87  character :: itransp
88  ! error management
89  integer :: err_act
90  character(len=*), parameter :: name='qrm_analyse'
91 
92  call qrm_err_act_save(err_act)
93 
94  __qrm_prnt_dbg('("Entering the analysis driver")')
95 
96  call _qrm_check_spmat(qrm_mat, qrm_analyse_)
97  __qrm_check_ret(name,'qrm_check_spmat',9999)
98 
99  call qrm_adata_destroy(qrm_mat%adata)
100  __qrm_check_ret(name,'qrm_adata_destroy',9999)
101 
102  ncsing=0
103  nrsing=0
104 
105  nullify(cperm_p,rperm_p,cperm,rperm,tmp)
106 
107  if(present(transp)) then
108  itransp = transp
109  else
110  itransp = 'n'
111  end if
112  ! in case transp=='t' switch temporarily the row and column indices as well as m and n
113  if(qrm_str_tolower(itransp) .eq. 't') then
114  tmp => qrm_mat%irn
115  qrm_mat%irn => qrm_mat%jcn
116  qrm_mat%jcn => tmp
117  i = qrm_mat%m
118  qrm_mat%m = qrm_mat%n
119  qrm_mat%n = i
120  end if
121 
122  call _qrm_compute_graph(qrm_mat, graph)
123  __qrm_check_ret(name,'qrm_compute_graph',9999)
124 
125  call qrm_get(qrm_mat, 'qrm_sing', i)
126  sing = i .eq. qrm_yes_
127 
128  if(sing) then
129  ! do singleton detection and prune graph
130  call _qrm_detect_singletons(graph, scol, srow, mrperm, mcperm, nrsing, ncsing)
131  __qrm_check_ret(name,'qrm_detect_singletons',9999)
132 
133  __qrm_prnt_dbg('("row/column singletons ",i10,x,i10)')nrsing,ncsing
134 
135  else
136  ncsing = 0
137  nrsing = 0
138  end if
139 
140  ! the cperm array will contain
141  call qrm_aalloc(graph%adata%cperm, graph%n) ! FIXME: can't rememebr why n+1, check
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
146 
147  graph%adata%ncsing = ncsing
148  graph%adata%nrsing = nrsing
149 
150  ! shortcut: in case a simple permutation is enough to have a triangular
151  ! matrix
152  if ( (ncsing .eq. qrm_mat%n) .or. (nrsing .eq. qrm_mat%m) ) goto 9998 ! FIXME: goto 10 is not the right thing to do
153 
154  ! at this point we the graph and the singletons. Time to go for the ordering
155  ! the ordering will be computed inside cperm(ncsing+1:n) and will be mapped to
156  ! the columns of A only later using mcperm
157  call _qrm_do_ordering(graph, cperm, qrm_mat%cperm_in)
158  __qrm_check_ret(name,'qrm_do_ordering',9998)
159 
160  ! build the elimination tree
161  call qrm_aalloc(parent, graph%n)
162  __qrm_check_ret(name,'qrm_aalloc',9998)
163 
164  call _qrm_elim_tree(graph, cperm, parent)
165  __qrm_check_ret(name,'qrm_elim_tree',9998)
166 
167  ! compute a postorder traversal of the tree
168  call qrm_postorder(parent, graph%n, cperm)
169  __qrm_check_ret(name,'qrm_postorder',9998)
170 
171  call qrm_aalloc(rc, graph%n)
172  __qrm_check_ret(name,'qrm_aalloc',9998)
173 
174  ! do the symbolic facto
175  call _qrm_rowcount(graph, parent, cperm, rc)
176  __qrm_check_ret(name,'qrm_rowcount',9998)
177 
178  call qrm_aalloc(nvar, graph%n)
179  __qrm_check_ret(name,'qrm_aalloc',9998)
180 
181  ! amalgamate the tree
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)
184  ! call qrm_prnt_array(nvar(1:qrm_mat%n), 'nvar')
185 
186  call qrm_aalloc(stair, graph%n)
187  __qrm_check_ret(name,'qrm_aalloc',9998)
188 
189  ! Compute the row permutation to put the matrix in staircase form
190  call _qrm_rowperm(graph, cperm, rperm, nvar, stair)
191  __qrm_check_ret(name,'qrm_rowperm',9998)
192 
193  ! Rework all the data and compress it to nnodes size instead of n
194  call qrm_compress_data(graph%adata, cperm, parent, rc, stair, graph%n)
195  __qrm_check_ret(name,'qrm_compress_data',9998)
196 
197  ! Do the symbolic facto: estimate fronts size, flops, parallelism etc.
198  call _qrm_symbolic(graph)
199  __qrm_check_ret(name,'qrm_symbolic',9998)
200 
201  ! reorder the tree in order to minimize something
202  call qrm_reorder_tree(graph%adata)
203 
204  if (ncsing .gt. 0) then
205  ! the result of the ordering, symbolic facto and amalgamation
206  ! has to be mapped back to the original matrix
207  call _qrm_attach_singletons(graph, scol, srow, mrperm, mcperm, nrsing, ncsing)
208  __qrm_check_ret(name,'qrm_attach_singletons',9998)
209 
210  end if
211 
212  call qrm_adata_move(graph%adata, qrm_mat%adata)
213  qrm_mat%gstats(1:3) = graph%gstats(1:3)
214 
215  ! compute the number of rows in each front and some global stats
216  ! call qrm_prnt_array(qrm_mat%adata%cperm(1:qrm_mat%n), 'cp',40)
217  ! call qrm_prnt_array(qrm_mat%adata%cp_ptr(1:qrm_mat%adata%nnodes+1), 'pt')
218  ! call qrm_prnt_array(qrm_mat%adata%rperm(1:qrm_mat%m), 'rp',30)
219 
220 
221 9998 continue
222 
223  ! in case transp=='t' restore the row and column indices as well as m and n
224  if(qrm_str_tolower(itransp) .eq. 't') then
225  tmp => qrm_mat%irn
226  qrm_mat%irn => qrm_mat%jcn
227  qrm_mat%jcn => tmp
228  i = qrm_mat%m
229  qrm_mat%m = qrm_mat%n
230  qrm_mat%n = i
231  end if
232 
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)
241  call qrm_spmat_destroy(graph, all=.true.)
242  __qrm_check_ret(name,'qrm_dealloc',9999)
243 
244 
245  qrm_mat%adata%ok = .true.
246  call qrm_err_act_restore(err_act)
247  return
248 
249 9999 continue ! error management
250  call qrm_err_act_restore(err_act)
251  if(err_act .eq. qrm_abort_) then
252  call qrm_err_check()
253  end if
254 
255  return
256 
257 end subroutine _qrm_analyse
258 
259 
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...
Definition: qrm_rowperm.F90:68
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.
Definition: qrm_analyse.F90:64
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.
int i
Definition: secs.c:40
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.