QR_MUMPS
qrm_elim_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 !! ##############################################################################################
33 
34 
35 #include "qrm_common.h"
36 
37 
40 
60 subroutine _qrm_elim_tree(graph, cperm, parent)
61 
62  use _qrm_spmat_mod
63  use qrm_mem_mod
64  implicit none
65 
66  type(_qrm_spmat_type), intent(in) :: graph
67  integer, intent(in) :: cperm(:)
68  integer, allocatable :: parent(:)
69 
70  integer, allocatable :: ancestor(:), prev_col(:)
71  integer :: iidx, jidx, i, j, k
72  ! error management
73  integer :: err_act
74  character(len=*), parameter :: name='qrm_elim_tree'
75 
76  call qrm_err_act_save(err_act)
77 
78  call qrm_aalloc(ancestor, graph%n)
79  call qrm_aalloc(prev_col, graph%m)
80  ! realloc parent in case it was not allocated before
81  call qrm_arealloc(parent, graph%n)
82  __qrm_check_ret(name,'qrm_aalloc',9999)
83 
84  ! ancestor is used for doing path compression
85  ancestor=0
86  prev_col=0
87  parent =0
88 
89  do jidx = 1, graph%n
90  j = cperm(jidx)
91  do iidx = graph%jptr(j), graph%jptr(j+1)-1
92  i = graph%irn(iidx)
93  k = prev_col(i)
94  if(k .ne. 0) call add_node(ancestor, parent, k, j)
95  prev_col(i) = j
96  end do
97  end do
98 
99  call qrm_adealloc(ancestor)
100  call qrm_adealloc(prev_col)
101  __qrm_check_ret(name,'qrm_adealloc',9999)
102 
103  call qrm_err_act_restore(err_act)
104  return
105 
106 9999 continue ! error management
107  call qrm_err_act_restore(err_act)
108  if(err_act .eq. qrm_abort_) then
109  call qrm_err_check()
110  end if
111 
112  return
113 
114 
115 contains
116 
117  subroutine add_node(ancestor, parent, k, j)
118  implicit none
119  integer, allocatable :: parent(:)
120  integer, allocatable :: ancestor(:)
121  integer :: k, j
122 
123  integer :: r
124 
125  continue
126  do
127  r = ancestor(k)
128 
129  ancestor(k)=j
130 
131  if (r .eq. j) then
132  ! node was already added, do nothing
133  return
134  end if
135 
136  ancestor(k)=j
137 
138  if(r .eq. 0) then
139  ! we reached the top of the subtree
140  ! add j on top
141  parent(k) = j
142  return
143  end if
144 
145  k = r
146 
147  end do
148 
149 
150  end subroutine add_node
151 
152 end subroutine _qrm_elim_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.
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 type defines the data structure used to store a matrix.
This module contains the definition of the basic sparse matrix type and of the associated methods...
Generic interface for the qrm_arealloc_i qrm_arealloc_s qrm_arealloc_d qrm_arealloc_c qrm_arealloc_z...
subroutine _qrm_elim_tree(graph, cperm, parent)
This subroutine builds the elimination tree for A'A.
subroutine add_node(ancestor, parent, k, j)
This module implements the memory handling routines. Pretty mucch allocations and deallocations...
Definition: qrm_mem_mod.F90:38