QR_MUMPS
 All Classes Files Functions Variables Enumerations Enumerator Pages
dqrm_attach_singletons.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 
40 
82 subroutine dqrm_attach_singletons(qrm_mat, scol, srow, mrperm, mcperm, nrsing, ncsing)
83 
84  use dqrm_spmat_mod
85  use qrm_adata_mod
86  implicit none
87 
88  type(dqrm_spmat_type) :: qrm_mat
89  integer :: scol(:), srow(:), mcperm(:), mrperm(:)
90  integer :: nrsing, ncsing
91 
92 
93  integer, allocatable :: tmp(:)
94  integer :: m, n, i, nnodes
95 
96  ! restore the values of m and n (to be equal to the original matrix
97  ! dimensions)
98  m = qrm_mat%m+nrsing
99  n = qrm_mat%n+ncsing
100 
101  ! restore the values of nnodes. One extra node will be attached to
102  ! the tree corresponding to the set of singletons
103  nnodes = qrm_mat%adata%nnodes+1
104 
105 
106  ! attach the singleton columns to the column permutation
107  call qrm_aalloc(tmp, n)
108  tmp(1:ncsing) = scol
109  do i=1, qrm_mat%n
110  tmp(ncsing+i) = mcperm(qrm_mat%adata%cperm(i))
111  end do
112  call qrm_adealloc(qrm_mat%adata%cperm)
113  call move_alloc(tmp, qrm_mat%adata%cperm)
114 
115 
116 
117  ! attach the singleton rows to the row permutation
118  call qrm_aalloc(tmp, m)
119  tmp(1:nrsing) = srow
120  do i=1, qrm_mat%m
121  tmp(nrsing+i) = mrperm(qrm_mat%adata%rperm(i))
122  end do
123  call qrm_adealloc(qrm_mat%adata%rperm)
124  call move_alloc(tmp, qrm_mat%adata%rperm)
125 
126 
127 
128  ! attach a new node to the cp_ptr array
129  call qrm_aalloc(tmp, nnodes+1)
130  tmp(1)=1
131  tmp(2:nnodes+1) = qrm_mat%adata%cp_ptr(1:qrm_mat%adata%nnodes+1)+ncsing
132  call qrm_adealloc(qrm_mat%adata%cp_ptr)
133  call move_alloc(tmp, qrm_mat%adata%cp_ptr)
134 
135 
136 
137  ! attach a new node to the parent array
138  call qrm_aalloc(tmp, nnodes)
139  tmp(1)=0
140  do i=1, qrm_mat%adata%nnodes
141  if(qrm_mat%adata%parent(i) .eq. 0) then
142  tmp(i+1)=0
143  else
144  tmp(i+1) = qrm_mat%adata%parent(i)+1
145  end if
146  end do
147  call qrm_adealloc(qrm_mat%adata%parent)
148  call move_alloc(tmp, qrm_mat%adata%parent)
149 
150 
151 
152  ! attach a new node to the rc array
153  call qrm_aalloc(tmp, nnodes)
154  tmp(1)=-1
155  tmp(2:nnodes) = qrm_mat%adata%rc(1:qrm_mat%adata%nnodes)
156  call qrm_adealloc(qrm_mat%adata%rc)
157  call move_alloc(tmp, qrm_mat%adata%rc)
158 
159 
160 
161  ! attach a new node to the stair array
162  call qrm_aalloc(tmp, nnodes)
163  tmp(1)=nrsing
164  tmp(2:nnodes) = qrm_mat%adata%stair(1:qrm_mat%adata%nnodes)+nrsing
165  call qrm_adealloc(qrm_mat%adata%stair)
166  call move_alloc(tmp, qrm_mat%adata%stair)
167 
168 
169 
170  ! attach a new node to the childptr array and fix the child array
171  call qrm_aalloc(tmp, nnodes+1)
172  tmp(1)=1
173  tmp(2:nnodes+1) = qrm_mat%adata%childptr(1:qrm_mat%adata%nnodes+1)
174  qrm_mat%adata%child = qrm_mat%adata%child+1
175  call qrm_adealloc(qrm_mat%adata%childptr)
176  call move_alloc(tmp, qrm_mat%adata%childptr)
177 
178 
179 
180  ! attach a new node to the nfrows array
181  call qrm_aalloc(tmp, nnodes)
182  tmp(1)=nrsing
183  tmp(2:nnodes) = qrm_mat%adata%nfrows(1:qrm_mat%adata%nnodes)
184  call qrm_adealloc(qrm_mat%adata%nfrows)
185  call move_alloc(tmp, qrm_mat%adata%nfrows)
186 
187 
188 
189  ! attach a new node to the fcol and fcol_ptr arrays
190  call qrm_aalloc(tmp, nnodes+1)
191  tmp(1)=1
192  tmp(2:nnodes+1) = qrm_mat%adata%fcol_ptr(1:qrm_mat%adata%nnodes+1)
193  call qrm_adealloc(qrm_mat%adata%fcol_ptr)
194  call move_alloc(tmp, qrm_mat%adata%fcol_ptr)
195  do i=1, qrm_mat%adata%fcol_ptr(nnodes+1)-1
196  qrm_mat%adata%fcol(i) = mcperm(qrm_mat%adata%fcol(i))
197  end do
198 
199 
200 
201  ! fix the leaves and small arrays
202  qrm_mat%adata%leaves = qrm_mat%adata%leaves+1
203  call qrm_aalloc(tmp, nnodes)
204  tmp(1)=0
205  tmp(2:nnodes) = qrm_mat%adata%small(1:qrm_mat%adata%nnodes)
206  call qrm_adealloc(qrm_mat%adata%small)
207  call move_alloc(tmp, qrm_mat%adata%small)
208 
209 
210 
211  ! recompute icperm
212  call qrm_arealloc(qrm_mat%adata%icperm, n)
213  do i=1, n
214  qrm_mat%adata%icperm(qrm_mat%adata%cperm(i)) = i
215  end do
216 
217 
218  qrm_mat%adata%nnodes = nnodes
219  qrm_mat%adata%ncsing = ncsing
220  qrm_mat%adata%nrsing = nrsing
221 
222 
223  return
224 
225 
226 end subroutine dqrm_attach_singletons
This module contains the definition of the analysis data type.
This module contains the definition of the basic sparse matrix type and of the associated methods...
subroutine dqrm_attach_singletons(qrm_mat, scol, srow, mrperm, mcperm, nrsing, ncsing)
This subroutine merges the results of the singletons detection into the results of the analysis phase...
This type defines the data structure used to store a matrix.
int i
Definition: secs.c:40