QR_MUMPS
 All Classes Files Functions Variables Enumerations Enumerator Pages
qrm_compress_data.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 !! ##############################################################################################
35 
36 
37 #include "qrm_common.h"
38 
39 
43 
44 !! The arrays on input are all of size n (where n is the column size of
45 !! the original matrix). The purpose of this subroutine is to store the
46 !! same information into arrays of size nsteps, where nsteps is the
47 !! number of nodes in the assembly tree. Only porder will remain of size n
48 !!
49 !! @param[out] adata on output the following fields will be
50 !! modified:
51 !! - cp_ptr
52 !! - rc
53 !! - parent
54 !! - nnodes
55 !! - stair
56 !! - child. for each node, it contains the list of its children
57 !! - childptr. pointers to the list of children. chil(childptr(i),...,childptr(i+1)-1)
58 !! contains all the children of node i
59 !! - icperm. the inverse column permutation is built
60 !!
61 !! @param[in] porder contains the postorder. porder(i)=k means that the i-th column
62 !! in the computed ordering is column k in the original matrix.
63 !! Inthis postorder principal variables always come before the
64 !! corresponding subordinate variables.
65 !!
66 !! @param[in] parent contains the assembly tree. parent(i)=k:
67 !! - k>0: means that i is a principal variable in a node and k
68 !! is the principal variable in the father's node
69 !! - k<0: means that i is a subordinate variable in a supernode
70 !! whose principal variable is k.
71 !!
72 !! @param[in] rc this array contains the rowcount, i.e., rc(i)=k means that
73 !! in the R factor the rows corresponding to the node whose
74 !! principal variable is i have k nonzeroes. rc(i)=0 for all
75 !! subordinate varibales and rc(i)=-1 for all the column singletons
76 !!
77 !! @param[in] stair stair(i) contains the number of rows in the step related to
78 !! node i. see @link ::_qrm_rowperm_ @endlink
79 !!
80 !! @param[in] n the number of columns in the original matrix (i.e. the size of the
81 !! input arrays)
82 !!
83 subroutine qrm_compress_data(adata, porder, parent, rc, stair, n)
84 
85  use qrm_mem_mod
86  use qrm_adata_mod
87  use qrm_common_mod, protect => qrm_compress_data
88  implicit none
89 
90  type(qrm_adata_type) :: adata
91  integer :: porder(:), parent(:), rc(:), stair(:)
92  integer :: n
93 
94  integer :: i, pnt, svar, f, ss
95  integer, allocatable :: work(:)
96  ! error management
97  integer :: err_act
98  character(len=*), parameter :: name='qrm_compress_data'
99 
100  call qrm_err_act_save(err_act)
101 
102  adata%nnodes = 0
103  do i=1, n
104  if(rc(i) .ne. 0) adata%nnodes = adata%nnodes+1
105  end do
106 
107  call qrm_aalloc(adata%cp_ptr, adata%nnodes+1)
108  call qrm_aalloc(work, n)
109  __qrm_check_ret(name,'qrm_aalloc',9999)
110 
111  work = 0
112  adata%cp_ptr=0
113  ! build pointers and inverse mapping
114  ! work(i)=k means that principal variable i is in node number k
115  adata%cp_ptr(1) = 1
116  work(porder(1)) = 1
117  pnt = 2
118  do i=2, n
119  svar=porder(i)
120  if (rc(svar) .ne. 0) then
121  adata%cp_ptr(pnt) = i
122  work(svar) = pnt
123  pnt = pnt+1
124  end if
125  end do
126  adata%cp_ptr(adata%nnodes+1) = n+1
127 
128  ! build adata%parent
129  ! build adata%rc
130  ! build adata%child and adata%childptr
131  call qrm_aalloc(adata%parent, adata%nnodes)
132  call qrm_aalloc(adata%rc, adata%nnodes)
133  call qrm_aalloc(adata%stair, adata%nnodes)
134  call qrm_aalloc(adata%child, adata%nnodes)
135  call qrm_aalloc(adata%childptr, adata%nnodes+1)
136  call qrm_aalloc(adata%icperm, n)
137  __qrm_check_ret(name,'qrm_aalloc',9999)
138 
139  adata%childptr = 0
140  ss = 0
141  do i=1, adata%nnodes
142  svar = porder(adata%cp_ptr(i))
143  f = parent(svar)
144  adata%rc(i) = rc(svar)
145  ! adata%nfrows(i) = stair(svar)-ss
146  adata%stair(i) = stair(svar)
147  ss = adata%stair(i)
148  if(f .eq. 0) then
149  adata%parent(i) = 0
150  else
151  adata%parent(i) = work(f)
152  adata%childptr(work(f)+1) = adata%childptr(work(f)+1)+1
153  end if
154  end do
155 
156  adata%childptr(1) = 1
157  do i=2, adata%nnodes+1
158  adata%childptr(i) = adata%childptr(i-1)+adata%childptr(i)
159  end do
160 
161 
162  adata%child=0
163  work(1:adata%nnodes)=0
164  do i=1, adata%nnodes
165  f = adata%parent(i)
166  if (f .ne. 0) then
167  adata%child(adata%childptr(f)+work(f)) = i
168  work(f) = work(f)+1
169  end if
170  end do
171 
172  do i=1, n
173  adata%icperm(adata%cperm(i)) = i
174  end do
175 
176 
177  call qrm_adealloc(work)
178  __qrm_check_ret(name,'qrm_adealloc',9999)
179 
180  call qrm_err_act_restore(err_act)
181  return
182 
183 9999 continue ! error management
184  call qrm_err_act_restore(err_act)
185  if(err_act .eq. qrm_abort_) then
186  call qrm_err_check()
187  end if
188 
189  return
190 
191 end subroutine qrm_compress_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.
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.
This module contains the definition of the analysis data type.
subroutine qrm_err_check()
This subroutine checks the errors stack. If something is found all the entries in the stack are poppe...
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
Generic interface for the ::qrm_compress_data routine.
This module implements the memory handling routines. Pretty mucch allocations and deallocations...
Definition: qrm_mem_mod.F90:38
subroutine qrm_compress_data(adata, porder, parent, rc, stair, n)
This routine compresses the results of a number of operations in the analysis phase. Basically, the input data is of size n and the output of size adatannodes (which is the number of nodes in the elimination tree).
int i
Definition: secs.c:40
subroutine qrm_err_act_save(err_act)
Saves a copy of the qrm_err_act variable.