QR_MUMPS
qrm_reorder_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 
38 
51 subroutine qrm_reorder_tree(adata)
52 
53  use qrm_adata_mod
54  use qrm_mem_mod
55  use qrm_sort_mod
56  implicit none
57 
58  type(qrm_adata_type) :: adata
59 
60  integer, allocatable :: peaks(:), children_peaks(:), aux(:), &
61  & roots(:), child(:), tnch(:)
62  integer :: root, maxch, nroots, nleaves, i, nl, leaf, f, c
63  integer :: nch, nnodes, node, peak
64 
65  call qrm_aalloc(peaks, adata%nnodes)
66 
67  nroots = 0
68  maxch = 0
69  nleaves = 0
70  do node = 1, adata%nnodes
71  if(adata%parent(node) .eq. 0) then
72  nroots = nroots+1
73  end if
74 
75  nch = adata%childptr(node+1)-adata%childptr(node)
76  if(nch .gt. maxch) maxch = nch
77  if(nch .eq. 0 ) nleaves = nleaves+1
78  peaks(node) = -nch
79  end do
80 
81  ! uncomment to bypass the reordertree
82  ! adata%nleaves = 0
83  ! call qrm_aalloc(adata%leaves, nleaves)
84  ! do node = 1, adata%nnodes
85  ! nch = adata%childptr(node+1)-adata%childptr(node)
86  ! if(nch .eq. 0) then
87  ! adata%nleaves = adata%nleaves+1
88  ! adata%leaves(adata%nleaves) = node
89  ! end if
90  ! end do
91  ! return
92  ! until here
93 
94 
95 
96  maxch = max(maxch, nroots)
97  call qrm_aalloc(children_peaks, maxch)
98  call qrm_aalloc(aux, maxch+2, lbnd=0)
99  call qrm_aalloc(roots, nroots)
100 
101  nnodes = adata%nnodes
102  node = 1
103  nl = nleaves
104  nroots = 0
105 
106  main: do leaf=1, adata%nnodes
107 
108  nch = adata%childptr(leaf+1)-adata%childptr(leaf)
109  if(nch .eq. 0) then
110  node = leaf
111 
112  ! start going up the subtree
113  do
114  if(nnodes .eq. 0) exit main
115  ! peaks(node).eq.0 only occurs if node is a leaf or
116  ! if its peak can be computed
117  if(peaks(node) .ne. 0) exit
118 
119  peak=0
120  ! sort node's children (if any) and compute its peak
121  nch = 0
122  if(adata%childptr(node) .ne. adata%childptr(node+1)) then
123  do i = adata%childptr(node), adata%childptr(node+1)-1
124  nch = nch+1
125  children_peaks(nch) = peaks(adata%child(i))
126  end do
127  ! sort
128  call qrm_mergesort(nch, children_peaks(1:nch), aux(0:nch+1), order=-1)
129 
130  call qrm_mergeswap(nch, aux(0:nch+1), children_peaks(1:nch), &
131  & adata%child(adata%childptr(node):adata%childptr(node+1)-1))
132  ! compute node's peak
133  do i=1, nch
134  peak = max(peak, i-1+children_peaks(i))
135  end do
136  end if
137 
138  peaks(node) = max(peak,nch+1)
139  nnodes = nnodes -1
140 
141  f = adata%parent(node)
142  if (f .ne. 0) then
143  peaks(f) = peaks(f)+1
144  if(peaks(f) .eq. 0) then
145  node = f
146  else
147  cycle main
148  end if
149  else
150  nroots = nroots+1
151  roots(nroots) = node
152  end if
153  end do
154 
155  end if
156 
157  end do main
158 
159  ! sort roots
160  do i=1, nroots
161  ! write(*,'("Root peak: ",i10)')peaks(roots(i))
162  children_peaks(i) = peaks(roots(i))
163  end do
164 
165  call qrm_mergesort(nroots, children_peaks(1:nroots), aux(0:nroots+1), order=-1)
166  call qrm_mergeswap(nroots, aux(0:nroots+1), children_peaks(1:nroots), &
167  & roots(1:nroots))
168 
169  call qrm_adealloc(children_peaks)
170  call qrm_adealloc(aux)
171  call qrm_aalloc(adata%leaves, nleaves)
172  call qrm_aalloc(tnch, adata%nnodes)
173  call move_alloc(peaks, child)
174  child = 0
175 
176  ! compute the list of leaves to start facto. Even if the whole tree
177  ! is reordered, only the leaves in l0 are returned; this may have to
178  ! be improved. TODO
179 
180  tnch = 0
181  do node=1, adata%nnodes
182  do i=adata%childptr(node), adata%childptr(node+1)-1
183  c = adata%child(i)
184  if(adata%small(c) .eq. 0) tnch(node) = tnch(node)+1
185  end do
186  end do
187 
188  nleaves = 0
189  roots_loop: do root=1, nroots
190  node = roots(root)
191  ! go down the subtree
192 
193  nodes_loop: do
194  if(tnch(node) .eq. 0) then
195  nleaves = nleaves+1
196  adata%leaves(nleaves) = node
197  else
198 
199  nch=adata%childptr(node+1) - adata%childptr(node)
200  do
201  if(child(node) .ge. nch) exit
202  child(node) = child(node)+1
203  c = adata%child(adata%childptr(node)+child(node)-1)
204  if(adata%small(c) .ne. 1) then
205  node = c
206  cycle nodes_loop
207  end if
208  end do
209  end if
210  node = adata%parent(node)
211  if (node .eq. 0) cycle roots_loop
212  cycle nodes_loop
213 
214  end do nodes_loop
215  end do roots_loop
216 
217  adata%nleaves = nleaves
218 
219 
220  call qrm_adealloc(roots)
221  call qrm_adealloc(child)
222  call qrm_adealloc(tnch)
223 
224 
225  return
226 
227 end subroutine qrm_reorder_tree
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 routines for sorting.
This module contains the definition of the analysis data type.
The main data type for the analysis phase.
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
subroutine qrm_reorder_tree(adata)
This subroutine reorders the assembly tree in order to reduce the tasks search space.
This module implements the memory handling routines. Pretty mucch allocations and deallocations...
Definition: qrm_mem_mod.F90:38