QR_MUMPS
 All Classes Files Functions Variables Enumerations Enumerator Pages
qrm_postorder.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 
49 subroutine qrm_postorder(parent, n, porder, weight)
50 
51  use qrm_mem_mod
52  implicit none
53 
54  integer :: n
55  integer :: parent(:), porder(:)
56  integer, optional :: weight(:)
57 
58  integer, allocatable :: son(:), brother(:), stack(:)
59  integer :: i, father, br, head, hp, pp, t, w, next
60  ! error management
61  integer :: err_act
62  character(len=*), parameter :: name='qrm_analyse'
63 
64  call qrm_err_act_save(err_act)
65 
66  call qrm_aalloc(son, n)
67  call qrm_aalloc(brother, n)
68  call qrm_aalloc(stack, n)
69  __qrm_check_ret(name,'qrm_aalloc',9999)
70 
71  son = 0
72 
73  ! build a tree that can be traversed top-to-bottom
74  if(present(weight)) then
75  ! use stack as a workspace
76  stack = 0
77  do i=1, n
78  w = weight(i)
79  brother(i) = stack(w)
80  stack(w) = i
81  end do
82 
83  do w=n, 1, -1
84  i = stack(w)
85  do
86  if(i .eq. 0) exit
87  next = brother(i)
88  father = parent(i)
89  if(father .ne. 0) then
90  brother(i) = son(father)
91  son(father) = i
92  end if
93  i = next
94  end do
95  end do
96  else
97  do i=n, 1, -1
98  father = parent(i)
99  if (father .ne. 0) then
100  br = son(father)
101  brother(i) = br
102  son(father) = i
103  end if
104  end do
105  end if
106 
107 
108  head = 0
109  hp = 0
110  pp = 1
111  ! the tree is processed in dfs. Starting from a root, we go down
112  ! and put all the encountered nodes on a stack. When we reach the bottom,
113  ! we pop the leaf from the stack, we put it in the postorder and we go down again
114  ! along the branch starting from the brother of the node just popped
115  do t=1, n
116  if (parent(t) .ne. 0) cycle
117  ! at this point t is the root of a tree
118  hp = hp+1
119  stack(hp) = t
120  head = t
121  do
122  if(son(head) .eq. 0) then
123  ! we reached the bottom
124  porder(pp) = head
125  pp = pp+1
126  hp = hp-1
127  if (parent(head) .ne. 0) then
128  son(parent(head)) = brother(head)
129  end if
130  if (hp .eq. 0) then
131  ! tree is exhausted
132  exit
133  end if
134  head = stack(hp)
135  else
136  ! go down one more level
137  hp = hp+1
138  stack(hp) = son(head)
139  head = son(head)
140  end if
141  end do
142  end do
143 
144  call qrm_adealloc(son)
145  call qrm_adealloc(brother)
146  call qrm_adealloc(stack)
147  __qrm_check_ret(name,'qrm_adealloc',9999)
148 
149  call qrm_err_act_restore(err_act)
150  return
151 
152 9999 continue ! error management
153  call qrm_err_act_restore(err_act)
154  if(err_act .eq. qrm_abort_) then
155  call qrm_err_check()
156  end if
157 
158  return
159 
160 end subroutine qrm_postorder
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_err_act_restore(err_act)
Restores the value of the qrm_err_act variable.
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
This module implements the memory handling routines. Pretty mucch allocations and deallocations...
Definition: qrm_mem_mod.F90:38
subroutine qrm_postorder(parent, n, porder, weight)
This subroutine computes a postorder by traversing a tree in dfs.
int i
Definition: secs.c:40
subroutine qrm_err_act_save(err_act)
Saves a copy of the qrm_err_act variable.