QR_MUMPS
 All Classes Files Functions Variables Enumerations Enumerator Pages
qrm_amalg_tree.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 !! ##############################################################################################
34 
35 
36 #include "qrm_common.h"
37 
39 
75 subroutine qrm_amalg_tree(n, parent, rowcount, porder, nvar, min_var, fill_thresh)
76 
77  use qrm_error_mod
78  use qrm_mem_mod
79  use qrm_common_mod, savesym => qrm_amalg_tree
80  implicit none
81 
82  integer :: n, min_var
83  integer :: parent(:), rowcount(:), porder(:), nvar(:)
84  real(kind(1.d0)) :: fill_thresh
85 
86  integer, allocatable :: snodes_ptr(:) ! pointers to the start of a supernode in porder
87  integer :: nsnodes, i, j, k, f
88  integer, allocatable :: snodes_map(:), snodes_parent(:), kids(:), &
89  & snodes_merged(:), snodes_nvar(:), snodes_rc(:)
90  real(kind(1.d0)) :: ratio
91  integer(kind=8), allocatable :: snodes_fill(:)
92  integer(kind=8) :: totfill, msize, fill
93  integer :: pv, fpv, rc, frc, col, fcol, tcol, s, ns
94  logical :: merge
95  ! error management
96  integer :: err_act
97  character(len=*), parameter :: name='qrm_amalg_tree'
98 
99  call qrm_err_act_save(err_act)
100 
101  ! sort nodes children in increasing row-count
102  call qrm_postorder(parent, n, porder, rowcount)
103  __qrm_check_ret(name,'qrm_postorder',9999)
104 
105  ! fundamental supernodes detection
106  call qrm_aalloc(snodes_ptr, n+1)
107  call qrm_aalloc(kids, n)
108 
109  kids = 0
110  ! count the number of children for each snode
111  do i=1, n
112  f = parent(i)
113  if(f .ne. 0) kids(f) = kids(f) + 1
114  end do
115 
116  ! detect snodes
117  nsnodes = 1
118  snodes_ptr(1) = 1
119  do k=2, n
120  i = porder(k)
121  j = porder(k-1)
122  if ( (parent(j) .ne. i) .or. &
123  & (rowcount(j) .ne. rowcount(i)+1) .or. &
124  & (kids(i) .gt. 1)) then
125  nsnodes = nsnodes+1
126  snodes_ptr(nsnodes) = k
127  end if
128  end do
129  snodes_ptr(nsnodes+1)=n+1
130 
131  ! compress the tree
132  ! note that to reduce memory consumption we could
133  ! do the amalgamantion on the uncompressed tree TODO
134  call move_alloc(kids, snodes_map)
135  call qrm_aalloc(snodes_parent, nsnodes)
136  call qrm_aalloc(snodes_nvar, nsnodes)
137  call qrm_aalloc(snodes_rc, nsnodes)
138 
139  snodes_nvar = 0
140 
141  do i=nsnodes, 1, -1
142  pv = porder(snodes_ptr(i))
143  snodes_nvar(i) = snodes_ptr(i+1)-snodes_ptr(i)
144  snodes_rc(i) = rowcount(pv)
145  do j=snodes_ptr(i), snodes_ptr(i+1)-1
146  ! map(i) = j means that j belongs to snode i
147  snodes_map(porder(j)) = i
148  end do
149  f = parent(porder(j-1)) ! the father of the last node in the snode
150  if(f .eq. 0) then
151  snodes_parent(i) = 0
152  else
153  snodes_parent(i) = snodes_map(f)
154  end if
155  end do
156 
157 
158  allocate(snodes_fill(nsnodes))
159  snodes_fill = 0
160  call move_alloc(snodes_map, snodes_merged)
161  snodes_merged = 0
162 
163 
164  snodes_loop: do s = nsnodes-1, 1, -1
165 
166  ! detect the real father of s
167  f = snodes_parent(s)
168  if(f .ne. 0) then
169  do
170  if(snodes_merged(f) .eq. 0) exit
171  f = snodes_merged(f)
172  end do
173 
174  ! path compression on the chain of amalgamated ancestors
175  i = snodes_parent(s)
176  do
177  if(i .eq. f) exit
178  j = snodes_merged(i)
179  snodes_merged(i) = f
180  i = j
181  end do
182  end if
183 
184  if(f .ne. s+1) then
185  ! don't merge if snode i+1 is not snode i's father
186  merge = .false.
187  cycle snodes_loop
188  end if
189 
190  rc = snodes_rc(s) ! the row count of the current snode
191  frc = snodes_rc(f) ! the row count of the father snode
192 
193  col = snodes_nvar(s) ! the variables in the current snode
194  fcol = snodes_nvar(f) ! the variables in the parent snode
195  tcol = col+fcol
196  fill = col*(frc+col-rc) ! the fill added in the case where the nodes are merged
197  totfill = snodes_fill(f) + fill
198 
199  if(tcol .lt. min_var) then
200  merge = .true.
201  else if (fill .eq. 0) then
202  merge = tcol.lt.min_var*50
203  else
204  msize = (tcol)*(tcol+1)/2 + tcol*(frc-fcol)
205 
206  ratio = real(totfill, kind(1.d0))/real(msize,kind(1.d0))
207 
208  merge = (tcol .lt. min_var*4 .and. ratio .lt. fill_thresh*16) .or.&
209  & (tcol .lt. min_var*12 .and. ratio .lt. fill_thresh*2) .or. &
210  & (ratio .lt. fill_thresh)
211 
212  end if
213 
214  if(merge) then
215  ! the two nodes ca be merged
216  snodes_fill(s) = totfill
217  snodes_nvar(s) = tcol
218  snodes_nvar(f) = 0
219  snodes_merged(f) = s
220  snodes_rc(s) = col+frc
221  snodes_rc(f) = 0
222  end if
223 
224 
225  end do snodes_loop
226 
227  nvar=0
228  ! rebuild snodes_ptr
229  ns = 1
230  snodes_ptr(1) = 1
231  do i=1, nsnodes
232  pv = porder(snodes_ptr(i))
233  if(snodes_merged(i) .eq. 0) then
234  ns = ns+1
235  rowcount(pv) = snodes_rc(i)
236  nvar(pv) = snodes_nvar(i)
237  snodes_ptr(ns) = snodes_ptr(ns-1)+snodes_nvar(i)
238  if(snodes_ptr(ns).gt.n) exit
239  end if
240  end do
241  ns = ns-1
242 
243  call move_alloc(snodes_merged, snodes_map )
244  do i=ns, 1, -1
245  pv = porder(snodes_ptr(i))
246  snodes_map(pv) = pv
247  f = parent(porder(snodes_ptr(i+1)-1)) ! the father of the last node in the snode
248  do j=snodes_ptr(i)+1, snodes_ptr(i+1)-1
249  ! map(i) = j means that j is the principal variable of the snode where i belongs
250  snodes_map(porder(j)) = pv
251  rowcount(porder(j)) = 0
252  parent(porder(j)) = -pv
253  end do
254  if(f .eq. 0) then
255  parent(pv) = 0
256  else
257  parent(pv) = snodes_map(f)
258  end if
259  end do
260 
261  call qrm_adealloc(snodes_map)
262  call qrm_adealloc(snodes_nvar)
263  call qrm_adealloc(snodes_rc)
264  call qrm_adealloc(snodes_ptr)
265  call qrm_adealloc(snodes_parent)
266 
267  call qrm_err_act_restore(err_act)
268  return
269 
270 9999 continue ! error management
271  call qrm_err_act_restore(err_act)
272  if(err_act .eq. qrm_abort_) then
273  call qrm_err_check()
274  end if
275 
276  return
277 
278 end subroutine qrm_amalg_tree
This module contains all the error management routines and data.
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 qrm_amalg_tree(n, parent, rowcount, porder, nvar, min_var, fill_thresh)
This subroutine performs amalgamation on the assembly tree.
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.
* s
Definition: secs.c:43
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.
Definition: qrm_mem_mod.F90:78
Generic interface for the ::qrm_amalg_tree routine.
Generic interface for the ::qrm_postorder routine.
This module implements the memory handling routines. Pretty mucch allocations and deallocations...
Definition: qrm_mem_mod.F90:38
int i
Definition: secs.c:40
subroutine qrm_err_act_save(err_act)
Saves a copy of the qrm_err_act variable.